TOI-2374 b and TOI-3071 b: two metal-rich sub-Saturns well within the Neptunian desert

We report the discovery of two transiting planets detected by the Transiting Exoplanet Survey Satellite (TESS), TOI-2374 b and TOI-3071 b, orbiting a K5V and an F8V star, respectively, with periods of 4.31 and 1.27 days, respectively. We confirm and characterize these two planets with a variety of ground-based and follow-up observations, including photometry, precise radial velocity monitoring and high-resolution imaging. The planetary and orbital parameters were derived from a joint analysis of the radial velocities and photometric data. We found that the two planets have masses of $(57 \pm 4)$ $M_\oplus$ or $(0.18 \pm 0.01)$ $M_J$, and $(68 \pm 4)$ $M_\oplus$ or $(0.21 \pm 0.01)$ $M_J$, respectively, and they have radii of $(6.8 \pm 0.3)$ $R_\oplus$ or $(0.61 \pm 0.03)$ $R_J$ and $(7.2 \pm 0.5)$ $R_\oplus$ or $(0.64 \pm 0.05)$ $R_J$, respectively. These parameters correspond to sub-Saturns within the Neptunian desert, both planets being hot and highly irradiated, with $T_{\rm eq} \approx 745$ $K$ and $T_{\rm eq} \approx 1812$ $K$, respectively, assuming a Bond albedo of 0.5. TOI-3071 b has the hottest equilibrium temperature of all known planets with masses between $10$ and $300$ $M_\oplus$ and radii less than $1.5$ $R_J$. By applying gas giant evolution models we found that both planets, especially TOI-3071 b, are very metal-rich. This challenges standard formation models which generally predict lower heavy-element masses for planets with similar characteristics. We studied the evolution of the planets' atmospheres under photoevaporation and concluded that both are stable against evaporation due to their large masses and likely high metallicities in their gaseous envelopes.


INTRODUCTION
Since the discovery of 51 Pegasi b (Mayor & Queloz 1995), more than 5500 exoplanets have been confirmed, thanks largely to detections from ground based RV surveys such as HARPS (High Accuracy Radial velocity Planet Searcher; Mayor et al. 2003a) , transiting surveys such as WASP (Wide Angle Search for Planets; Pollacco et al. 2006) and space telescopes such as CoRot (Convection, Rotation, and planetary Transits; Baglin et al. 2008), Kepler (Borucki et al. 2010) and, most recently, the Transiting Exoplanet Survey Satellite (TESS, Ricker et al. 2015).TESS was launched in 2018 with the objective of known as the "hot Neptune desert", "subJovian desert", and "evaporation desert"), which roughly encompasses Neptunian-sized planets (approximately 2 R ⊕ <   < 9 R ⊕ and 0.03 M J <   < 0.8 M J ) with periods up to ∼ 5 days (Szabó & Kiss 2011;Mazeh et al. 2016;Owen & Wu 2017;Owen & Lai 2018;Deeg et al. 2023).This should not be due to an observational bias, as planets with short period and intermediate size are easily discovered by transit surveys like Kepler.Several theories have been put forward to explain the existence of the desert and its boundaries in the parameter space (e.g., Owen & Lai 2018;Vissapragada et al. 2022): the lower boundary could be caused by photoevaporation stripping away their gaseous H/He envelopes, therefore reducing their radii/masses and leaving behind a dense core; while the upper boundary seems to be stable against photoevaporation, and may instead be understood as a "tidal disruption barrier": in order to form planets below and left of the boundary through high-eccentricity migration, they would have to come so close to their host star and no longer be able to succesfully circularize and stabilize.
Many of these planets are readily observable with missions such as TESS and JWST, and also with ground-based spectroscopic instruments such as HARPS (Mayor et al. 2003b) and ESPRESSO (Pepe et al. 2021) due to their close orbits and short periods.The HARPS-NOMADS programme (PI Armstrong, 1108.C-0697) aims to significantly increase the number of planet confirmations in this parameter space with precise masses and radii.With these parameters we can constrain the densities and thus the internal structures of these planets.The few planets found in this desert are likely to have undergone unusual formation and/or evolutionary processes compared to those in more populated parameter spaces so this analysis could lead towards an improved understanding of the formation and evolution of the Neptunian desert.
Residing in the Neptunian desert with masses between 0.1 and 0.4   , hot Saturns are a group whose study can help us further understand the divergent formation pathways of small planets and gas giants: on the one hand, they could be the smallest planets formed via runaway gas accretion.They can therefore provide information on the limits of the core accretion mechanism (Petigura et al. 2018) .On the other hand, in accordance with the trend that associates the presence of short-period planets and of large planets with higher host star metallicity (Mulders et al. 2016;Dong et al. 2018), Petigura et al. (2018) found that stars hosting hot sub-Saturns (planets with radii between 4 and 8  ⊕ and periods between 1 and 10 days) have the highest mean stellar metallicity among all planet hosts.This evidence suggests some kind of mechanism connected to high stellar metallicity that leads to these short-period large planets.Lastly, because their lower surface gravities lead to larger atmospheric scale heights than those of typical hot Jupiters, hot Saturns are some of the best targets for transmission spectroscopy observations (Wakeford et al. 2017).
We present here the detection and characterisation of TOI-2374 b and TOI-3071 b, two hot sub-Saturns transiting a K5V and an F8V star, respectively.The observations leading to the detection and confirmation of the two planets include photometry from TESS and ground-based telescopes, spectroscopy from HARPS and PFS and high-spatial resolution imaging from SOAR.The details of these observations are outlined in Section 2. In Section 3 we present the analysis of these observations, including the determination of stellar parameters and the joint modelling used to constrain the planetary parameters.In Section 4 we present the results and discuss the nature of the two targets, including the position of the planets in the Neptune desert and in mass-radius parameter space.This section also includes an analysis of the internal structure and evaporation history of the planet.Finally, in Section 5 we report the conclusions of our work.

