The TESS-Keck Survey. XVIII. A Sub-Neptune and Spurious Long-period Signal in the TOI-1751 System

We present and confirm TOI-1751 b, a transiting sub-Neptune orbiting a slightly evolved, solar-type, metal-poor star (T eff = 5996 ± 110 K, log(g)=4.2±0.1 , V = 9.3 mag, [Fe/H] = −0.40 ± 0.06 dex) every 37.47 days. We use TESS photometry to measure a planet radius of 2.77−0.07+0.15R⊕ . We also use both Keck/HIRES and APF/Levy radial velocities (RV) to derive a planet mass of 14.5−3.14+3.15M⊕ , and thus a planet density of 3.6 ± 0.9 g cm−3. There is also a long-period (∼400 days) signal that is observed in only the Keck/HIRES data. We conclude that this long-period signal is not planetary in nature and is likely due to the window function of the Keck/HIRES observations. This highlights the role of complementary observations from multiple observatories to identify and exclude aliases in RV data. Finally, we investigate the potential compositions of this planet, including rocky and water-rich solutions, as well as theoretical irradiated ocean models. TOI-1751 b is a warm sub-Neptune with an equilibrium temperature of ∼820 K. As TOI-1751 is a metal-poor star, TOI-1751 b may have formed in a water-enriched formation environment. We thus favor a volatile-rich interior composition for this planet.


INTRODUCTION
The Transiting Exoplanet Survey Satellite (TESS , Ricker et al. 2014) has discovered 415 confirmed exoplanets to date 1 , and has identified thousands of planet candidates.These planets join the extensive population of over 5500 known planets.With this sample, we can search for demographic trends in planet radii, masses, compositions, and occurrence, and begin to probe the mechanisms of planet formation and evolution that govern the planets in our galaxy.
The most commonly detected exoplanets are those smaller than Neptune (1 − 4 R ⊕ , commonly referred to as sub-Neptunes and/or super-Earths 1 ).TESS has discovered over 200 such planets, and ≳ 75% of all confirmed planets with measured radii are smaller than 4 R ⊕ 1 .Population-level trends offer windows into these planets' evolution and formation.
One such feature is the radius gap, a valley in the distribution of planet radii near 1.8 R ⊕ (Fulton et al. 2017).This feature may be explained by two distinct populations: primarily rocky planets with H/He envelopes, and planets with volatile-rich interiors and more voluminous envelopes (e.g.Zeng et al. 2019;Izidoro et al. 2022;Venturini & Helled 2017;Venturini et al. 2020).Lee et al. (2022) suggest the radius gap is imprinted on the planet population at early times through planet for-mation in gas-poor disks.Alternatively, Luque & Pallé (2022) posit that the two populations are distinct in density, rather than radius, around M dwarfs.In this framework, compositional differences between the populations are set by the materials available in their formation environments, suggesting distinct formation locations within their protoplanetary disks.Water-rich planets may form beyond the snow line and then migrate inwards (Léger et al. 2004;Bitsch et al. 2019).On the other hand, rocky planets may form in situ from ice-poor pebbles, and accrete primordial atmospheres from their natal disks (Lee & Chiang 2016;Chiang & Laughlin 2013).
Other mechanisms involving the sculpting of a single underlying population to create the radius gap have also been suggested.These include atmospheric massloss through XUV-driven photoevaporation (Owen & Wu 2013;Rogers et al. 2021), core-powered mass loss (Ginzburg et al. 2018), and atmospheric mass-loss and growth through giant impacts (Wyatt et al. 2020).
In this paper, we present and confirm TOI-1751 b, a sub-Neptune with a period of 37.47 days, orbiting a slightly evolved G0 star.TOI-1751 b is at the upper radius boundary of the sub-Neptune population, where planet occurrence begins to decrease.While much work has been done to investigate the mechanism(s) sculpting the radius gap, there are few models to explain this "occurrence cliff" (Dattilo et al. 2023).By characterizing planets in this regime, we may refine models that encapsulate both the radius gap and occurrence cliff.Furthermore, due to the lower transit probability for planets with longer periods, only about 10% of sub-Neptunes confirmed by TESS have orbital periods greater than 25 d.TOI-1751 b thus resides in a population that is challenging to probe with transit surveys.This target was observed by the TESS -Keck Survey (TKS), a collaboration spanning several institutions that pools time on the Keck-I telescope on Maunakea.It was initially selected for radial velocity (RV) observations under several science cases: searching for distant giants, probing planets across the radius gap, and analyzing the diversity of gaseous envelopes (see Chontos et al. 2022 for a comprehensive of the TKS science cases).In Section 2 we discuss the data collected for this target.In Section 3 we characterize the stellar host.Section 4 presents our photometric and RV analyses and resulting planet parameters.In Section 5 we put this system in context and discuss possible planetary compositions, and conclude in Section 6.

Photometric Observations
The TESS mission observed TOI-1751 (TIC 287080092, HD 146757) for a total of 27 sectors (15, 17-26, 40-41, and 47-59) between 15 August 2019 and 23 December 2022 at 120 s cadence.The Science Processing Operations Center (SPOC, Jenkins et al. 2016) conducted a transit search of Sectors 15, 17, 18 and 19 on 24 January 2020 with an adaptive, noise-compensating matched filter (Jenkins 2002;Jenkins et al. 2010Jenkins et al. , 2020)).This produced a threshold-crossing event (TCE) to which an initial limb-darkened transit model was fit (Li et al. 2019).Diagnostic tests were then conducted to investigate the planetary nature of the signal (Twicken et al. 2018).The transit signal passed all of the diagnostic tests.The TESS Science Office (TSO) reviewed the vetting information and issued an alert on 27 February 2020 (Guerrero et al. 2021).The signal was repeatedly recovered as additional observations were made in sectors 20-26, 40-41 and 47-59.The final transit search located the host star within 5.8 ± 3.6 ′′ of the source of the transit signal using a difference image centroiding test.We note that observations taken in sectors 15 and 17-26 were impacted by a bias in the sky background correction algorithm, which tended to overestimate the sky background flux.However, the impact on the derived planetary radius of TOI-1751 b is below 0.4% in all affected sectors and typically between 0.1% -0.2%, so this is not a dominant error source in our analysis.

Imaging Observations
Eclipsing binaries with small projected separations from a putative planet host star can create a falsepositive transit signal.This effect is particularly important to consider for TESS Objects of Interest (TOIs), due to the TESS mission's larger pixel size (21 ′′ × 21 ′′ ) compared to that of Kepler (4 ′′ ×4 ′′ ).The SPOC transit search of TOI-1751 was able to constrain the location of the transit event to within a TESS pixel (see Section 2.1), suggesting a blended eclipsing binary was not the cause of the transit events.However, the flux from companion stars can lead to underestimated planetary radii, overestimated bulk densities, and erroneous stellar parameters (Ciardi et al. 2015;Furlan et al. 2017).The contaminating flux from a nearby star may also inhibit the detection of shallow transits (Lester et al. 2021).We thus obtained high-resolution imaging observations of TOI-1751 to search for nearby companion stars.

Lucky Imaging
We obtained two observations of TOI-1751 with the AstraLux instrument (Hormuth et al. 2008) installed at the 2.2-m telescope in the Calar Alto Observatory (Almería, Spain) under average weather and atmospheric conditions (seeing around 1 ′′ ) on 22 March and 14 September 2021 with the SDSS z filter.AstraLux uses the lucky imaging technique to obtain thousands of short exposure frames and selects a few percent of these frames with the best Strehl ratio (Strehl 1902).We obtained 62,200 frames with an exposure time of 10 ms each and selected the best 10% for a final effective exposure time of 62.2 s.We used the final stacked image to obtain the contrast curve by using the astrasens code (Lillo-Box et al. 2012, 2014).The result provides a contrast of ∆z = 6.3 mag for separations above 0.5 ′′ on the 22 March 2021.We found no additional sources in the field of view of the instrument (3 ′′ × 3 ′′ ) within these sensitivity limits.

Speckle Imaging
TOI-1751 was observed on 6 June 2020 using the 'Alopeke speckle instrument on the Gemini North 8-m telescope (Scott et al. 2021).'Alopeke provides simultaneous speckle imaging in two bands (562 nm and 832 nm).Three sets of 1000 × 0.06 s exposures were collected and subjected to Fourier analysis in the standard reduction pipeline (see Howell et al. 2011).The Fourier transform of the summed autocorrelation of each set of images is used to make a fringe image of the target, which is then used to reconstruct the image.We find no companions fainter than the target star by 4.58 magnitudes at 562 nm and by 6.73 magnitudes at 832 nm at separations of 0.5 ′′ (i.e.57 AU) or greater.

Adaptive Optics Imaging
We used the Shane Adaptive optics infraRed Camera-Spectrograph (ShARCS, Kupke et al. 2012;Gavel et al. 2014) mounted on the 3-m Shane Telescope at Lick Observatory to collect AO imaging of TOI-1751.We observed the target using the K s and J band filters.We conducted observations using a four-point dither pattern with a spacing of 4 ′′ on each side.We analyzed the data using the Stellar Image Maturation via Efficient Reduction (SImMER) package (Hirsch et al. 2019;Savel et al. 2020).We find no stellar companions within 5 magnitudes at separations ≥ 1 ′′ in K s band.

Reconnaissance Spectra
We obtained one reconnaissance spectrum with an exposure time of 750 s on 2 March 2020 using the Tillinghast Reflector Echelle Spectrograph (TRES, Fűrész 2008) as part of the TESS Follow-up Observing Program (TFOP) SubGroup 2 (SG2) Reconnaissance Spectroscopy program.TRES is a fiber-fed echelle spectrograph on the 1.5 m Tillinghast Reflector at the Fred Lawrence Whipple Observatory (FLWO) in Arizona, USA, operating between 390 and 910 nm.The spectrograph has a resolving power of R ∼ 44,000.We used the Stellar Parameter Classification (SPC, Buchhave et al. 2012) tool to derive stellar parameters.SPC cross-correlates a ∼ 310 Å region of the observed spectrum surrounding the Mg b lines against a library grid of synthetic spectra calculated using the Kurucz (1992) atmospheric models to derive effective temperature (T eff ), surface gravity, (log(g)), rotational velocity, (vsin(i)), and metallicity, ([m/H]).Metallicity is derived using all available metal lines and is therefore reported as [m/H].
We also obtained follow-up spectra for stellar classification with Las Cumbres Observatory's Network of Robotic Echelle Spectrographs (LCO/NRES, Siverd et al. 2018).NRES is a network of four identical spectrographs located at different observatories spanning a wide longitudinal range.Each instrument is a fiber-fed echelle spectrograph operating between 380-860 nm with a resolution of R ∼ 53,000.Specifically, we obtained spectra for TOI-1751 at our NRES facility at the Wise Observatory, Israel, on 4 March 2020 at 23:00 UTC.The exposure time was 1800 s and the final SNR was 36 at 5500 Å.The wavelength calibrated spectrum was obtained through the standard BANZAI-NRES pipeline (Brandt et al. 2020).The resulting stellar parameters from these observations are reported in Table 1.

Keck/HIRES
The HIgh Resolution Echelle Spectrometer (HIRES) on the Keck-I telescope on Maunakea operates between 360 and 800 nm.We used Keck/HIRES to collect highresolution spectra of TOI-1751 in order to derive precision RVs.We collected 71 spectra between July 2020 and July 2023.We observed TOI-1751 using the red cross-disperser, B5 decker (3.5 ′′ × 0.861 ′′ , R = 50,000) or C2 decker (14 ′′ × 0.861 ′′ , R = 50,000), and a median exposure time of 382 s.We took 69 RV observations with a warm (50 • C) iodine cell in the light path for wavelength calibrations as per Butler et al. (1996).We also took two further higher-resolution spectra without the iodine cell in the light path ("iodine-out") in March and September 2020 in order to obtain a spectral template, one using the B3 decker (14 ′′ × 0.574 ′′ , R = 72,000) and one using the C2 decker.The spectra were reduced using the standard procedures described in Howard et al. (2010).The RVs, RV errors, and Mount Wilson S-Index (a proxy for stellar activity derived from Ca II H & K lines) are reported in Table 4, and the RV data are shown in Fig. 1.We also derived stellar parameters for TOI-1751 using these data (see Section 3, Table 1).