TESS photometry
TESS observed TOI-2374, with TESS Input Catalog (TIC) ID 439366538, in sector 1 (UT 2018 July 25-August 22) with a 30 minute cadence and in sector 28 (UT 2020 July 30-August 26) with a 10 minute cadence, both on camera 1 CCD 4. TOI-2374 b was detected in the Full-frame images (FFIs) by the MIT's Quick-Look Pipeline (QLP, Huang et al. 2020) and alerted as a TESS object of interest (TOI) on 28 October 20202 (Guerrero et al. 2021).The detection gave a period of 4.31362 ± 0.00001 d, an epoch of 2458326.555387days in BJD, and a depth of 9.490 ± 0.008 partsper-thousand (ppt).
TOI-3071 (TIC ID 452006073) was monitored by TESS in sectors 10 (UT 2019 March 26-April 22) and 37 (UT 2021 April 2-28) with camera 2 CCD 2 with a 30 and 10 minute cadence FFI, respectively.The QLP detected the candidate and it was alerted as a TOI on 4 June 20213 .The detection gave a period of 1.2671 ± 0.0003 d, an epoch of 2459331.9673± 0.0017 days in BJD, and a depth of 2.580 ± 0.002 ppt.It was subsequently observed in sector 64 (UT 2023 April 6-May 4) on camera 2 CCD 1 as a 2-minute cadence target.The Science Processing Operations Center (SPOC) at NASA Ames Research Center conducted a transit with a noise-compensating matched filter (Jenkins 2002;Jenkins et al. 2010Jenkins et al. , 2020) ) and detected the transit signature of TOI-3071 b.An initial limb-darkened model fit was performed (Li et al. 2019) and a suite of diagnostic tests were conducted help make or break the planetary nature of the signal (Twicken et al. 2018).The candidate transit signature passed all these diagnostic tests including the difference image centroiding test, which constrained the location of the host star to within 2.1 +/-2.6 arc sec of the transit source.
For observations of TOI-2374 from TESS sectors 1 and 28 and observations of TOI-3071 from sectors 10 and 37, we downloaded the light curves computed by the QLP and used the Kepler Spline Simple Aperture Photometry (KSPSAP) flux, which is detrended by applying a high-pass filter to remove low-frequency variability from stellar activity or instrumental noise (Huang et al. 2020).We then used the lightkurve python package (Lightkurve Collaboration et al. 2018) to remove all data points flagged as being affected by excess noise (corresponding to stray light from Earth or the Moon in the camera field of view or scattered light).For TOI-3071 observations from sector 64, we downloaded the photometry provided by the SPOC pipeline, and used the Presearch Data Conditioning Simple Aperture Photometry (PDCSAP), from which common trends and artefacts, as well as an estimated contamination from blended sources, have been removed by the SPOC Presearch Data Conditioning (PDC) algorithm (Twicken et al. 2010;Smith et al. 2012;Stumpe et al. 2012Stumpe et al. , 2014)).No further detrending of the light curves was deemed necessary as they are relatively flat across the whole time series.The mediannormalised light curves for both targets were used in the joint model outlined in section 3.3.
TIC contamination ratios of 6.65 and 0.62 are reported for TOI-2374 and TOI-3071 respectively.These values are computed as the nominal flux from the contaminants divided by the flux from the source.The large values reported for both targets are most likely due to close sources as can be seen in the TESS Target Pixel Files images shown in Figure 1 (created with tpfplotter4 , Aller et al. 2020).
The TPFs show several sources of potential contamination around the target stars (particularly TOI-2374, which has a very close source two magnitudes brighter).Lightcurves from QLP and SPOC data processing pipelines are approximately corrected for blending from nearby stars.The SPOC pipeline computes a crowding metric using a series of assumptions described in Section 2.3.11 of Stumpe et al. (2012).The QLP uses a different set of assumptions described in step 7 in Huang et al. (2020).To further measure the extent of this blending and contribute to the correction, follow up ground-based photometry was necessary; the details of these observations are described below.Additionally, the TPF for TOI-2374 shows that both the target and the bright companion have the same proper motion direction so they might be bounded companions.

Light Curve Follow-up Observations
The TESS pixel scale is ∼ 21 ′′ pixel −1 and photometric apertures typically extend out to roughly 1 arcminute, generally causing multiple stars to blend in the TESS photometric aperture.This is precisely the case for TOI-2374 and TOI-3071 (Fig. 1).To determine the true source of the TESS detection, we acquired ground-based time-series follow-up photometry of the field around both TOI-2374 and TOI-3071 as part of the TESS Follow-up Observing Program (TFOP; Collins et al. 2018) 5 .The follow-up light curves are also used to confirm the transit depth and thus the TESS photometric deblending factor and refine the TESS ephemeris.We used the TESS Transit Finder, which is a customized version of the Tapir software package (Jensen 2013), to schedule our transit observations.The observations are summarised below and confirm the transit-like events detected by TESS to be occurring within the small TOI-2374 and TOI-3071 follow-up photometric apertures.

Brierfield Private Observatory
We observed one full transit window of TOI-2374 b on UTC 2021 May 18 in B band using the Brierfield Observatory, located near Bowral, New S. Wales, Australia.The 0.36 m telescope is equipped with a 4096×4096 Moravian 16803 camera.The image scale after binning 2×2 is 1.47 ′′ pixel −1 , resulting in a 50 ′ × 50 ′ field of view.The photometric data were extracted using the AstroImageJ (AIJ) software package (Collins et al. 2017).Circular photometric apertures with radius 7.4 ′′ were used as they excluded all of the flux from the nearest known neighbor.

PEST
We observed a egress window of TOI-3071 b in Rc band on UTC 2021 June 04 using the Perth Exoplanet Survey Telescope (PEST), located near Perth, Australia.The 0.3 m telescope is equipped with a 5544×3694 QHY183M camera.Images are binned 2 × 2 in software giving an image scale of 0.7 ′′ pixel −1 resulting in a 32 ′ ×21 ′ field of view.A custom pipeline based on C-Munipack6 was used to calibrate the images and extract the differential photometry.
We used circular photometric apertures with radius 6.4 ′′ .The target star aperture included the flux from the nearest known neighbor in the Gaia DR3 catalog (Gaia DR3 5342880462297608192), which is ∼ 3.8 ′′ east of TOI-3071.The PEST observation was not included in the modeling due to the lack of ingress coverage and low signal to noise detection.

Spectroscopy
In this section we describe in detail the high-precision RV measurements used in this paper.These observations are presented in figures 11 and 12.In section 3.2 we perform a further exploration of this data.

HARPS RVs
We obtained radial velocity (RV) measurements of TOI-2374 and TOI-3071 with the High Accuracy Radial velocity Planet Searcher (HARPS) spectrograph mounted on the ESO 3.6 m telescope at the La Silla Observatory in Chile (Mayor et al. 2003a).Under the HARPS-NOMADS large programme (ID 1108.C-0697, PI: Armstrong), a total of 21 spectra were obtained between 9 April 2022 and 7 July 2022 for TOI-2374 and 14 spectra between 19 March and 2 April 2022 for TOI-3071.The instrument (with resolving power  = 115 000) was used in high-accuracy mode (HAM).The exposure time for the TOI-2374 observations was 1800 s.This resulted in a signal to noise ratio (SNR) of between 22 and 35 per pixel.The exposure time for the TOI-3071 observations was 2400 s, resulting in a SNR of between 10 and 25 per pixel.The data were reduced using the standard offline HARPS data reduction pipeline, and a K5 template (TOI-2374) or G2 template (TOI-3071) was used to form a weighted cross-correlation function (CCF) to determine the RV values (Baranne et al. 1996;Pepe et al. 2002).The line bisector (BIS) and full-width at half-maximum (FWHM) were measured using methods by Boisse et al. (2011).Three bad data points were removed from the TOI-2374 dataset.One of them (BJD = 2459711.8)was reduced incorrectly and the other two (BJD = 2459712.8and BJD = 2459716.8)contained large errors due to bad weather during the observations.The achieved instrumental precision for the rest of the RV measurements is ≈ 4 m s −1 .The observations are presented in tables A1 and A2.

PFS RVs
We observed TOI-2374 with the Planet Finder Spectrograph (PFS; Crane et al. 2006Crane et al. , 2008Crane et al. , 2010) ) on the Magellan/Clay telescope at Las Campanas Observatory in Chile.We obtained seven observations of this target between 2021 Sep 14 and 2021 Nov 17 through an iodine cell, with relatively short exposure times of 5 minutes.These observations were made in 3x3 binning mode to reduce readout noise.We also obtained an iodine-free template spectrum at higher signal to noise, which we used to extract relative RVs.The RVs were extracted using the pipeline described in Butler et al. (1996), achieving a typical instrumental RV precision of ≈ 3 m s −1 .We also quantified potential variations in the spectral line profiles using bisector inverse slopes (BIS) measurements.These were derived using the procedures described in Hartman et al. (2019), applied to the iodine-free orders of each spectrum.The BIS measurements are shifted such that the median bisector span for each order is zero.All measurements are presented in Table A3.

TRES
We obtained two reconnaissance spectra on UT 2021 July 21 and UT 2021 August 5 of TOI-2374 using the Tillinghast Reflector Echelle Spectrograph (TRES, Fűrész 2008) located at the Fred Lawrence Whipple Observatory (FLWO) in Arizona, USA.TRES has a resolving power of ∼44,000 and operates in the wavelength range 390-910 nm.The spectra were extracted as described in Buchhave et al. (2010).