APF/Levy
The Automated Planet Finder (Radovan et al. 2014) is a robotic 2.4-m telescope that hosts the Levy spectrograph, a high-resolution (R=100,000) echelle spectrograph that operates between 374 and 980 nm (Vogt et al. 2014a).We obtained 129 observations of TOI-1751 using APF/Levy between April 2020 and January 2023, with multiple observations per night on several nights.These observations had a median exposure time of 2100 s.Similarly to the Keck/HIRES observations, spectra were taken through a warm iodine cell for wavelength calibration, and RVs were derived using the methods described in Howard et al. (2010).As with other faint targets on APF (V ∼ 9), we use the Keck template splined onto the APF wavelength solution to calculate precise RVs.We report the APF/Levy RVs in Table 4, and show the data in Fig 1 .3. STELLAR PARAMETERS

Spectroscopically-derived Parameters
We used the SpecMatch-Emp2 algorithm (Yee et al. 2017) to characterize TOI-1751.SpecMatch-Emp uses a library of high resolution (R∼55,000) and high signalto-noise (> 100) Keck/HIRES spectra to measure stellar parameters from optical spectra.This method achieves accuracies of 100 K in effective temperature, 15% in stellar radius, and 0.09 dex in metallicity for FGK stars.The algorithm takes the observed stellar spectrum, shifts it to the rest frame wavelength scale, identifies the most similar library spectra, and interpolates between them to derive parameters for the target star.We measure a stellar mass of 0.90 ± 0.03M ⊙ , a stellar radius of 1.17 ± 0.18R ⊙ , and an effective temperature of 5996 ± 110 K.We also used the SpecMatch-Synthetic 3 code (Petigura et al. 2017) to measure stellar parameters from the iodine-out Keck/HIRES spectrum by interpolating between a grid of model spectra from Coelho et al. (2005).
Additionally, we report stellar parameters derived from reconnaissance spectra (described in Section 2.3.1).However, we do not adopt these parameters in our analysis due to the the higher resolution of Keck/HIRES.Some of these parameters (e.g.vsin(i)) are discrepant to greater than 2σ from those measured using Keck/HIRES 3 github.com/petigura/specmatch-syndata, likely due to the low-resolution and low-cadence (often 1 observation per instrument) nature of these observations.Our adopted parameters are consistent to 1σ with those derived homogeneously for the TKS sample in MacDougall et al. (2023).
We find that TOI-1751 is a slightly evolved (log g = 4.24 ± 0.10), metal-poor ([Fe/H] = −0.40 ± 0.06), solartype star.The spectroscopic v sin i and R ⋆ measured using SpecMatch-Emp imply a projected rotation period for the star of P rot / sin i = 46 ± 37 d.We derive additional estimates of projected rotation period from the logR ′ HK observations, and using the SED of TOI-1751 (see Section 3.2).We report all calculated stellar parameters in Table 1.

Photometrically-derived Parameters
We also performed an analysis of the broadband spectral energy distribution (SED) of the star together with the Gaia DR3 parallax (Gaia Collaboration et al. 2023), in order to determine an empirical measurement of the stellar radius (Stassun & Torres 2016;Stassun et al. 2017Stassun et al. , 2018)).We pulled the JHK S magnitudes from 2MASS (Skrutskie et al. 2006), the W1-W4 magnitudes from WISE (Wright et al. 2010), the G BP G RP magnitudes from Gaia, and the FUV and NUV magnitudes from GALEX (Martin et al. 2005).Together, the available photometry spans the full stellar SED over the wavelength range 0.2-20 µm.
We note that the stellar mass derived from the spectroscopic data (0.90 ± 0.03M ⊙ ) is inconsistent to 1.4 σ with that reported in the TESS Input Catalog (1.15 ± 0.17M ⊙ , TIC v8.2, Paegert et al. 2022), possibly due to the slightly evolved nature of the star.
We also do not see any evidence of stellar-rotationrelated signals in the Lomb-Scargle periodogram (Lomb 1976;Scargle 1982) of the TESS light curve.The strongest power in the periodogram is at 4.46 d, but given the measured v sin(i) and R * of this target (1.27 ± 1.00 km/s, 1.17 ± 0.18R * ), this is unlikely to correspond to stellar rotation, and may be due to TESS momentum dumps.We attribute the lack of photometric rotation signals to the known challenges with identifying rotation periods longer than the TESS orbital period of 13.7 d, or roughly half a sector length (Canto Martins et al. 2020;Fetherolf et al. 2023).

Constraining Stellar Age and Companions
To further investigate the rotation period and age of TOI-1751, we use the BANYAN Σ (Bayesian Analysis for Nearby Young AssociatioNs Σ, Gagné et al. 2018) analysis tool to investigate whether TOI-1751 is a member of any of the 27 known young stellar associations within 150 pc.We find that TOI-1751 is highly likely (99.9%) to not be a member of these associations, pointing to an older age.This is in line with the reported surface gravity, suggesting that TOI-1751 is nearing the end of its time on the main sequence.
We use the empirical activity-age relations of Mamajek & Hillenbrand (2008) along with the log R ′ HK of TOI-1751 to find an age of 10.0±1.5 Gyr.This includes both the measurement uncertainty on log R ′ HK and the systematic error reported for the Mamajek & Hillenbrand (2008) relation.The same empirical relations of Mamajek & Hillenbrand (2008) predict a rotation period from the log R ′ HK of 23.3 ± 3.2d.This is consistent with the P rot / sin i value determined using the Torres et al. (2010) empirical relations.
Additionally, Gaia DR3 (Gaia Collaboration et al. 2023) does not contain any proper motion companions within 100 ′′ of TOI-1751.We also used tpfplotter (Lillo-Box et al. 2014) to search for sources in GAIA DR3 within the pipeline aperture mask used to generate the light curve.We found that there were no contaminating sources up to 8 magnitudes fainter than TOI-1751 within the aperture mask, and thus likely no substantial contaminating flux.Furthermore, GAIA DR3 reports a Renormalized Unit Weight Error (RUWE) of 0.90, indicating that this target does not have a detectable companion (Lindgren 2018).

Galactic Context
TOI-1751 is a high proper motion star, and we find that is also metal-poor ([Fe/H] = −0.40 ± 0.06).These attributes indicate that it may be a member of the thick disk of the Milky Way (Bensby 2004;Bensby & Feltzing 2010).Thick disk stars are kinematically hotter, more depleted in metals, more enriched in alpha-elements, and older than those in the thin disk.The formation of the thick disk is still unclear (van der Kruit & Freeman 2011), with possible mechanisms including the merger of the Milky Way with a dwarf galaxy (Quinn et al. 1993, e.g.), and radial mixing of gas and stars (e.g.Schönrich & Binney 2009;Loebman et al. 2011).Planets orbiting thick-disk stars (e.g.Kepler-444, Campante et al. 2015) show that planet formation has occurred for ≳ 11 Gyr.Carrillo et al. (2020) calculated the galactic velocity of TIC stars, and report 1 2 ≈ 110km/s for this target.This result is consistent with that calculated using the methods in Rodriguez (2016), and may suggest thick-disk membership (Nissen 2004).TOI-1751 is, however, still consistent with the low-metallicity tail of the thin-disk distribution, and is 1.6 times more likely to belong to the thin disk than the thick disk (Carrillo et al. 2020).Therefore, we do not conclusively report its thick disk membership, and defer a detailed discussion of the star's alpha abundances to a future work (Polanski et al., in prep).We used the lightkurve package (Lightkurve Collaboration et al. 2018) to download the 2-minute SPOC TESS light curves for TOI-1751.We normalized and stitched together the light curves from each TESS Sector (see Fig. 1) and computed a box least-squares (BLS, Kovács et al. 2002) periodogram.We recovered a significant periodic signal at 37.468 d with an associated Signal Detection Efficiency (SDE) of 24.8.This signal corresponds to the planet candidate TOI-1751.01.
We visually inspected the light curve to confirm the time of the first transit of TOI-1751.01.Next, we reduced the computation time for the subsequent transit fit using the BLS period to trim data points falling outside of a two-day window on either side of each transit center.We also removed any observations flagged with a quality flag greater than 0 to exclude scattered light, cosmic rays, and additional anomalous events.
Next, we flattened each transit individually by first fitting a second-order polynomial to the out-of-transit baseline flux spanning 12 to 48 hours before and after each transit midpoint.We then divided out the best-fit polynomial from the Presearch Data Conditioning Simple Aperture Photometry (PDCSAP, Smith et al. 2012;Stumpe et al. 2012Stumpe et al. , 2014) ) flux within the full four-day window for each transit.To ensure robust flattening, we required each four-day window to be at least 80% complete (i.e., at least 80% of possible cadences contained data).This resulted in 13 complete flattened transits.
To constrain the size and orbital properties of TOI-1751 b, we modeled the transit photometry using the exoplanet package (Foreman-Mackey et al. 2021).To remain agnostic to the stellar properties, we defined a transit model in terms of the planet-to-star radius ratio (R p /R * ), time of first transit (T 0 ), orbital period (P), semi-major axis in units of stellar radii (a/R * ), impact parameter (b), and quadratic limb-darkening coefficients (q 1 , q 2 , using the parameterization of Kipping 2013).Our model also included a mean baseline flux term (⟨F ⟩) and a photometric jitter term (s), which was added in quadrature to the reported flux uncertainties.Finally, we assumed a circular orbit for the planet, setting the eccentricity equal to zero.
We optimized the model parameters using Bayesian inference, implementing a Hamiltonian Monte Carlo (HMC) No U-Turn Sampler (NUTS, Hoffman & Gelman 2011) with PyMC3 (Salvatier et al. 2016) to sample the posterior probability distributions.The prior distributions we selected for each parameter are given in Table 2.We set the target acceptance rate to 0.95 (to account for the higher acceptance fractions returned by HMC samplers compared with Metropolis-Hastings samplers) and initialized the sampler by adapting a dense mass matrix from the sample covariances.We then ran the sampler using a total of 4 chains, each one drawing 20,000 samples after discarding 5,000 burn-in steps.To check for convergence, we computed the Gelman-Rubin Diagnostic and visually inspected the sampler trace plot for each parameter.The median and 68% confidence range for each parameter are given in Table 2.

RV Analysis
We analyze the combined Keck/HIRES and APF/Levy RV observations in order to confirm TOI-1751.01.The Lomb-Scargle periodogram of the combined Keck/HIRES and APF/Levy RV data (see Fig. 2) shows a significant peak with < 0.1% false alarm probability at 37.4 d, which we confirm as the exoplanet TOI-1751 b.We also note a peak at ∼ 400 d.This second, longer-period signal may indicate the presence of a non-transiting distant giant planet in the system.
We used the RadVel package (Fulton et al. 2018) which fits Keplerian models using maximum posterior probability optimization to model the RV data.In our models, we allowed several combinations of the orbital period (P), time of conjunction (t c ), mass (M pl ), argument of periastron (ω), and eccentricity (e) of the planet(s) to vary.We performed several fits: including or excluding a linear trend ( γ), including or excluding a second planet at ∼ 400 d, and circular or eccentric orbits.In all cases, we fixed the orbital period (P b ) and time of conjunction (t c,b ) using the precise constraints from our initial photometric analysis (see Section 4.1).
We find that the mass of TOI-1751 b is consistent to 1σ across all models.Furthermore, the mass is consistent to 1σ when fitting the Keck/HIRES and APF/Levy data separately.The preferred model using the full data set is a circular 1-planet model, which returns a minimum mass of 14.5 +3.15 −3.14 M ⊕ for TOI-1751 b.The ∆BIC between the circular model and eccentric model is less than 1, so there is no clear evidence for eccentricity (Kass & Raftery 1995).We report the circular posteriors in Table 2.