High resolution imaging
High-angular resolution imaging is needed to search for nearby sources that can contaminate the TESS photometry, resulting in an underestimated planetary radius, or be the source of astrophysical false positives, such as background eclipsing binaries (Daemgen et al. 2009;Lillo-Box et al. 2012, for example).We searched for stellar companions to TOI-2374 with speckle imaging on the 4.1-m Southern Astrophysical Research (SOAR) telescope (Tokovinin 2018) on 3 December 2020 UT, observing in Cousins I-band, a similar visible bandpass as TESS.This observation was sensitive to a 4.7-magnitude fainter star at an angular distance of 1 arcsec from the target.More details of the observations within the SOAR TESS survey are avail-able in Ziegler et al. (2020).The 5 detection sensitivity and speckle auto-correlation functions from the observation are shown in Figure 2. No nearby stars were detected within 3 ′′ of TOI-2374 in the SOAR observations.TOI-3071 was similarly observed by SOAR on 20 March 2022, sensitive to a 5.5-magnitude fainter star at a 1 arcsec separation from the target.The 5 detection sensitivity and speckle auto-correlation functions from this observation are shown in Figure 3.We detected no nearby stars within 3 ′′ of TOI-3071 in the SOAR observations.

Stellar parameters
Because the final physical parameters of the planets depend directly on the values of the stellar parameters, here we perform several independent methods to measure and derive a range of stellar parameters for TOI-3071 and TOI-2374.

Spectroscopic parameters
The stellar spectroscopic parameters ( eff , log , microturbulence, [Fe/H]) were estimated using the ARES+MOOG methodology from the respective combined HARPS spectrum of each star.The methodology is described in detail in Sousa et al. (2021); Sousa (2014); Santos et al. (2013).We used the latest version of ARES7 (Sousa et al. 2007(Sousa et al. , 2015) ) to consistently measure the equivalent widths (EW) on the list of iron lines.For TOI-3071 we used the list of lines presented in Sousa et al. (2008), while for TOI-2374 we used the list of lines presented in Tsantaki et al. (2013) which is more appropriate for stars cooler than 5200 K.In the analysis we use a minimization process to find the ionization and excitation equilibrium to converge for the best set of spectroscopic parameters.This process makes use of a grid of Kurucz model atmospheres (Kurucz 1993) and the radiative transfer code MOOG (Sneden 1973).We also derived a more accurate trigonometric surface gravity using recent GAIA data following the same procedure as described in Sousa et al. (2021).We estimated rotational velocity  sin  values by performing spectral synthesis (with the same code and model atmospheres as for stellar parameters) of two iron lines in the 6705Å region.We fix the macroturbulence velocity with the empirical formula by Doyle et al. (2014), which provides v  =4.6 km s −1 for TOI-3071.Cool stars are outside of such calibration, hence we adopted a value v  =2 km s −1 for TOI-2374.Both stars seem to be slow rotators with  sin  = 2.0 km s −1 for TOI-3071 and  sin  < 0.5 km s −1 for TOI-2374.The stellar parameters are presented in tables 1 and 2.
Independently, we derived the stellar atmospheric parameters for TOI-2374 from the TRES spectra using the Stellar Parameter Classification (SPC, Buchhave et al. 2012) tool.SPC cross correlates an observed spectrum against a grid of synthetic spectra based on Kurucz atmospheric models (Kurucz 1992).With this method we obtained  eff = 4937 ± 50 K, log  = 4.61 ± 0.10 and [Fe/H] = 0.24 ± 0.08, which are in agreement (within errors) with the results derived from the HARPS spectra.

Abundances
Stellar abundances of the elements were derived using the same tools and models as for stellar parameter determination under the assumption of local thermodynamic equilibrium.For the derivation of chemical abundances of refractory elements we closely followed the methods described in e.g.Adibekyan et al. (2012Adibekyan et al. ( , 2015)); Delgado Mena et al. (2014Mena et al. ( , 2017)).Abundances of the volatile elements, C and O, were derived following the method of Delgado Mena et al. (2021) and Bertran de Lis et al. (2015).Nine out of the 14 individual spectra of TOI-3071 were contaminated around the 6300.3Å oxygen line and had to be discarded before being coadded, leading to a lower S/N spectrum and thus a higher error in oxygen abundance.These two elements could not be derived for TOI-2374 because their lines are very weak in cool stars spectra.All the [X/H] ratios are obtained by doing a differential analysis with respect to a high S/N solar (Vesta) spectrum from HARPS.The abundances of the elements are shown in table B1 in the appendix.

Physical parameters
From the parameters obtained in the spectroscopic analyisis we estimated the stellar radius and mass using the calibrations presented in Torres et al. (2010) (see tables 1 and 2).In addition, we performed an analysis of the broadband spectral energy distribution (SED) of the star together with the Gaia DR3 parallax, in order to determine an empirical measurement of the stellar radius (Stassun & Torres 2016;Stassun et al. 2017Stassun et al. , 2018)).We pulled the   magnitudes from 2MASS, the W1-W3 magnitudes from WISE, and the  BP  RP magnitudes from Gaia.We also utilized the absolute flux-calibrated Gaia spectrum where available.Together, the available photometry spans the full stellar SED over the wavelength range 0.4-10 m (see figures 4 and 5).We performed a fit using PHOENIX stellar atmosphere models (Husser, T.-O. et al. 2013), adopting from the spectroscopic analysis the effective temperature ( eff ), metallicity ([Fe/H]), and surface gravity (log ).We fitted for the extinction   ments of the Mount-Wilson S-index.With these values we obtained estimates of the rotational period of the stars  rot .For both stars,  rot is consistent with the maximum rotation period calculated with the  sin  estimates,  rot /sin  > 70 d for TOI-2374 and 33 ± 8 d for TOI-3071.We used Eq. 3 and Eq.16 in Barnes (2007) to calculate the gyrochronological stellar ages from the rotation periods, giving 4.34 ± 0.67 Gyr and 4.9 ± 0.8 Gyr for TOI-2374 and TOI-3071, respectively.Moreover, we obtained two other independent estimates of the stellar ages : first, we used the web interface PARAM 1.38 , which estimates the basic intrinsic parameters of stars given their photometric and spectroscopic data following the method described in da Silva et al. (2006).This resulted in stellar age estimates of 4.8 ± 4.2 Gyr and 0.8 ± 0.6 Gyr for TOI-2374 and TOI-3071, respectively.As expected, the age determination based on evolutionary models is extremely uncertain for main-sequence stars.The age of TOI-2374 formally agrees with the one from gyrochronology.On the other hand, PARAM provides a much younger age for TOI-3071.We then used the chemical abundances of some elements to derive ages through the so-called chemical clocks (i.e.certain chemical abundance ratios which have a strong correlation for age).We applied the 3D formulas described in Table 10 of Delgado Mena et al. (2019), which also consider the variation in age produced by the effective temperature and iron abundance.
and [Sr/Si] were used from which we obtain a weighted average age of 1.1 ± 0.2 Gyr for TOI-3071.This age is again much younger than the age from gyrochornology.The ages from chemical clocks for cool stars such as TOI-2374 must be taken with caution since the calibrations are obtained using hotter stars.Therefore, we choose to apply the 2D formulas considering the abundance ratios and the metallicity (see Table 8 B2.We include in Table 1 the geometric average of the three stellar age estimates for TOI-2374.The error is given by the deviation of the three values.In Table 2 we report two values for the age of TOI-3071, corresponding to the two disagreen results: the gyrochronological age and the weighted average of the ages derived with PARAM and with the chemical clocks.

RV analysis
We carried out a preliminary stage of RV data exploration in two parts: periodogram analysis through the DACE platform9 and testing of linear RV models with different stellar activity indicators as model covariates.These examinations showed no indication of stellar activity in either of the systems.However, the time span is short for TOI-3071.
Despite having a relatively small number of points, we were able to find the planets of both systems independently of photometry.-3071).The RV periodograms present peaks above the 0.01 per cent False Alarm Probability at the expected planetary periods.After removing the best fit model for a Keplerian orbit, the periodograms of the RV residuals show no further significant peaks, therefore ruling out the detection of a second planet in each system from these time series.
When looking at periodograms of TOI-3071 stellar activity indicators, we found that the bisector span, contrast and Ca-index show power at low frequency.However, as the time span is relatively short (13.9 days), we cannot constrain the period.We also found one peak in the   periodogram (and, to a lesser extent, in the Na periodogram) that is similar to the planetary period but is not large enough to be deemed significant.No other stellar activity indicator shows signs of significant frequencies in this system.As to TOI-2374, there is a significant peak around 1 day in the Na periodogram, which corresponds to a 1 d alias of a low frequency signal, which disappears when we substract a long-term drift from the timeseries.Other than that, this star does not show signs of any other activity.
We computed the Pearson's R statistic for the RV data and the stellar activity indicators and found no significant correlation between those variables for any system.As for the correlations between the stellar activity indicators and the RV residuals after removing the best  fit for a Keplerian orbit, the highest value obtained was R∼ 0.61 between TOI-2374 RV residuals and its FWHM, but it strongly depends on a single point (BJD = 2459750.72),without which the correlation disappeared and the Keplerian model stays essentially unchanged.The same situation is seen between the TOI-3071 RV residuals and its contrast, which shared a R of ∼ 0.61, but strongly depends on four points.No other indicator shows significant correlation with the RV residuals for none of the systems.
To further examine the effects of stellar activity on the RV data, we tested if the RV curves of our joint model described in section 3.3 should include linear terms corresponding to any of the stellar activity indicators.For this we built simple linear models where the dependent variable was the radial velocity and the covariates were each of the available stellar activity indicators, placed individually in each linear model.We performed an F-test for each model compared to a constant RV model and none of them obtained a p-value smaller than 0.01.We then incorporated into each model a circular orbit with a period fixed to the planetary period found and we repeated the F-test comparing them to the model that only included the circular orbit.Again, none of the stellar activity indicators obtained a p-value smaller than 0.01.This shows that it is not necessary to incorporate linear terms with these covariates to explain the RV data.