Evidence Against a Distant Giant Planet
There is much interest in probing the outer regions of planetary systems, and investigating the occurrence of small, inner planets that may depend on the presence of distant giants (e.g.Van Zandt et al. 2023;Zhu & Wu 2018).We find that the ∼ 400 d periodic signal seen in the combined data set is only present in the Keck/HIRES RVs, and is not seen in the APF/Levy data (see Fig. 2).This casts doubt on a planetary interpretation of the 400 d signal.
In order to investigate whether a distant giant planet would be detected in the APF/Levy data, we simulated observations using parameters for TOI-1751 b and the putative outer planet determined from a two-planet Keck/HIRES-only fit.We used M = 14.9 ± 3.8M ⊕ for the inner planet and M = 88±15M ⊕ (∼0.9 M Saturn ) for the outer planet to simulate 129 observations, using real APF/Levy observational time stamps and errors.All of the simulated observations show a long-period peak greater than that in the real APF/Levy data (see Fig. 3).This suggests that if such a distant giant planet were present in the system, its signal would have been captured by the APF observations.We performed a similar exercise for simulated APF data assuming lower masses (10, 35, and 50 M ⊕ ) for a putative distant giant.None of these scenarios reproduced the periodogram of the APF data.We used Radvel to model these sets of simulated observations.The measured planet mass for the outer planet is inconsistent with the injected values.As such, we did not find a planetary solution for the 400 d signal that is consistent with observations of TOI-1751 b, and conclude that this long-period signal is not planetary in nature.
A possible cause of this spurious signal is that the stellar template used in reducing the RV observations was sub-optimal.This would cause a correlation between the barycentric velocity of the observations and the RV measurements.However, we would expect such an effect to impact the APF/Levy observations as well.We also computed the spectral window function of our observations, i.e. the Lomb-Scargle periodogram of observation times.There is substantial power at ∼ 200 d in the window function of the Keck/HIRES observations, indicating that the periodic signal in the RV data may be a harmonic of a signal in the window function due to patterns in observation times.Given that the ∼ 400 d periodic signal is only present in the Keck/HIRES data, that it would be detectable in the APF/Levy data assuming it was astrophysical, and that there is a corresponding peak in the Keck/HIRES window function, we conclude that this signal is not planetary in nature.This highlights the valuable role that the APF, and other comparable telescopes, can play in large RV surveys.High-cadence APF/Levy observations of TOI-1751, despite being lower precision than the Keck/HIRES data, allowed us to identify and remove an alias in our data, and to gain a clearer picture of this system.Using the results from the previous two sections as a baseline, we performed a joint fit of the 13 flattened TESS transits described in Section 4.1 and the APF and HIRES RV measurements described in Section 4.2.

Joint Photometry and RV Analysis
As in Section 4.1, we used the exoplanet package (Foreman-Mackey et al. 2021).We define the combined photometry and RV model by a set of 16 free parameters: R p /R * , T 0 , P , b, q 1 , q 2 , and ⟨F ⟩ (as in the transitonly model), as well as stellar mass and radius, M * and R * , and planet mass M p .We also included three jitter terms (one per instrument), which were added to the respective instrument uncertainties in quadrature, as well as RV offsets, γ, for both APF and HIRES.Finally, we added an RV acceleration term, γ.We also derived the following quantities from the fit parameters: the stellar density (ρ * ), planet radius (R p ), planet bulk density (ρ p ), RV semi-amplitude (K), transit duration (t dur ), and semi-major axis (a).As stated in Section 4.2, the ∆BIC between circular and eccentric is less than 1, indicating no clear evidence for eccentricity, so we adopt the circular model posteriors.We report both circular and eccentric model posteriors in Rable 2 for completeness.
Again, we used the HMC NUTS implemented in PyMC3 (Salvatier et al. 2016) to optimize the model parameters and sample their posterior probability distributions.The priors we adopted for each of these parameters are shown in Table 2.We followed the procedure described in Section 4.1 to initialize the sampler and sample the posteriors.The full posterior distributions for the joint transit and RV fit are shown in Figure 8 in the Appendix, and the posteriors for the auxiliary parameters are shown in Figure 9. Finally, phase-folded plots of the data and models are shown in Figures 4   and 5.The median posterior values for each parameter, along with 68% confidence intervals, are shown in Table 2.We confirm the planetary nature of TOI-1751 b, and measure its radius (2.77 +0.15  −0.07 R ⊕ ), mass (14.5 +3.15 −3.14 M ⊕ ), and bulk density (3.6 ± 0.9 g cm −3 ).TOI-1751 b has a longer orbital period (37.47d) than ∼ 90% of all confirmed planets with measured radii consistent to 1σ with TOI-1751 b and with well-constrained masses ( M pl σ M pl > 3).Additionally, it has a larger mass than 94% of this sample, (see Fig. 6).It is also one of only 5 sub-Neptunes with periods longer than 30 days orbiting bright (V < 10) stars, representing one of the best cases for investigating the warm sub-Neptune population.
Although TOI-1751 is likely a thin-disk star, metalpoor Neptune-hosting stars are typically enriched in alpha-elements (Adibekyan et al. 2012).This suggests that planets orbiting these stars had fewer metals, and more water, available when forming.We thus favour a volatile-rich interior composition for TOI-1751 b (see Section 5.2).

Interior Composition
The number of exoplanets with radii greater than ∼ 3R ⊕ declines sharply, with planets between 2.7 and 3.0 R ⊕ around FGK stars within 100 d being up to 10 times more common than planets between 3.3 and 3.7 R ⊕ (Fulton & Petigura 2018;Hsu et al. 2019).Currently, the only proposed mechanism for this phenomenon is the "fugacity crisis": the increased solubility of hydrogen in magma at high pressures (Kite et al. 2019).At ∼ 3R ⊕ , the pressures at the base of a planetary atmosphere are sufficient to sequester atmospheric hydrogen in magma.The associated loss of H 2 (and thus, volume) from the atmosphere during formation may thus decrease the occurrence of ≳ 3R ⊕ planets.This makes an interesting case for studying the interior composition of planets with a radius just below ∼ 3R ⊕ such as TOI-1751 b.These planets may have sequestered hydrogen from their atmospheres, which may be observable as low atmospheric mass fractions.
The radius of TOI-1751 b (2.77 +0.15 −0.07 R ⊕ ) implies that it is unlikely to be a purely rocky planet (Rogers 2015).Furthermore, the instellation (∼ 75S ⊕ ) and escape velocity (∼ 25km/s) of TOI-1751 b place it firmly on the right of the cosmic shoreline, implying that this planet has an atmosphere (Zahnle & Catling 2017).
However, mass and radius measurements alone are not sufficient to uniquely constrain planet composition (Adams et al. 2008;Rogers & Seager 2010).The planet's bulk density is consistent with both a rocky interior with a few percent by mass H/He atmosphere and a volatilerich interior with a H/He atmosphere (see Fig. 6).This degeneracy between silicate-and water-rich models is compounded by the unconstrained albedo and heat redistribution patterns, and is an ongoing challenge in characterizing sub-Neptunes (e.g.Valencia et al. 2013;Nixon & Madhusudhan 2021;Luque & Pallé 2022).
As shown in Fig. 6, if TOI-1751 b has a rocky interior, it will have a higher H/He atmospheric mass fraction compared to if it has a volatile-rich interior.A volatile-rich planet has a lower bulk density than that of a planet with an iron-rich core, and thus needs a less massive atmosphere to make up the bulk density.The models described in Zeng et al. (2019) are also consistent with those described in Lopez & Fortney (2013), which predict a 0.5-3% H/He atmosphere by mass for TOI-1751 b.
The relatively low atmospheric mass fraction of TOI-1751 b may have been sculpted by mass loss during evolution.Ongoing mass loss may be probed through future spectroscopic observations of Lyman α or meta-stable helium (e.g.Spake et al. 2018;Kulow et al. 2014), although the planet's long transit duration will make these observations resource-intensive.Although this planet has a low Transmission Spectroscopy Metric (TSM = 20, Kempton et al. 2018), it may yet prove an interesting target for comparative atmospheric studies with next-generation instruments.
The mass, radius, and instellation of TOI-1751 b are also consistent with models of putative irradiated ocean worlds (Aguichine et al. 2021).These models are composed of refractory layers (iron core and rocky mantle), a hydrosphere with an equation of state that extends to the plasma regime, and a steam atmosphere.Within this framework, TOI-1751 b is consistent with a ∼ 20 − 50% water mass fraction for irradiation temperatures between 600K and 1000K (T irr = 700K for TOI-1751 b, see Eqn. 9 of Aguichine et al. 2021).The instellation of TOI-1751 b at its current location also lies above the water vapour runaway greenhouse threshold (Zahnle & Catling 2017), suggesting that any water in this planet's atmosphere may exist as steam.
We used MAGRATHEA4 (Huang et al. 2022)    layers: an iron core, a silicate mantle, a steam/water/ice hydrosphere, and a H/He atmosphere.We note that these models do not encapsulate the true complexity of planet interiors, which may also include mixing between rock and water layers (Kovačević et al. 2022;Vazan et al. 2022).Ice and rock may remain mixed in planetary interiors for billions of years if no significant mass loss occurs.The timescale for a planet to lose its atmosphere via XUV-driven photoevaporation scales with M 2 p (or, more precisely, M 2 core ) and scales inversely with instellation (see Eqn. 4 in Lopez & Fortney 2013).A massive sub-Neptune-sized planet such as TOI-1751 b may have been able to resist photoevaporative massloss and thereby retain a small but non-negligible H/He atmosphere, and thus, an ice-rock mixture, in its interior over its lifetime.
We used two MAGRATHEA modes to model TOI-1751 b: the default model and a composition solver, which we discuss in Sections 5.2.1 and 5.2.2, respectively.

Solving for Planet Radius
Given the inherent degeneracies present in modelling sub-Neptune interiors, we consider several compositional "families" of solutions based on Solar System objects and confirmed exoplanets in our analysis, differing in their core, mantle, water, and atmospheric mass fractions (CMF, MMF, WMF, and AMF, respectively).We summarize the families used in Table 3.
The default MAGRATHEA mode calculates the planet radius for a given planet mass and distribution of mass between layers.For each step in enclosed mass, the pressure, density, and temperature are calculated, moving from M = 0 to M = M p .This mode takes as input Models of TOI-1751 b assuming an Earth-like rock/iron ratio, no atmosphere (i.e., AMF = 0%), a supercritical and condensed water layer, and a water envelope.By allowing the planet's envelope to be made of water, the WMF required to reproduce the planet's radius is substantially reduced.Grey x's mark mass and radius samples which are too low density to be modeled with AMF=0%.Right: Models of TOI-1751 b assuming an Earth-like mantle/core ratio, a WMF of 20%, and a non-ideal H/He atmosphere.The non-ideal atmosphere is less compressible than an ideal H/He envelope, and thus a smaller AMF is needed to reproduce the measured planet radius.Iso-composition curves for WMF and AMF are shown with dashed lines.
the temperature discontinuities at each boundary layer (if any), as well as the surface temperature.We assume thermal equilibrium between each boundary layer and a surface temperature equal to the calculated equilibrium temperature (820 K).We use the default phase diagrams and equations of state for the iron core and magnesium silicate mantle.We use an updated version of MAGRATHEA which includes tabulated equations of state for water in liquid, solid, vapor, and supercritical forms in the hydrosphere (Haldemann et al. (2020)) and for a H/He atmosphere with solar composition (Chabrier & Debras (2021)).
We find that Uranus-like, Neptune-like, and Earth-like compositions do not reproduce the measured radius of TOI-1751 b.We find that a WMF between that of a Kepler-11c-like and Ganymede-like model most closely reproduced the measured radius of TOI-1751 b, though these objects are very different in scale and environment.
We also manually refined the mass distribution between layers in order to optimize the model.In this, we set an atmospheric mass fraction of 1%, based on other models of planet composition (Zeng et al. 2019;Lopez & Fortney 2013).The model which most closely reproduces the radius of TOI-1751 b includes a CMF of 36%, MMF of 53%, WMF of 10%, and an AMF of 1%.However, we stress that this is a possible composition rather than a prescriptive model of the interior of this planet.