Joint modelling
For each of the TOI-2374 and the TOI-3071 planetary systems, we used the exoplanet package (Foreman-Mackey et al. 2020) to model the RV data and light curves simultaneously.exoplanet utilises the probabilistic programming package PyMC3 (Salvatier et al. 2016) and the light curve modelling package Starry (Luger et al. 2019).The radial velocity mean for each spectroscopic instrument was subtracted from the data.QLP and ground-based light curves were already normalized from the pipelines to have an out-of-transit flux of one.The TESS-SPOC light curves were normalized by dividing them by the median of the flux.All timestamps were converted to the time system used by TESS, i.e.BJD -2457000 (BJD-TDB).
To model the planetary transits, we used a limb-darkened transit model utilising the quadratic limb-darkening parameterisation in Kipping (2013) and a Keplerian orbit model, which was parameterised for the planet in terms of the orbital period , the epoch  0 , the eccentricity , the argument of periastron , the impact parameter , the stellar radius  ★ and the stellar mass  ★ .These parameters were then input into light curve models created with Starry, alongside parameter   (planetary radius), hyperparameter  exp (the exposure time of the instrument) and the time series of the data .The RV model was computed from the Keplerian orbit model by adding the semi-amplitude  of the RV signal as a parameter.Based on the analysis described in the previous section, we modeled only one planet per system and found no need to include a Gaussian process as a model for the noise given that periodogram analysis of RV and LC data showed no indication of stellar activity or signals from stellar rotation at sufficiently significant frequencies.Instead, we assumed the errors are normal and included an additional white noise term for every instrument, which was added as a free nonnegative jitter term in quadrature with the reported errors.
All prior distributions set on the parameters fit in this model are given in tables C1 and C2.We put Gaussian priors on the stellar radius and mass informed by the values reported in section 3.1.3.We used values given by the QLP detection for the epoch and period to inform our Gaussian priors for log- and  0 in both planets, and depth values from the same pipeline to inform the Gaussian prior on log-, where  is the planet-to-star radius ratio  =    ★ , with a mean of 0.5 log-depth.We inflated the widths of these priors to  = 1.Because it is well known that sampling directly the eccentricity and the argument of periastron can be problematic for most MCMC samplers (Parmentier & Crossfield 2018), we sampled for

√
sin  and

√
cos  with a uniform prior within a unit disk.This leads to a uniform prior on  as noted by Anderson et al. (2011).For the impact parameter we chose a uniform prior between 0 and 1 + .For each lightcurve, we put a Gaussian prior on the transit normalisation factor  0 (the light curve flux level out of transit) with a mean of 1 and standard deviation of 0.1.For the jitter parameter of each instrument we used a wide Gaussian prior, the mean of which was the log of the error median on each light curve.We used independent limb-darkening parameters for each filter, parameterised following Kipping (2013).For each spectroscopic instrument, we introduced a constant offset term to the reported RVs.We also sampled for a  linear and quadratic trend in the RV data.We used uniform priors for log- and for the log-jitter of each RV instrument.
We added one dilution factor for each TESS sector as in Almenara et al. (2022). 10Here, we use this parameter to complement the point of the PyMC3 sampler, which draws samples from the posterior using a variant of Hamiltonian Monte Carlo, the No-U-Turn Sampler (NUTS).We allowed for 10000 burn-in samples which were discarded, and then 30000 steps with 10 chains.The MCMC chains for all model parameteres have Gelman-Rubin statistics (rank normalized split-R) (Vehtari et al. 2021) very close to unity (< 1.01), effective sample sizes (bulk-ESS) larger than 9100 and tail-ESS larger than 4600, so we have a good level of confidence that the chains mixed well and we have a stable estimate of uncertainty.The median model and the interval between the 16th and 84th percentiles are shown in figures 8-12.We present the parameters' marginal posterior sample median for the TOI-2374 and TOI-3071 systems from the joint fits in table 3.

RESULTS & DISCUSSION
The combined analysis of high resolution spectroscopy and spaceand ground-based photometry allowed the determination of the masses and radii of planetary companions for the two target systems.We determined that both planets are sub-Saturns with mass and radius precision of 4−7% (see Table 3).From this we inferred bulk densities of 0.98 +0.15  −0.13 g cm −3 for TOI-2374 b and 1.02 +0.26 −0.22 g cm −3 for TOI-3071 b.These parameters place both planets within in the Neptunian desert, with TOI-3071 b being much deeper in than TOI-2374 b (Fig. 13).We found the eccentricity of both targets to be consistent with 0. The 95% highest density interval (HDI) for the eccentricity ranges from 0 to 0.13 and 0.09 for TOI-2374 b and TOI-3071 b, respectively.We did not find any statistically significant long-term trends in the RV measurements of either star.
We note that the MCMC posterior distributions for the TESS light curves dilution correction factors are narrower than the priors (see figures C1 and C2 and tables C1 and C2 in the Appendix).For TOI-2374, the posterior modes of the QLP dilution correction factors are systematically negative, hinting to an overestimation of the dilution factor estimated by this pipeline.On the other hand, in TOI-3071, the posterior mean of the QLP dilution correction factor is positive, suggesting that for this star, the correction procedure underestimates the dilution.Finally, our model also hints to an overestimation of the dilution factor in the SPOC lightcurve of TOI-3071, as the posterior mean of  is negative.Exploring the reasons of these behaviours is outside the scope of this article.Additionaly, in all cases the 83% HDI contains the null correction factor value.In any case, our model accounts for the uncertainty in the dilution factors estimated by the pipelines on the determination of the planetary radii.