Solving for Mass Fractions
The second MAGRATHEA input method used, a composition solver, calculates a planet's atmospheric and water mass fractions given a fixed planet mass and radius.Again, we assumed thermal equilibrium between each boundary layer and calculated the null albedo surface temperature to be 690 K.We drew 1000 samples from the planet mass and radius posteriors measured using RV and photometric data, and used these as inputs for the composition solver.Our base model used an interior that has an Earth-like core-to-mantle ratio (32.5% core to 67.5% mantle).We then found the WMF (including a water-vapor atmosphere) required to match the radius of TOI-1751 b.As shown in Fig. 7, we found WMFs from 21 to 100% with a median WMF of 49%.This median WMF may suggest the planet formed beyond the snow line, much further out than its present location.
However, if we include a H/He atmosphere above supercritical water, the planet is consistent with a lower WMF.We now assumed an Earth-like core to mantle ratio and a WMF of 20%, set by the lowest WMF found in the previous section.We then redistributed mass from the interior to form a H/He atmosphere.The atmosphere has an isothermal temperature profile with temperature increasing adiabatically at pressures above 100 bar (similar to the models of Nixon & Madhusudhan 2021).The AMF under this model ranges from 0.04-3.65%with a median of 0.72% (see Fig. 7), much smaller than the atmospheric mass fractions of Neptune and Uranus.These findings are consistent with Section 5.2.1 and demonstrate the range of possible interior solutions for TOI-1751 b.

CONCLUSIONS
In this work, we study the solar-type star TOI-1751 using TESS photometry, and Keck/HIRES and APF/Levy RV observations.We confirm the planetary nature of TOI-1751 b, a transiting sub-Neptune, and measure its radius (2.77 +0.15  −0.07 R ⊕ ) and mass (14.5 +3.15 −3.14 M ⊕ ).From these measurements, we infer a bulk density of 3.6 ± 0.9 g cm −3 .This points to several possible compositions: a rocky planet with a H/He atmosphere, a sub-Neptune with a volatile-rich interior, or an irradiated ocean planet with a volatile-rich interior and steam atmosphere.
TOI-1751 is a metal-poor ([Fe/H] = −0.40 ± 0.06) star, and as such, may have had a water-enriched protoplanetary disk (Adibekyan et al. 2012).We thus favor a volatile-rich (i.e.metal-poor) interior for this planet.We find its bulk density (3.6±0.9 g cm −3 ) can be matched by a volatile-rich interior with an atmosphere mass fraction of ∼ 1%.TOI-1751 b is a relatively longperiod (37.47 d) but highly irradiated (T eq = 820K) planet, giving us an insight into the small warm Neptune population with precisely measured masses and radii.
Facility: Keck:I (HIRES, Vogt et al. 1994), APF (Vogt et al. 2014b), TESS (Ricker et al. 2014) Software: This research made use of exoplanet (Foreman-Mackey et al. 2021;Foreman-Mackey et al. 2021) and its dependencies (Agol et al. 2020;Kumar et al. 2019;Astropy Collaboration et al. 2013, 2018, 2022;Salvatier et al. 2016;Theano Development Team 2016).This research also made use of Lightkurve, a Python package for Kepler and TESS data analysis (Lightkurve Collaboration et al. 2018), and its dependencies (Astropy Collaboration et al. 2013Collaboration et al. , 2018;;Ginsburg et al. 2019;Brasseur et al. 2019).Additionally, this work made use of the following Python packages: numpy (Harris et al. 2020), pandas (pandas development team 2020), matplotlib (Hunter 2007), corner (Foreman-Mackey 2016), RadVel (Fulton et al. 2018), SpecMatch-Emp (Yee et al. 2017), SpecMatch-Synth (Petigura et al. 2017), and CMasher (Van der Velden 2020).Some of the data presented in this paper were obtained at the W. M. Keck Observatory, which is operated as a scientific partnership among the California Institute of Technology, the University of California and the National Aeronautics and Space Administration.The Observatory was made possible by the generous financial support of the W. M. Keck Foundation.The authors wish to recognize and acknowledge the very significant cultural role and reverence that the summit of Maunakea has always had within the indigenous Hawaiian community.We are most fortunate to have the opportunity to conduct observations from this sacred mountain which is now colonized land.
This paper made use of data collected by the TESS mission and are publicly available from the Mikulski Archive for Space Telescopes (MAST) operated by the Space Telescope Science Institute (STScI).Funding for the TESS mission is provided by NASA's Science Mission Directorate.We acknowledge the use of public TESS data from pipelines at the TESS Science Office and at the TESS Science Processing Operations Center.Resources supporting this work were provided by the NASA High-End Computing (HEC) Program through the NASA Advanced Supercomputing (NAS) Division at Ames Research Center for the production of the SPOC data products.
We thank Artem Aguichine for useful discussions on interpreting composition models, and Andreia Carrillo for insight on the galactic context of TOI-1751.We also acknowledge Erik Petigura and BJ Fulton for contributions to the construction, design, and team management of the TESS -Keck Survey.
M. Rice acknowledges support from Heising-Simons Foundation Grant.

Figure 1 .
Figure 1.Top: RV measurements from APF/Levy (green squares) and Keck/HIRES (blue diamonds).For clarity, the RV offsets determined from the joint RV and transit fit (see Section 4.3) for each instrument have been subtracted from the data, and the respective jitter terms have been combined with the instrumental uncertainties to produce the error bars.Bottom: Flattened and normalized TESS PDCSAP flux light curve.

Figure 2 .
Figure 2. Lomb-Scargle periodogram of a) all RV observations, b) Keck/HIRES RVs, c) APF/Levy RVs, d) Keck/HIRES window function, e) APF/Levy window function,and f) S-values.The vertical green line marks the 37d period, and the gray-shaded region marks the 300-500d range.The long-period peak seen in the periodogram of the full RV data set is not present in the APF data, and is at a period similar to that of a significant peak in the HIRES window function.The dark red dashed line represents the 0.1% false-alarm level.Data consisting of Gaussian noise with no periodic signal would produce peaks of this height (or above) in < 0.1% of samples.The red dashed line is the 1% false-alarm level.

Figure 3 .
Figure 3. Lomb-Scargle periodogram of APF/Levy RV data (blue), and 100 simulations of APF/Levy observations of a 2-planet TOI-1751 system.The period of TOI-1751 b is indicated by a vertical green line.The peaks in the periodogram of the simulated data are consistently higher than those of the real data at ∼ 400 d.

Figure 4 .
Figure 4. Top: Phase-folded TESS light curve of TOI-1751 b.The black points with error bars show the flattened PDC-SAP flux from 13 transits, where the error bars are calculated by adding the jitter and instrument uncertainty in quadrature.The orange line shows the best-fit (median) transit model from the joint RV and transit analysis.For clarity, the data are also shown in 1-hour bins (white circles).Bottom: Residual flux between the best-fit (median) transit model and the data.

Figure 5 .
Figure 5. Top: RV measurements from APF (blue diamonds) and HIRES (green squares) and the best-fit (median) RV model (orange line) from the joint RV and transit analysis.For clarity, the RV offset for each instruments has been subtracted from the data, and the respective jitter terms have been combined with the instrumental uncertainties to produce the error bars.The orange shaded region shows the 68% confidence range for the model.Middle: Phase-folded radial velocities of TOI-1751 b, with the instrumental offsets and linear trend removed.The APF data (blue triangles) and HIRES data (green squares) are shown with error bars and the best-fit (median) RV model for the planet (orange line).The jitter term for each instrument has been added to the RV uncertainty in quadrature to calculate the error bars.The orange shaded region indicates the model 68% confidence range.For clarity, repeated orbital phase coverage is shown in the gray shaded regions, and we have also binned the data to ∼4 d intervals (white circles with error bars).Bottom: Residual radial velocities between the best-fit (median) model and the data.

2 .
are shown with the 68% confidence interval.** The quadratic limb darkening coefficients are implemented using the Kipping (2013) parameterization.† The RV trend coefficients are for a linear function: vr(t) = γt + γ. ‡ Each jitter term was added to the respective data uncertainties in quadrature, such that the effective error was √ σ 2 + s Note-U (a, b) is a uniform distribution that is nonzero only between a and b; N (µ, σ) is a normal distribution with mean µ and standard deviation σ; T (µ, σ, a, b) is a truncated normal distribution with mean µ and standard deviation σ, bounded between a and b; and L(µ, σ) is a log-normal distribution with log mean µ and log standard deviation σ.
Compositional "families" of planet interiors that we used as inputs for the default MAGRATHEA mode to model TOI-1751 b in Section 5.2.1.For each model we include the core mass fraction (CMF), mantle mass fraction (MMF), water mass fraction (WMF), and atmospheric mass fraction (AMF).

Figure 6 .
Figure 6.Mass vs. radius diagram showing confirmed planets with masses and radii measured to > 3σ precision (black points) and TOI-1751 b (orange).Each panel includes different theoretical composition curves fromZeng et al. (2019) .Earth (E), Venus (V), Uranus (U), and Neptune (N) are also shown for context, with the precise masses and radii for these planets lying in the center of the letter symbol.TOI-1751 b is consistent with both rocky and water-rich models.Left: The Earth-like composition is assumed to be 32.5% Fe/Ni-metal plus 67.5% MgSiO3-rock.Curves for 100% H2O, and 100% MgSiO3 are also shown.Middle: Composition curves assuming an Earth-like planet with the addition of a H2/He atmosphere (made up of a mixture of 75% H2 and 25% He).These curves are evaluated along interior adiabats at different internal specific entropies, labelled by the temperature of the corresponding specific entropy at 100-bar level in the gas envelope (2000K is sparsely dashed, 1000K is solid, 700K is dashed, 500K is dotted, 300K is thin solid).Right: Composition curves correspond to an Earth-like planet with an H/He isothermal envelope at various surface temperatures atop an ice-VII (a cubic crystalline form of ice) layer.

Figure 7 .
Figure 7. Mass vs. radius posterior samples from the joint transit and RV fit for TOI-1751 b modeled with MAGRATHEA.Left:Models of TOI-1751 b assuming an Earth-like rock/iron ratio, no atmosphere (i.e., AMF = 0%), a supercritical and condensed water layer, and a water envelope.By allowing the planet's envelope to be made of water, the WMF required to reproduce the planet's radius is substantially reduced.Grey x's mark mass and radius samples which are too low density to be modeled with AMF=0%.Right: Models of TOI-1751 b assuming an Earth-like mantle/core ratio, a WMF of 20%, and a non-ideal H/He atmosphere.The non-ideal atmosphere is less compressible than an ideal H/He envelope, and thus a smaller AMF is needed to reproduce the measured planet radius.Iso-composition curves for WMF and AMF are shown with dashed lines.
#2023-4478.J. Akana Murphy is supported by the National Science Foundation (NSF) Graduate Research Fellowship Program (GRFP) under Grant No. DGE-1842400.A. Desai would like to acknowledge the UC Berkeley Physics Undergraduate Research Scholars Program, and the UC Berkeley Mathematical & Physical Sciences Scholars for funding their work.E. Turtelboom acknowledges support from a David & Lucile Packard Foundation grant.C. Harada acknowledges support from the National Science Foundation Graduate Research Fellowship Program under Grant No. DGE 2146752.D. Huber acknowledges support from the Alfred P. Sloan Foundation, the National Aeronautics and Space Administration (80NSSC21K0652), and the Australian Research Council (FT200100871).M. Hill would like to acknowledge NASA support via the FINESST Planetary Science Division, NASA award number 80NSSC21K1536.

Table 2 .
Prior distributions and posterior quantiles derived from model fits.

Table 3 .
to investigate the composition of TOI-1751 b.MAGRATHEA is a 1D structure code which assumes fully differentiated planet

Table 4 .
TOI-1751 RV and S Value