Comparison to other exoplanets
In order to understand how TOI-2374 b and TOI-3071 b fit into the landscape of known planets, we compare them to planets of similar size, mass, and period ( < 20 ; 40 <   < 100  ⊕ ;  < 2   −3 ) with mass precision better than 30% and radius precision better than 20% (top panel of Fig. 14).This includes all "hot" giant planets below Saturn's mass.We use the parameters from the planetary systems table on the NASA Exoplanet Archive11 accessed on 28 May 2024.The empirical R-M relation for volatile-rich planets ( < 3.3   −3 ) found by Otegi et al. (2020) is also shown in the diagram.The bottom panel of Fig. 14 shows the equilibrium temperature ( eq ) vs radius plane for the same subset of exoplanets, but extending the mass range to 10 <   < 300  ⊕ . eq is calculated from the planetary orbital semi-major axes  and the host-star effective temperature  ★ and radius  ★ assuming full day-night heat distribution, according to where  is the Bond albedo of the planet.For this analysis we assumed  = 0.5, similar to that of Jupiter (Li et al. 2018), to calculate  eq for TOI-2374 b, TOI-3071 b and every other planet in the dataset (even if they already had a reported  eq ).Both planets are hot and highly irradiated, with  eq ≈ 744 K for TOI-2374 b and  eq ≈ 1812 K for TOI-3071 b.TOI-2374 b is slightly below the usual inflated hot-Jupiter cutoff of  * ≈ 2 × 10 8 The error bars represent the reported uncertainty and the empirically derived per-instrument jitter, added in quadrature.These RV measurements are also listed in Table A1.Top: RV measurements plotted against time.Bottom: phase-folded RV measurements.erg/s/cm 2 (Miller & Fortney 2011), and it is therefore not expected to be anomalously inflated.On the other hand, TOI-3071 b is well inside the hot-Jupiter instellation regime.However, previous studies suggest that planets with masses below about 0.4  show a weaker response to high instellation fluxes compared to more massive planets, and are less likely to be strongly inflated (Thorngren & Fortney 2018;Sestovic et al. 2018).Therefore, standard evolution models of giant planets may reflect the radius more accurately than for higher mass planets, for which the inflation is known to be underestimated.From Fig. 14 it is clear that TOI-3071 b has a similar radius to previously detected giant exoplanets with similar masses, despite having a much higher equilibrium temperature.Besides hydrogen and helium, the planet is therefore expected to contain a large amount of heavier elements.

Inferred heavy-element masses
Both planets have available masses, radii and ages, which make them viable targets for the estimation of their heavy-element masses (e.g., Miller & Fortney 2011;Thorngren et al. 2016;Müller & Helled 2023b).The bulk heavy-element mass of a giant planet is a key property since it can be used to test planet formation models and provide additional constraints on the possible formation pathways (e.g., Guillot et al. 2006; Johansen & Lambrechts 2017; Hasegawa et al. 2018).Unfortunately, it cannot be measured directly and must be inferred from a combination of measurements and giant-planet evolution models.To estimate the composition of the two planets, we modeled their evolution with a modified version of the stellar evolution code Modules for Experiments in Stellar Astrophysics (MESA; Paxton et al. (2011Paxton et al. ( , 2013Paxton et al. ( , 2015Paxton et al. ( , 2018Paxton et al. ( , 2019)); Jermyn et al. (2023)).Müller et al. (2020b,a) modified the equation of state of MESA to make it suitable to model giant planets that contain a significant amount of heavy elements.Here, we used a recently updated version (Müller & Helled 2024), which uses a new hydrogen-helium equation of state that includes non-ideal interactions in the hydrogen-helium mixture (Chabrier & Debras 2021).These effects have been shown to lead to a significant change in the cooling of planets (Howard & Guillot 2023).We assumed a core-envelope structure and a 50-50 water-rock composition for the heavy elements.The radiative opacity was from Freedman et al. (2014), and to account for the instellation fluxes we used a gray atmospheric model developed for irradiated planets (Guillot 2010;Parmentier & Guillot 2014).
We first calculated the planetary evolution assuming that they have the same composition as their stars.Both planets are much smaller compared to these predictions, demonstrating that they must have super-stellar metallicities (see solid lines in Fig. D1 in Appendix).
To quantify the metal enrichment of the planets, we used the standard Monte-Carlo approach (see, e.g., Müller & Helled 2023a, for a review): sample planets were created by drawing from the observation prior distributions for the mass, radius, and stellar age.For the mass and radius, we used normal distributions informed by the parameters' marginal posterior sample median and standard deviation from the joint fit (Sect.3.3), and for the age we used a uniform distribution based on the estimates from the stellar parameters analysis (Sect.3.1.3).The evolution models were then used to calculate the cooling for different heavy-element masses   to find the one that agrees with the observations.This process is repeated until there is enough data to estimate the posterior distribution for the heavy-element mass.For TOI-3071b we used two different stellar age priors (from chemical clocks and gyrochronology) and then combined their posteriors for the final result assuming that the priors are equally likely.The MCMC approach required the calculation of a large number of evolution models.Since this would have been too slow, we followed a similar approach as outlined in Müller & Helled (2021).For each planet, we calculated a grid of evolution models (in the planet mass-metallicity space) and then used 2d piecewise monotonic cubic interpolation to generate new cooling tracks.
The resulting posterior distributions of the inferred heavy-element masses are shown in Fig. 15.We find heavy-element masses (mass fractions) of   = 33.3± 3.8 ⊕ ( = 0.59 ± 0.05) for TOI-2374 b and   = 45.1 ± 5.5 ⊕ ( = 0.66 ± 0.07) for TOI-3071 b.For comparison, Saturn, the most similar planet in the solar system, has an estimated bulk heavy-element mass (mass fraction) of   ≃ 30 ⊕ ( ≃ 0.3).The exact value is unknown and depends on several assumptions when modeling Saturn's interior.Both targets are metalrich, with TOI-3071b being more enriched both in absolute and relative terms (see Figure D2 for posterior distributions of the heavyelement masses in mass fractions).TOI-3071b in particular has an extremely high bulk metallicity and is about 60-70% heavy elements (water and rocks).
The fact that these planets have such high metallicities challenges standard formation models which generally tend to predict lower heavy-element masses (Pollack et al. 1996;Helled et al. 2014;Johansen & Lambrechts 2017;Bitsch et al. 2018).Interestingly, both planets correspond to masses below Saturn's mass, which was recently suggested to be the transition mass to become giant planets (Helled 2023).In that view, the formation of a giant planet can last for a few million years where gas accretion is being delayed due to continuous accretion of planetesimals (Alibert et al. 2018;Helled et al. 2022b).Such a formation scenario also suggests that planets below about Saturn's mass are expected to be metal-rich in composition, consistent with our estimates for TOI-2374 b and TOI-3071 b.An alternative explanation is a post-formation enrichment due to planetesimal accretion or collisions between planetary embryos (Mousis et al. 2009;Shibata et al. 2020;Ginzburg & Chiang 2020).However, such collisions have difficulty enriching planets with tens of Earth masses of heavy elements.
Based on the models from Santos et al. (2015Santos et al. ( , 2017)); Adibekyan et al. (2021), we used the stellar abundances to calculate the summed mass fraction of all heavy elements available as planetary building blocks.The values are   −2374 = 0.018±0.03and   −3071 = 0.024 ± 0.03, respectively.For comparison, the same models yield  ⊙ = 0.013 ± 0.01 for solar abundances.Both mass fractions are above the solar value, which hints at a possible connection between the stellar abundances and the high inferred bulk metallicities.
Finally, we note that measurements of the atmospheric chemical composition could be used to further constrain the planetary for-

Internal structure
To model the evaporation histories of the two planets, we assumed a simpler two-layer internal structure with a rocky core surrounded by a pure H/He envelope.This model can be described with four quantities: the rocky core mass M core and core radius R core , and the gaseous envelope thickness R env and envelope mass fraction f env , defined as the ratio of envelope mass to planet mass M env /M p .In order to determine these parameters, we adopted the empirical massradius relation by Otegi et al. (2020) for the rocky core, and the envelope structure model by Chen & Rogers (2016) for the gaseous atmosphere, which is based on MESA simulations.Uisng this model, we thus determined envelope mass fractions of 38 ± 5 % for TOI-2374 b and 26 ± 6 % for TOI-3071 b.These correspond to core (heavy-element) masses of M core = 34.7 ± 3.7 for TOI-2374 b, and M core = 50.5 ± 4.9 for TOI-3071 b, both consistent with the heavyelement masses determined in Sect.4.2 within 1.

Evolution modelling
Simulating the evolution of a planet's atmosphere under photoevaporation requires knowledge of its X-ray irradiation history.We can estimate the X-ray and EUV emission (together, XUV) of a star from its spin period, as the two are related via rotation-activity relations (e.g.Wright et al. 2011;Pizzolato et al. 2003).We adopted the stellar rotation evolution models by Johnstone et al. (2021), which simulate the evolution of the angular momentum of stars constrained by observations of young star clusters.Moreover, to estimate the EUV component, we adopted the scaling relation by King & Wheatley (2021), which are based on Solar TIMED/SEE data.We plot the evolution of the spin period and corresponding XUV emission of TOI-2374 and TOI-3071 on Figures 16 and 17, respectively.
Regarding TOI-2374, we found the rotational models by Johnstone et al. (2021) are consistent with the measured age and spin period of the star.
For TOI-3071, we considered the two age estimates separately: an age of 1 Gyr from the chemical clocks and PARAM isochrones, and an age of 5 Gyr from gyrochronology (see Sect. 3.1.3).
We found the choice of age had a strong impact on the X-ray emission history of the star.The younger age is more consistent with a short spin period of 2-4 d, whereas the older age from gyrochronology is more consistent with a slower spin period of 20-30 d, as the star has had more time to spin down.The spin period of the star, however, is constrained to 26 d, and so adopting the younger age would imply that the star must be unusually slowly rotating and X-ray faint (Fig. 17, red line).
For that reason, we adopted two scenarios for the X-ray emission history of TOI-3071: (1) the high activity scenario where we adopted the predicted evolution from the models by Johnstone et al. (2021), consistent with the older gyrochronological age (see Fig. 17, black line); and (2) the low activity scenario, where we adopted a fainter emission history model consistent with the younger age from the 3D chemical clocks (Fig. 17, red line).
We then simulated the past and future evaporation histories of the planets with the photoevolver code (Fernández Fernández et al. 2023).The simulation is built upon three ingredients: (1) the XUV emission history of the star, which determines the X-ray irradiation on the planet, (2) a mass loss formulation, which determines the amount of gas lost due to XUV irradiation, and (3) an envelope structure model, which relates the mass of the gaseous atmosphere to its thickness.Therefore, on each time step, the X-ray luminosity of the star at that age is used to calculate the mass loss rate, which is then used to recalculate the mass and size of the planet.We adopted the models by Johnstone et al. (2021) for the stellar XUV emission history (described above), and the model by Chen & Rogers (2016) for the structure of the gaseous atmosphere, described in Section 4.3.1.Regarding the mass loss formulation, we adopted the model by Kubyshkina et al. (2018).We evolved the planets' atmospheres using the 4th order Runge-Kutta algorithm (Dormand & Prince 1980) as the integration method with a variable time step no larger than 1 Myr.We ran the simulation backwards to the age of 10 Myr and forwards to 10 Gyr.In the case of the F-type TOI-3071, we stop the simulation at 6 Gyr, the approximate age at which the star leaves the main sequence (Choi et al. 2016).
Moreover, for each of the two planets, we also explored a range of possible evaporation histories taking into account the uncertainties on the planet's mass, radius, and XUV irradiation.To maximise the evaporation rate, we picked a radius on the higher end of its 1 uncertainty and a mass on the lower end of its 1 uncertainty, which minimises the density and maximises the escape rates.We also picked the highest XUV emission history allowed by the model's 2 spread, as shown on Figures 16 and 17.Likewise, to minimise the evaporation rate, we picked the opposite parameters: values for the mass and radius within their 1 errors that maximise density under the faintest XUV irradiation history allowed by the model's 2 spread.

Simulation results
The evaporation histories of TOI-2374 b and TOI-3071 b are plotted on Figure 18 (left and right panels, respectively).
We found that TOI-2374 b is stable against evaporation, as its envelope mass fraction changes only by a few percent across its lifetime.The evolution of its radius is thus largely driven by thermal contraction.We also constrained its current mass loss rates to 0.4 − 0.8 × 10 8 g s −1 .
On the other hand, TOI-3071 b is more susceptible to escape despite its higher mass.This is due to being twice as close to its star compared to TOI-2374 b, as well as the higher X-ray luminosity of its F-type host star.Under the higher X-ray irradiation history consistent with the gyrochronological age (5 Gyr), our simulations show that the planet is reaching the point of losing its atmosphere completely when its star starts evolving off the main sequence (see Fig. 18), with an envelope mass fraction of 10-20% at 6 Gyr down from 40-60% initially.Under this X-ray irradiation history, the planet currently experiences mass loss rates of 7.5 − 40.7 × 10 8 g s −1 .
However, under the fainter irradiation history consistent with the much younger age derived from chemical clocks (1 Gyr), TOI-3071 b experiences lower escape rates throughout its life and is stable against evaporation.We constrained its current mass loss rates to 0.9 − 2.6 × 10 8 g s −1 , one order of magnitude lower compared to the high activity scenario.
In Sect.4.2, we presented more detailed internal structure models to estimate the bulk metallicity of the two planets.These models cannot predict the atmospheric composition, which can only be determined by observations.However, since both are metal-rich, they are also likely to have high-metallicity atmospheres.Such highmetallicity atmospheres experience much lower escape rates as heavier species require more energy to remove (Wilson et al. 2022;Owen & Jackson 2012).Since we only considered a pure H/He composition in our atmospheric evolution models, our results in Fig. 18 constitute an upper bound to the evaporation histories of the planets.
For TOI-2374,b, the atmosphere is stable against evaporation throughout its life, even when assuming a pure hydrogen envelope.TOI-3071,b experiences moderate escape rates and is significantly affected by evaporation only as its host star transitions out of the main sequence, risking total atmospheric loss.However, heavy-mass planets like TOI-3071 b could present a metallicity gradient in their gaseous envelopes, with a relatively pure H/He upper atmosphere that progresses towards metal-rich layers deeper in the atmosphere.If the planet experienced stronger escape rates earlier in its life, either from photoevaporation or Roche-lobe overflow, the pure H/He upper layers of the envelope would be stripped off first, exposing the deeper metal-rich layers which would then curtail atmospheric escape (Owen & Murray-Clay 2018;Hoyer et al. 2023).
Therefore, the presence of heavy species in the hydrogen-rich envelope of TOI-3071 b at the present time would lead to diminished escape rates, making its atmosphere more resistant to evaporation and thus likely stable.As a result, we find the two planets are not significantly affected by photoevaporation and are thus able to hold onto their gaseous envelopes throughout their lifetimes.

Prospects for atmospheric characterization
The transmission spectroscopy metric (TSM) is defined as (2) where   is the apparent magnitude of the host star in the J band and the scale factor is 1.15 for 4 <   < 10  ⊕ (Kempton et al. 2018).
We calculated the targets' TSM using the equilibrium temperatures for albedos of  = 0.5 and obtained a value of 84 ± 13 for TOI-2374 b and 36 ± 9 for TOI-3071 b.The former is slightly below the threshold of 90 set by Kempton et al. (2018) for sub-jovians to be considered high quality atmospheric characterization targets.However, if we assume lower albedos of 0.3 and 0.0, we obtain TSM values of 91 ± 15 and 100 ± 16 for TOI-2374 b, respectively.For TOI-3071 b the TSM is always lower than 50.As mentioned above, atmospheric measurements of these targets could provide key information to further constrain the planetary formation, evolution and interior structure.

CONCLUSION
We have presented the discovery and characterisation of two new transiting hot Saturns orbiting the K-type and F-type stars TOI-2374 and TOI-3071, respectively.Our analysis includes photometry from TESS, follow up ground-based photometry and spectroscopy from HARPS and PFS.The combination of this data in a joint light curve and radial velocity model allowed the confirmation of the planets and the determination of their orbital and physical properties.
Both planets are hot sub-Saturns in close, virtually circular orbits.They have similar masses and radii and are located within the Neptunian desert.Due to its closer orbit, TOI-3071 b is much deeper in the desert than TOI-2374 b.This, and the fact that TOI-3071 is a much hotter star and emits higher XUV flux, leads to TOI-3071 b having a much higher equilibrium temperature, higher even than that of all previously detected giant exoplanets with similar masses or radii.Both planets are very metal rich, with TOI-3071 b being more enriched both in absolute and relative terms.By studying the evolution of the planets' atmospheres under photoevaporation we find that both are stable against evaporation.
Both planets are amenable to further follow-up with JWST, as despite the relatively low TSM of TOI-3071 b, it is a fairly bright (11.87 T mag) target undergoing atmospheric evaporation.Meanwhile, TOI-2374 b has a higher TSM and lies close to the Continuous Viewing Zone of both TESS and JWST for ease of observation scheduling.Measuring the atmospheric chemical compositions would prove very useful to further constrain the formation, evolution and interior of the planets.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).
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.
Funding for the TESS mission is provided by NASA's Science Mission Directorate.KAC and CNW acknowledge support from the TESS mission via subaward s3449 from MIT.
This research was funded in part by the UKRI, (Grants ST/X001121/1, EP/X027562/1).For the purpose of open access, the author has applied a Creative Commons Attribution (CC BY) licence to any Author Accepted Manuscript version arising from this submission.

DATA AVAILABILITY
The TESS data are available from the Mikulski Archive for Space Telescopes (MAST), at https://heasarc.gsfc.nasa.gov/docs/tess/data-access.html.The other photometry from LCOGT, Briefrield Private Observatory and PEST, and the high-resolution imaging data, are available for public download from the ExoFOP-TESS archive at https://exofop.ipac.caltech.edu/tess/target.php?id=439366538 for TOI-2374 and at https://exofop.ipac.caltech.edu/tess/target.php?id=452006073 for TOI-3071.The PFS RV data is also available here.The HARPS RV data is shown in tables A1 and A2 and the full HARPS data products can be found on the ESO archive.The model code underlying this article will be shared on reasonable request to the corresponding author.The shaded regions are cooling curves using the best-fit inferred compositions within 1.The error bars show the observational data.Both planets are much smaller than the predictions from the evolution models, suggesting super-stellar bulk metallicities. 13

Figure 1 .
Figure 1.Target Pixel Files (TPF) for TOI-2374 from TESS S28 (top) and TOI-3071 from S64 (bottom).Targets are marked as a white cross.Other Gaia DR3 sources within a limit of 4 Gaia magnitudes difference are marked as red circles, and are numbered in distance order from the targets.The aperture mask is outlined and shaded in red.This figure was created with tpfplotter (Aller et al. 2020).
We observed two full transit windows of TOI-2374 b on UTC 2021 May 10 and 2021 May 18 in Sloan  ′ and Sloan  ′ bands, respectively, using the Las Cumbres Observatory Global Telescope (LCOGT; Brown et al. 2013) 1.0 m network nodes at South Africa Astronomical Observatory (SAAO) and Siding Spring Observatory (SSO).The images were calibrated by the standard LCOGT BANZAI pipeline (McCully et al. 2018) and differential photometric data were extracted using AstroImageJ (Collins et al. 2017).We used circular photometric apertures with radius 5.1 ′′ for the  ′ band observations and 5.4 ′′ for the  ′ band observations.The target star apertures excluded all of the flux from the nearest known neighbor in the Gaia DR3 catalog (Gaia DR3 6828814283414902784), which is ∼ 22 ′′ northeast of TOI-2374.Each light curve is included in the global modelling described in section 3.3.

Figure 2 .Figure 3 .
Figure 2. SOAR observation of TOI-2374.Black circles correspond to measured data points and the black lines show the fit in two different separation regimes.The 5 detection sensitivity is plotted with the speckle imaging auto-correlation functions inset.
, limited to the maximum line-of-sight value from the Galactic dust maps ofSchlegel et al. (1998).The resulting fits (figures 4 and 5) have   = 0.02 ± 0.02 and 0.20 ± 0.03 for TOI-2374 and TOI-3071, respectively, with a reduced  2 of 1.3 and 1.4, respectively.Integrating the (unreddened) model SED gives the bolometric flux at Earth,  bol = 4.981 ± 0.058 × 10 −10 erg s −1 cm −2 and 3.291 ± 0.038 × 10 −10 erg s −1 cm −2 , respectively.Taking the  bol together with the Gaia parallax directly gives the bolometric luminosity,  bol = 0.2867 ± 0.0034 L ⊙ and 2.539 ± 0.032 L ⊙ , respectively.The Stefan-Boltzmann relation then gives the stellar radius,  ★ = 0.774 ± 0.025 R ⊙ and 1.393 ± 0.030 R ⊙ , respectively.These estimates differ by 2.6 and 1.7 with the ones obtained from the spectroscopic analysis, respectively.Given that the two sets disagree by less than 3 and that the former analysis works with the full spectrum, we decided to use the ones estimated from the calibrations inTorres et al. (2010) as input for the joint model described in section 3.3.We used the relations stated inNoyes et al. (1984) to derive chromospheric activity index log  ′ HK values from the HARPS measure-

Figure 4 .
Figure 4. Spectral energy distribution of TOI-2374.Red symbols represent the observed photometric measurements, where the horizontal bars represent the effective width of the passband.Blue symbols are the model fluxes from the best-fit PHOENIX atmosphere model (black).

Figure 5 .
Figure 5. Spectral energy distribution of TOI-3071.Red symbols represent the observed photometric measurements, where the horizontal bars represent the effective width of the passband.Blue symbols are the model fluxes from the best-fit PHOENIX atmosphere model (black).The absolute flux-calibrated Gaia spectrum is shown as a grey swathe in the inset figure.

Figure 6 .
Figure 6.GLS Periodograms (Zechmeister & Kürster 2009) for the TOI-2374 HARPS data (top) and HARPS + PFS data (bottom).The highest peak in the RV periodogram, corresponding to the orbital period of TOI-2374 b is denoted by a dashed vertical red line.The 0.1, 0.05, and 0.01 per cent False Alarm Probabilities (FAP) are calculated using the approach of Baluev (2008) and are shown as dashed horizontal lines.In the HARPS + PFS figures, the second periodogram corresponds to the the radial velocity residuals after the best fit model for a Keplerian orbit has been removed.

Figure 7 .
Figure 7. GLS Periodograms(Zechmeister & Kürster 2009) for the TOI-3071 HARPS data.The highest peak in the RV periodogram, corresponding to the orbital period of TOI-3071 b is denoted by a dashed vertical red line.The 0.1, 0.05, and 0.01 per cent False Alarm Probabilities (FAP) are calculated using the approach ofBaluev (2008) and are shown as dashed horizontal lines.The second periodogram corresponds to the the radial velocity residuals after the best fit model for a Keplerian orbit has been removed.

Figure 8 .
Figure 8.Light curves of TOI-2374.TESS light curve for Sectors 1 and 28 (grey dots), with time given as Barycentric Julian Date (BJD).The top plots show the TESS QLP KSPSAP light curves from sector 1 (30 minute cadence) and sector 28 (10 min cadence).Overplotted in green is the best-fit solution to the global model resulting from the analysis described in section 3.3.The middle plots show this same data phase folded.The binned points indicate the mean of each bin, with error bars representing the standard error of the mean with the binned data shown in black.Residuals of this fit can be found in the bottom plots.

Figure 9 .
Figure 9.Light curves of TOI-2374 from ground-based telescopes.From left to right, the light curves correspond to the Brierfield Observatory, the LCO-SAAO, the LCO-SSO with ip filter, and the LCO-SSO with gp filter.The data are shown as grey dots, with binned values in black.Overplotted in red is the best fit model from the joint fit described in section 3.3.The bottom plots show the residuals of this fit for each lightcurve.

Figure 10 .
Figure 10.Light curves of TOI-3071.TESS light curve for Sectors 10, 37 and 64 (grey dots), with time given as Barycentric Julian Date (BJD).The top plots show the TESS QLP KSPSAP light curves from sector 10 (30 minute cadence) and sector 37 (10 minute cadence) and SPOC PDCSAP light curves from sector 64 (2 minute cadence).Overplotted in blue is the best-fit solution to the global model resulting from the analysis described in section 3.3.The middle plots show this same data phase folded.The binned points indicate the mean of each bin, with error bars representing the standard error of the mean with the binned data shown in black.Residuals of this fit can be found in the bottom plots.

Figure 11 .
Figure 11.Precise RV measurements of TOI-2374.HARPS measurements are shown as red circles and PFS measurements as blue diamonds.The model plotted is the MCMC median of the global model, corresponding to the orbital period of TOI-2374b.Empirically derived linear and quadratic trends and perinstrument offset  have been subtracted from the raw RV measurements.The error bars represent the reported uncertainty and the empirically derived perinstrument jitter, added in quadrature.These RV measurements are also listed in tables A2 (HARPS) and A3 (PFS).Top: RV measurements plotted against time.Bottom: phase-folded RV measurements.

Figure 12 .
Figure 12.Precise RV measurements of TOI-3071.The model plotted is the MCMC median of the global model, corresponding to the orbital period of TOI-3071b.Empirically derived linear and quadratic trends and perinstrument offset  have been subtracted from the raw RV measurements.The error bars represent the reported uncertainty and the empirically derived per-instrument jitter, added in quadrature.These RV measurements are also listed in TableA1.Top: RV measurements plotted against time.Bottom: phase-folded RV measurements.

Figure 13 .
Figure 13.TOI-2374 b (black marker) and TOI-3071 b (red marker) in the context of the Neptunian Desert regions in Mass vs Period (top panel) and Radius vs Period (bottom panel).The dashed blue lines show the delineation of the Neptunian desert from Mazeh et al. (2016), whereas the horizontal dotted blue lines show the updated lower limits for periods ≤ 2 days from Deeg et al. (2023).Known planets were sourced from the NASA exoplanet archive (https://exoplanetarchive.ipac.caltech.edu/)on 28 May 2024.

Figure 14 .
Figure 14.Top: mass-radius diagram of the confirmed planet population with  < 20 , masses 40 <   < 100  ⊕ and bulk densities  < 2   −3 , with radius determined to better than 20% precision and masses determined to better than 30% precision.The dashed lines show constant density values.For comparison, the position of Saturn in the plot is represented with a white diamond marker.The bulk densities for each planet are color coded.The R-M empirical relation found by Otegi et al. (2020) is shown by a red dashed line with the 1 regions colored.Bottom: equilibrium temperatures and radii of the same set of exoplanets, extending the mass range to 10 <   < 300  ⊕ .The bulk densities for each planet are color coded and the markers are sized according to their mass.

Figure 15 .
Figure15.Inferred heavy-element mass posteriors for the two planets.For TOI-3071b, the two curves show the results for the two different stellar age estimates (chemical clocks and gyrochronology), and the filled area is their combined posterior.Both planets are heavy-element rich, with a bulk heavyelement fraction of  ≳ 0.6.

Figure C1 .
Figure C1.Posterior distributions for the dilution correction factors of the TOI-2374 TESS light curves contrasted with the prior distributions.

Figure C2 .Figure D1 .
Figure C2.Posterior distributions for the dilution correction factors of the TOI-3071 light curves contrasted with the prior distributions.

Figure D2 .
Figure D2.Inferred heavy-element mass fraction posteriors for the two planets.For TOI-3071b, the two curves show the results for the two different stellar age estimates (chemical clocks and gyrochronology), and the filled area is their combined posterior.

Table 1 .
An overview of stellar properties for TOI-2374.

Table 2 .
An overview of stellar properties for TOI-3071.
in Delgado Mena et al. 2019) which provide a weighted average age estimate of 3.3 ± 0.7 Gyr.Ages from each individual clock are shown in Table 04434/2020); 2022.06962.PTDC (http://doi.org/10.54499/2022.06962.PTDC).KAC and CNW acknowledge support from the TESS mission via subaward s3449 from MIT. FB carried out this work within the framework of the National Centre of Competence in Research PlanetS supported by the Swiss National Science Foundation.-funded by the European Union (ERC, FIERCE, 101052347).Views and opinions expressed are however those of the author(s) only and do not necessarily reflect those of the European Union or the European Research Council.Neither the European Union nor the granting authority can be held responsible for them.This work was supported by FCT -Fundação para a Ciência e a Tecnologia through national funds and by FEDER through COMPETE2020 -Programa Operacional Competitividade e Internacionalização by these grants: UIDB/04434/2020; UIDP/04434/2020.This work makes use of observations from the LCOGT network.Part of the LCOGT telescope time was granted by NOIRLab through the Mid-Scale Innovations Program (MSIP).MSIP is funded by NSF.This research has made use of the Exoplanet Follow-up Observation Program (ExoFOP; DOI: 10.26134/ExoFOP5) website, which is operated by the California Institute of Technology, under contract with the National Aeronautics and Space Administration under the Exoplanet Exploration Program.This research has made use of the NASA Exoplanet Archive, which is operated by the California Institute of Technology, under contract with the National Aeronautics and Space Administration under the Exoplanet Exploration Program.

Table A1 .
The full HARPS data products can be found on the ESO archive