TOI-1695 b: A Water World Orbiting an Early-M Dwarf in the Planet Radius Valley

TESS is an ongoing NASA mission surveying the entire sky to search for nearby transiting planets (Ricker et al. 2015). The spacecraft orbits the Earth in an elliptical, 2:1 linear synchronous orbit with a period of 13.7 days. Annual observation cycles are split into sectors lasting two orbits, which is about 27 days. The detector consists of four contiguous CCD cameras, each covering a 24° × 24° field of view, making up a 24° by 96° strip aligned along ecliptic latitude lines. TESS uses the MIT Lincoln Laboratory CCID-80 detector with a depletion depth of 100 μm, allowing for sensitivity just redward of 1000 nm. At the blue end, the spectral response is limited by a longpass filter with a cut on wavelength of 600 nm (Sullivan et al. 2015). In years 1 and 2, the data were stored with a short cadence of 2 minutes and a long cadence of 30 minutes. In year 3, an additional ultrashort cadence of 20 s is available, and the long cadence was reduced to 10 minutes.

TOI-1695 was observed in five nonconsecutive TESS sectors between 2019 November 2 and 2022 June 13 UT at two-minute cadence. The full baseline of the TESS observations is 952 days. TESS observations occurred on CCD 2 on camera 3 in Sector 18 (2019 November 2–November 27 UT), on CCD 1 on camera 3 in Sector 19 (2019 November 27–December 24 UT), on CCD 4 on camera 4 in Sector 24 (2020 April 16–May 13 UT), on CCD 3 on camera 4 in Sector 25 (2020 May 13–June 8 UT), and on CCD 3 on camera 4 in Sector 52 (2022 May 18–June 13 UT).

The TESS images were processed by the NASA Ames Science Processing Operations Center (SPOC; Jenkins et al. 2016), which produced Simple Aperture Photometry (SAP; Twicken et al. 2010; Morris et al. 2020) and Presearch Data Conditioning Simple Aperture Photometry (PDCSAP) light curves per sector. The latter were corrected for systematic uncertainties exhibited by all sources within the field (Smith et al. 2012; Stumpe et al. 2012, 2014). The light curves are corrected for dilution during the SPOC processing with TOI-1695 suffering low levels of contamination as indicated by its average dilution correction factor across all five sectors of 0.983. We only consider reliable TESS measurements for which the measurements quality flag is equal to zero. TOI-1695's PDCSAP light curve is depicted in Figure 1 and shows no compelling signs of coherent photometric variability, in particular from rotation.

Zoom In Zoom Out Reset image size

Following PDCSAP light-curve construction, the SPOC conducted a transit search using the Transiting Planet Search (TPS) module (Jenkins 2002; Jenkins et al. 2010). A repeating transit-like signal with a reported period of 3.13 days was detected in all five sectors independently. A total of 38 transit events were observed over the five TESS sectors and are highlighted in Figure 1. The signal passed a set of internal data validation (DV) tests (Twicken et al. 2018) and was fit with a preliminary limb-darkened transit model. After a review of the diagnostic tests in the DV reports, the TESS Science Office classified the planet candidate as TOI-1695.01 (Guerrero et al. 2021). The SPOC reported a preliminary R_p /R_s value of 0.0316, which corresponds to a planetary radius of 1.78 R_⊕ using our adopted stellar radius of 0.515 R_⊙ (see Section 3).

Active regions on the surface of a rotating star will introduce a time-varying signal in both photometric and in precise radial velocity (RV) measurements. Our ability to construct complete models of each of these data sets therefore benefits from a priori knowledge of the stellar rotation period P_rot. Because the TESS light curve does not exhibit any signature of stellar rotation (Figure 1), we queried the All-Sky Automated Survey for Supernovae (ASAS-SN; Shappee et al. 2014; Kochanek et al. 2017) data archive to search for long-term photometric monitoring of TOI-1695. ASAS-SN is a global network of 24 telescopes, hosted by the Las Cumbres Observatory, whose ongoing goal is to monitor the entire sky on a continuous basis to search for transient phenomena.

Our data archive search revealed that TOI-1695 was monitored throughout the ASAS-SN campaign in the V band for more than four years from 2014 July 8 to 2018 November 29 UT (Figure 2). The nightly cadence of the light curve is sufficient to detect rotational variability on timescales that exceed a few days, such as what we expect for TOI-1695 given the lack of photometric variability in its TESS light curve. We computed the generalized Lomb–Scargle (GLS) periodogram (Zechmeister & Kurster 2009) of these data and uncovered a strong periodicity at approximately 48 days that is not seen in the light curve's window function (not shown in Figure 2). We interpret this periodic signal as the stellar rotation period and fit the ASAS-SN photometry with a sinusoidal function as shown in the bottom row of Figure 2. From this fit we measure a photometric amplitude of 7.7 ppt and P_rot = 47.7 ± 2.2 days. We note that the GLS periodogram of the light-curve residuals does not show any signature of a significant residual periodicity. We also note that this value of P_rot is consistent with the expectation from population studies of inactive early-M dwarf rotation periods (e.g., Newton et al. 2016). We will use this measurement of P_rot as a prior in our forthcoming data analyses.

Zoom In Zoom Out Reset image size

We obtained a pair of reconnaissance spectra of TOI-1695 on 2020 February 2 and 13 UT using the Tillinghast Reflector Échelle Spectrograph (TRES; Fűrész 2008). These observations were coordinated as part of the TESS Follow-up Observing Program (TFOP). TRES is an R = 44,000 fiber-fed optical spectrograph (310–910 nm) mounted on the 1.5 m Tillinghast Telescope at the Fred Lawrence Whipple Observatory on Mount Hopkins, Arizona. The spectra were reduced and extracted following the standard procedure (Buchhave et al. 2010) and were subsequently cross-correlated against a custom template of Barnard's star over a range of $v\sin i$ values (Winters et al. 2018). These spectra reveal that TOI-1695 is single lined, lacks a measurable signal from rotational broadening (i.e., $v\sin i$ < 3.4 km s⁻¹), and exhibits Hα in absorption. Taken together, these findings are consistent with the absence of strong photometric variability and confirm that TOI-1695 is likely a slowly rotating and chromospherically inactive star.

We also measured stellar RVs at opposing quadrature phases of −59.766 ± 0.081 km s⁻¹ and −59.892 ± 0.056 km s⁻¹, indicating that there is no significant RV variation that would have been produced if the system were a spectroscopic binary. As such, the transit-like signal TOI-1695.01 remains a viable planet candidate that we continue to vet observationally in the following subsections.

With a pixel scale of 21'' pixel⁻¹, and photometric apertures that typically extend out to roughly $1^{\prime}$, TESS commonly produces light curves with multiple blended sources. Indeed, there are 68 sources from Gaia DR3 within 25 of TOI-1695, with TOI-1695 being the brightest source (ΔG = 2.693) and requiring only a marginal dilution correction in the TESS bandpass (i.e., 0.983). To resolve this, we acquired seeing-limited ground-based transit follow-up photometry of TOI-1695.01. These observations were taken as part of TFOP Sub Group 1 (TFOP SG1; Collins 2019) ³⁴  to rule out or identify nearby eclipsing binaries (NEBs) and to check for a transit-like event on target using the greater spatial resolution compared to the TESS images. Our data also constrain the chromatic transit depth across complementary optical filter bands. We used the TESS Transit Finder to schedule our transit observations and the photometric data were extracted using AstroImageJ (Collins et al. 2017).

2.4.1. MLO

2.4.2. LCOGT

We observed four predicted transit windows of TOI-1695.01 with the Las Cumbres Observatory Global Telescope (LCOGT; Brown et al. 2013) 1.0 m network node at McDonald Observatory in Texas, USA, on 2020 August 24, 2020 October 23, 2020 November 14, and 2020 December 9 UT. All but one observation was conducted in Pan-STARRS z_s band, with one of two observations on 2020 November 14 UT being taken in the Sloan $g^{\prime}$ band. The 1 m telescopes are equipped with 4096 × 4096 SINISTRO cameras having an image scale of 0389 per pixel, resulting in a $26^{\prime} \times 26^{\prime}$ field of view. The images were calibrated by the standard LCOGT BANZAI pipeline (McCully et al. 2018).

The first LCOGT observation on 2020 August 24 UT should have been sensitive to the 1.2 ppt event if it occurred on TOI-1695, but the data likely ruled out an event on or off target that would have been deep enough to cause the event detected by the SPOC pipeline. However, like the MLO observations one predicted orbit earlier, these observations turned out to be out-of-transit relative to the later multisector SPOC ephemeris. The four remaining LCOGT observations were conducted according to the precise multisector SPOC ephemeris and achieved continuous coverage across the full transit events. We used AstroImageJ to extract the photometric data using circular photometric apertures with radii in the range 43 to 58. All of the TOI-1695 apertures exclude flux from the nearest known Gaia DR3 and TESS Input Catalog neighbor (TIC 629325854) 146 east. We perform a least squares fit to the individual light curves using a combined transit plus systematics model. We construct the latter using a linear combination of time, the FWHM of the point-spread function (PSF), and the sky background.

We detect the transit event within the TOI-1695 photometric apertures in the three z_s -band and one $g^{\prime}$-band light curves (see Figure 3). We find that in all but the z_s -band light curve taken on 2020 October 23 UT, the measured values of R_p /R_s ∈ [0.033, 0.036] are consistent with the value measured by TESS. One exceptional light curve produced an anomalously large R_p /R_s = 0.045, which we consider to be an outlier and not a true chromatic effect due to its inconsistency with the other z_s -band light curves. Our results confirm that the transit event of TOI-1695.01 occurs on target such that we are able to rule out NEBs and continue to interpret TOI-1695.01 as a viable planet candidate.

Zoom In Zoom Out Reset image size

We note that we do not include these observations in our global transit analysis (Section 4.1) because of the strong dependence of each light curve on the exact systematics model used. Plus, because of the multiyear TESS baseline, which extends well beyond our most recent ground-based observation, the ground-based light curves presented herein do not provide stronger constraints on the planet candidate's orbital period when compared to TESS alone.

2.5.1. Keck/NIRC2

Following Ciardi et al. (2015), we assessed the possible contamination of the TESS light curve by bound or unbound companions using high-resolution adaptive optics (AO) imaging by NIRC2 on Keck II (Wizinowich et al. 2000). We observed TOI-1695 on 2020 May 28 UT in the narrowband Br–γ filter with an integration time of 4 seconds with one coadd per frame for a total of 36 seconds on target. Our AO data were processed and analyzed following the standard procedure (Furlan et al. 2017), which includes the calculation of the 5σ contrast curve via source injection (see Figure 4).

Zoom In Zoom Out Reset image size

We detect no additional companions around TOI-1695 given the sensitivity of our data. We demonstrate sensitivities down to ∼3.5 mag at 006 (2.7 au) and ∼7 mag at 05 (22 au). These contrast limits indicate that there are likely no stellar companions down to M6–L9 spectral types within 30 of the position of TOI-1695. From these results we conclude that TOI-1695.01 remains a viable planet candidate and is worthy of the time investment to obtain precise RV measurements for planet confirmation and planetary mass measurement.

2.5.2. SAI Speckle Polarimetry

We obtained 49 spectra of TOI-1695 with the HARPS-N optical échelle spectrograph at the 3.6 m Telescopio Nazionale Galileo on La Palma in the Canary Islands. HARPS-N has a resolving power of R = 115,000 and is stabilized in pressure and temperature, thus enabling submeter per second instrumental stability (Cosentino et al. 2012). We observed TOI-1695 over 428 days between 2020 December 6 UT and 2022 February 7 UT as part of the HARPS-N Guaranteed Time Observations program. We fixed the exposure time to 1800 s throughout the campaign, which yielded a median signal-to-noise ratio (S/N) per order of 15.9 across all orders redward of aperture 18 (440–687 nm).

We reduced the spectra using version 3.7 of the HARPS-N Data Reduction Software (DRS; Lovis & Pepe 2007), which includes automated RV extraction using an M0 template with a cross-correlation function (CCF). We opted to conduct a separate RV extraction using the template-matching algorithm TERRA (Anglada-Escude & Butler 2012), which has been shown to outperform the CCF method for M dwarfs (e.g., Anglada-Escude & Butler 2012; Astudillo-Defru et al. 2015). TERRA constructs an empirical master spectral template by coadding the individual HARPS-N spectra after being translated into the barycentric frame. We ignore spectral regions in which the atmospheric transmission is <99% and only consider échelle orders redward of aperture 18, following the recommended procedure for M dwarfs (Anglada-Escude & Butler 2012). From the remaining spectral regions, we compute the RV shift of each spectrum via least squares matching to the master template. We obtain a median RV uncertainty and a raw RV rms of 1.81 and 6.32 m s⁻¹, compared to values of 6.11 m s⁻¹ and 8.02 m s⁻¹, respectively, from the DRS.

We also measure logR'_HK = −4.74 ± 0.41 from our master template. From this we derive an expected P_rot of ${27.8}_{-13.2}^{+26.8}$ days from the activity–rotation relation from Astudillo-Defru et al. (2017), which is consistent with the measured P_rot from ASAS-SN. TOI-1695 is in the unsaturated chromospheric activity regime, which is consistent with the absence of observed broadening in the HARPS-N spectrum ($v\sin i$ < 1.3 km s⁻¹) and with its Hα being seen in absorption.

The raw RVs from TERRA are shown in the upper row of Figure 5. The corresponding GLS periodogram of the raw RVs clearly exhibits a strong periodicity at the 3.13 day period of the planet candidate TOI-1695.01. We therefore proceed by interpreting TOI-1695.01 as a validated planet, which we will now refer to as TOI-1695 b for the remainder of this paper.

Zoom In Zoom Out Reset image size

We first model the raw TESS PDCSAP light curve (Figure 1, top row) in which the planet candidate TOI-1695.01 was originally detected. The PDCSAP light curve has already undergone systematics corrections via a linear combination of cotrending basis vectors; however, some low-amplitude and temporally correlated signals that are unrelated to planetary transits persist. We model this residual systematic noise in the PDCSAP curve using a GP simultaneously with our transit model. We employ the exoplanet software package (Foreman-Mackey et al. 2019) to sample the posterior of the joint GP and transit model parameters at each step in our Markov chain Monte Carlo (MCMC) simulation. The exoplanet package uses the STARRY package (Luger et al. 2019) to compute analytical transit models and celerite (Foreman-Mackey et al. 2017) to evaluate the likelihood of the GP model.

We adopt a covariance kernel of the form of a stochastically driven, damped, simple harmonic oscillator in Fourier space. The power spectral density of the kernel is

$$\begin{eqnarray}&&S(\omega )=\sqrt{\displaystyle \frac{2}{\pi }}\displaystyle \frac{{S}_{0}{\omega }_{0}^{4}}{{\left({\omega }^{2}-{\omega }_{0}^{2}\right)}^{2}+{\omega }_{0}^{2}{\omega }^{2}/{Q}^{2}}.\end{eqnarray} \tag{ 1 }$$

The kernel is parameterized by the undamped period of the oscillator ρ = 2π/ω₀, where ω₀ is the undamped angular frequency (ω is the angular frequency); the standard deviation of the process ${\sigma }_{\mathrm{TESS}}=\sqrt{{S}_{0}{\omega }_{0}Q};$ and the fixed quality factor Q = 1/$\sqrt{2}$. Our GP model is jointly fit with a transit model for TOI-1695 b with the following free parameters: stellar mass M_s , stellar radius R_s , a reparameterization of the quadratic limb-darkening coefficients q1 and q2 (Kipping 2013a), orbital period P, time of midtransit t₀, log transit depth lnδ, baseline flux f_0,TESS, impact parameter b, eccentricity e, and argument of periastron ω. e is reparameterized with two shape parameters and sampled from a beta distribution ${ \mathcal B }$ as described in Kipping (2013b). Our full TESS transit model therefore contains the following 13 parameters: {lnρ, lnσ_TESS, M_s , R_s , q1, q2, lnP, t₀, lnδ, f_0,TESS, b, e, ω}.

We execute an MCMC simulation to sample the joint posterior probability density function (PDF) of our full set of model parameters using the PyMC3 MCMC package (Salvatier et al. 2016) within exoplanet. The MCMC simulation is initialized with two simultaneous chains, each with 1500 tuning steps and 1000 draws in the final sample. Point estimates of the maximum a posteriori values from the marginalized posterior PDFs of the GP hyperparameters are selected to construct the GP predictive distribution, whose mean function we adopt as our detrending model of the PDCSAP light curve. This mean detrending function and the detrended light curve are both shown in Figure 1. Similarly, we recover the maximum a posteriori point estimates of the transit model parameters to construct the transit model shown in the bottom panel of Figure 1. Median maximum a posteriori values and uncertainty point estimates from the 16th and 84th percentiles for all model parameters are reported in Table 2. Recall that we do not consider ground-based transits jointly with our TESS transit analysis because they are susceptible to residual systematic uncertainties and our longest transit baseline is spanned by TESS Sectors 18 through 52.

Table 2. Point Estimates of the TOI-1695 Model Parameters

Parameter	Value
Transit Parameters
Baseline flux, f0,TESS [ppt]	${0.017}_{-0.040}^{+0.042}$
Limb-darkening coefficient, q1	${0.32}_{-0.48}^{+0.53}$
Limb-darkening coefficient, q2	${0.33}_{-0.48}^{+0.51}$
lnρ [days]	${1.00}_{-0.11}^{+0.12}$
lnσTESS	−0.91${}_{-0.07}^{+0.08}$
RV Parameters
Prot [days]	${48.9}_{-2.7}^{+2.2}$
σrot	${5.83}_{-0.92}^{+0.94}$
lnQ0	${1.4}_{-0.9}^{+1.0}$
lndQ	${0.049}_{-1.984}^{+2.026}$
f	${0.64}_{-0.31}^{+0.25}$
RV offset, v0 [m s−1]	−0.91${}_{-0.88}^{+0.91}$
Jitter, s [m s−1]	${2.32}_{-0.98}^{+0.94}$
TOI-1695 b Parameters
Orbital period, P [days]	${3.1342791}_{-0.0000063}^{+0.0000071}$
Time of midtransit, T0 [BJD—2,457,000]	${1791.52056}_{-0.00111}^{+0.00098}$
Transit duration, D [hr]	${1.311}_{-0.189}^{+0.324}$
Transit depth, δ [ppt]	1.236 ${}_{-0.081}^{+0.083}$
Semimajor axis, a [au]	${0.033548}_{-0.000268}^{+0.000260}$
Planet-to-star radius ratio, Rp /Rs	0.034 ± 0.002
Impact parameter, b	${0.69}_{-0.31}^{+0.11}$
Inclination, i [deg]	${87.2}_{-0.5}^{+1.3}$
Eccentricity, e	<0.097 a
Planet radius, Rp [R⊕]	${1.90}_{-0.14}^{+0.16}$
RV semiamplitude, K [m s−1]	4.39 ± 0.69
Planet mass, Mp [M⊕]	6.36 ± 1.00
Bulk density, ρp [g cm−3]	${5.0}_{-1.3}^{+1.8}$
Surface gravity, gp [m s−2]	${17.1}_{-3.4}^{+4.3}$
Escape velocity, vesc [km s−1]	20.5 ± 1.8
Instellation, F [F⊕]	39 ± 3
Equilibrium temperature, Teq [K]	698 ± 14 b

Notes.

^a 95% confidence interval. ^b Zero albedo assumed.

Download table as:  ASCIITypeset image

We impose strong priors on P and T₀ derived from our TESS analysis. The raw RVs and their GLS periodograms are shown in Figure 5. The periodic signal induced by the orbit of TOI-1695 b is clearly visible at 3.13 days. The rotation period P_rot = 48 days of TOI-1695 is well constrained from 4 yr of ASAS-SN data (see Figure 2), and a moderate signal is observed near this period in the GLS periodogram, consistent with the star being inactive. We simultaneously fit the observed data with a Keplerian orbit and a quasi-periodic GP regression model of stellar activity. We adopt a GP covariance kernel of the form of two stochastically driven, damped, simple harmonic oscillators in Fourier space, both described by Equation (1). The parameters of the two simple harmonic oscillator terms are

$$\begin{eqnarray}&&{Q}_{1}=1/2+{Q}_{0}+{dQ},\end{eqnarray} \tag{ 2 }$$

$$\begin{eqnarray}&&{Q}_{2}=1/2+{Q}_{0},\end{eqnarray} \tag{ 3 }$$

$$\begin{eqnarray}&&{\omega }_{1}=\displaystyle \frac{4\pi {Q}_{1}}{{P}_{\mathrm{rot}}\sqrt{4{Q}_{1}^{2}-1}},\end{eqnarray} \tag{ 4 }$$

$$\begin{eqnarray}&&{\omega }_{2}=\displaystyle \frac{8\pi {Q}_{1}}{{P}_{\mathrm{rot}}\sqrt{4{Q}_{1}^{2}-1}},\end{eqnarray} \tag{ 5 }$$

$$\begin{eqnarray}&&{S}_{1}=\displaystyle \frac{{\sigma }_{\mathrm{rot}}^{2}}{(1+f){\omega }_{1}{Q}_{1}},\end{eqnarray} \tag{ 6 }$$

$$\begin{eqnarray}&&{S}_{2}=\displaystyle \frac{f{\sigma }_{\mathrm{rot}}^{2}}{(1+f){\omega }_{2}{Q}_{2}}.\end{eqnarray} \tag{ 7 }$$

This kernel has two modes in Fourier space: one at period P and one at half the period. We parameterize the kernel by the primary period of variability P_rot; the standard deviation of the process σ_rot; the quality factor (minus one half) for the secondary oscillation Q₀, which keeps the system underdamped; the difference between the quality factors of the first and second modes dQ, which ensures that the primary mode always has higher quality; the fractional amplitude f of the secondary mode compared to the primary; and finally an additive scalar jitter term s is added to account for any excess noise in the activity model. The GP only models correlated noise, so the jitter term is added to capture the uncorrelated noise and is added in quadrature to the RV uncertainty.

The RV Keplerian model is parameterized by the orbital period P, the time of midtransit T₀, the log of the semiamplitude $\mathrm{ln}K$, an RV offset v₀, and the following reparameterizations of eccentricity e and argument of periastron ω: $h=\sqrt{e}$ cosω and $k=\sqrt{e}$ sinω. Hence, our full RV model consists of our GP and Keplerian models with 12 free parameters: {P_rot, σ_rot, lnQ₀, lndQ, f, lns, v₀, P, T₀, lnK, h, k}. The adopted model parameter priors are also included in Table 3. We fit the RV data with our full model using the exoplanet package, which is an extension of the PyMC3 inference engine. Point estimates of the model parameters are derived from their respective marginalized PDFs and are reported in Table 2. The point estimates reported represent each parameter's median maximum a posteriori value and uncertainties from its 16th and 84th percentiles.

Table 3. TESS Light-curve and RV Model Parameter Priors

Parameter	Fiducial Model Priors
Stellar Parameters
Ms [M⊙]	${ \mathcal N }$(0.513, 0.012)
Rs [R⊙]	${ \mathcal N }$(0.515, 0.015)
Light-curve Parameters
q1	${ \mathcal U }(0,1)$ a
q2	${ \mathcal U }(0,1)$ a
lnρ [days]	${ \mathcal U }(2,20)$
lnσTESS	${ \mathcal N }$(ln(std( f TESS)), 10) b
lnsTESS	${ \mathcal N }$(ln(std( f TESS)), 10) b
f0,TESS	${ \mathcal N }$(0, 2)
GP and RV Parameters
Prot [days]	${ \mathcal N }$(48, 3)
σrot	${ \mathcal N }$(std(data), 1)
lnQ0	$| { \mathcal N }(0,2)|$
lndQ	${ \mathcal N }(0,2)$
f	${ \mathcal U }(0,1)$
lnsRV [m s−1]	${ \mathcal U }(-3,1)$
γ [m s−1]	${ \mathcal N }$(0, 2)
Planetary Parameters
P [days]	${ \mathcal N }(3.1343,0.0001)$ c
T0 [BJD—2,457,000]	${ \mathcal N }(1791.518,0.001)$ c
lnK [m s−1]	${ \mathcal U }(-1,4)$
lnδ [ppt]	${ \mathcal U }(-2,2)$ d
b	${ \mathcal U }(0,1)$
e	${ \mathcal B }(0.867,3.03)$ e
	half ${ \mathcal N }(0,0.32)$ f
ω	${ \mathcal U }(-\pi ,\pi )$

Notes.

^a Reparameterization of the limb-darkening coefficients, u₁ and u₂ described by Kipping (2013a). ^b f _TESS is the normalized PDCSAP flux. ^c Priors for RV analysis were based on the TESS transit analysis. We did not impose any prior during the TESS analysis. ^d Informed by a BLS search and fed into R_p /R_s . ^e TESS transit analysis, Kipping (2013b). ^f RV analysis, Van Eylen et al. (2019). Other prior values and parameterizations were explored (including ecosω and esinω) and revealed negligible differences in point estimates.

Download table as:  ASCIITypeset image

Figure 5 shows the raw RVs, individual model components, and GLS periodograms. A noticeable but not significant periodicity at P_rot = 48 days emerges in the GLS periodogram of the RV activity signal, and also at P_rot/2 = 24 days to a lesser extent. The GLS periodogram corresponding to the orbital model solution after removing the mean GP activity model is clearly dominated by the 3.13 day periodicity with no other significant periodic signal. We measure an RV semiamplitude of K = 4.39 ± 0.69 m s⁻¹, which is clearly visible in the phase-folded RVs also shown in Figure 5. The RV residuals show the data minus the mean GP activity model and maximum a posteriori Keplerian solution, and have an rms of 3.47 m s⁻¹ and a reduced χ² of 1.39. The GLS periodogram of the residuals lack signals with significant power, which suggests that we do not have evidence for an additional planetary companion in the TOI-1695 system.

We also performed a search for additional transiting planets in the TESS light curve using the Transit Least Squares (TLS) algorithm (Hippke & Heller 2019). We conducted the search on the detrended TESS light curve following the removal of our maximum a posteriori transit model for TOI-1695 b. We ran separate TLS transit searches on each set of consecutive TESS sectors and over the range of orbital periods from 0.5–30 days. Our TLS search revealed no significant signals indicating that we do not have evidence for repeating transits from an additional planet in the system. We also conducted a close visual inspection of the TESS light curve, which revealed no obvious signature of any single-transit events. Our null detection of a second planet is consistent with the results of the SPOC transit search as indicated in the DV report.

We also attempted to measure the individual transit times of TOI-1695 b in each of the five TESS sectors to search for transit timing variations (TTVs). However, given the low S/N of any individual transit event in the TESS light curve, ³⁵  we were unable to recover 70% of the individual midtransit times confidently. For the remaining transit events for which our transit model did converge, we measured the individual T_0,i values with a typical uncertainty of ∼15 minutes. These measurements did not reveal any significant TTVs when compared to a linear ephemeris.

Our analysis of the PDCSAP light curve reveals that TOI-1695 b has an orbital period of $P={3.1342791}_{-0.0000063}^{+0.0000071}$ days. Using the stellar parameters presented in Table 1, we find that the semimajor axis of the planetary orbit is $a\,={0.033516}_{-0.000086}^{+0.000087}$ au where it receives an instellation flux of F = 39 ± 3 F_⊕. Assuming uniform heat redistribution and a Bond albedo of zero, the corresponding equilibrium temperature of TOI-1695 b is T_eq = 698 ± 14 K.

We also measure a planetary radius and mass of ${R}_{p}={1.90}_{-0.14}^{+0.16}\ {R}_{\oplus }$ and M_p = 6.36 ± 1.00 M_⊕, corresponding to 12.7σ and 6.4σ detections, respectively. The mass and radius measurements yield a 2.8σ bulk density measurement of ${\rho }_{p}={5.0}_{-1.3}^{+1.8}$ g cm⁻³.

Figure 7 compares the mass and radius of TOI-1695 b to the population of small M dwarf planets with masses measured to better than 3σ, which we retrieve from the NASA Exoplanet Archive. Planets are classified based on their bulk compositions inferred from their mass and radius measurements. Earth-like planets are defined as those consistent with an Earth-like composition curve, gas-rich planets cannot be explained by even 100% water composition and require an extended H/He envelope, and the remaining planets we broadly classify as "intermediate" given that their masses and radii are between those of the aforementioned groups. We find that the bulk composition of TOI-1695 b is inconsistent with that of the Earth at 1.7σ. Specifically, TOI-1695 b is underdense relative to an Earth-like composition of the same mass. The planet could therefore belong to the population of enveloped terrestrials whose rocky components resemble Earth but require extended gaseous envelopes to explain their masses and radii. Such is the expected composition of the majority of sub-Neptunes predicted by thermally driven mass-loss models (Owen & Wu 2017; Gupta & Schlichting 2019). In Section 5.2, we show that TOI-1695 b would require an envelope mass fraction of 0.06% to explain its mass and radius. Assuming an H/He envelope with solar metallicity (mean molecular mass = 2.35 u), we also show that the planet is unlikely to retain such an envelope over long timescales since it is susceptible to thermally driven hydrodynamic escape at 39 times Earth insolation. Thus, it is unlikely to be an enveloped terrestrial with a primordial H/He atmosphere. Instead, it is likely that TOI-1695 b possesses a high mean molecular weight envelope supplied by volatile delivery, or formation beyond the ice line followed by migration and volatile retention. We find that the mass and radius of TOI-1695 b is consistent with a water mass fraction of ${31}_{-22}^{+33}$% assuming a two layer model of MgSiO₃ and H₂O. This composition is consistent within 1σ of the subpopulation of water worlds with water mass fractions of 50% (see Section 5.4; Luque & Palle 2022).

Zoom In Zoom Out Reset image size

We ran photoevaporation simulations to assess the possibility that TOI-1695 b could be an Earth-like core enveloped by an H/He-dominated atmosphere. Our model assumes that extreme ultraviolet (EUV; 10 nm < λ < 130 nm)-driven hydrodynamic escape is the primary driver of mass loss. We opt to ignore X-ray irradiation since approximately 80%–95% of the high-energy flux is in the EUV and to remain consistent with the available synthetic spectra discussed below (Fontenla et al. 2016; Loyd et al. 2016; Peacock et al. 2019a,2019b, 2020). We also assume an H/He composition with solar metallicity (mean molecular mass = 2.35 u) and assume that hydrogen is atomic due to photodissociation.

In order to model the time evolution of the incident EUV flux, we constructed a power-law function of the EUV flux's temporal evolution. Our model is informed by semi-empirical spectra generated by the HAZMAT team for populations of M1 stars over a range of ages from 10 Myr to 5 Gyr (Figure 9; Peacock et al. 2020)

$$\begin{eqnarray}\begin{array}{c}{F}_{{\rm{EUV}}}=\left\{\begin{array}{cc}{F}_{{\rm{EUV}},0} & t\lt {t}_{{\rm{sat}}},\\ {F}_{{\rm{EUV}},0}{\left(\displaystyle \frac{t}{{t}_{{\rm{sat}}}}\right)}^{\beta } & t\gt {t}_{{\rm{sat}}},\end{array}\right.\end{array}\end{eqnarray} \tag{ 8 }$$

where t_sat is the saturation time marking the transition from a constant EUV flux phase to a power-law decay and F_EUV is the incident planetary EUV flux (related to the stellar EUV luminosity L_EUV by F_EUV = L_EUV/4π a², where a is the orbital semimajor axis). We found β = −1.23 and t_sat = 500 Myr to be consistent with the HAZMAT spectra as well as previously reported solar data (Ribas et al. 2005).

EUV-driven escape generates a mass flux, which we compute as a function of an efficiency factor , incident EUV flux F_EUV, and the planetary gravitational potential V_pot

$$\begin{eqnarray}&&\phi =\displaystyle \frac{\epsilon {F}_{\mathrm{EUV}}}{4{V}_{\mathrm{pot}}},\end{eqnarray} \tag{ 9 }$$

where ${V}_{\mathrm{pot}}={{GM}}_{p}/{R}_{p}^{2}$. encapsulates several heat-transfer processes, ultimately representing the fraction of incident radiation that drives escape. We chose a value of = 0.15, consistent with previously reported lower estimates (Watson et al. 1981; Schaefer et al. 2016).

The evolution of the planetary radius through time was determined by a combination of atmospheric escape described by ϕ and contraction due to cooling described by thermal evolution models (Lopez & Fortney 2014). Our simulations demonstrate that TOI-1695 b loses its entire atmosphere on rapid timescales of the order 1–100 Myr for a wide range of initial atmospheric mass fraction values (0.06%–5%; Figure 8). This finding is robust over a wide range of values for t_sat and . We find that the observed planetary radius is consistent with a 0.06% H/He-dominated atmosphere, which the planet could not have retained over reasonably observable periods. We conclude that a present-day H/He-dominated atmosphere is highly unlikely, barring a steady-state outgassing scenario.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Several physical mechanisms have been proposed to explain the emergence of the radius valley around FGK and late-K to mid-M dwarfs, each making distinct predictions of the slope of the radius valley in period–radius space parameterized as R_p,valley ∝ P^β . These include two thermally driven atmospheric mass-loss models: stellar X-ray and ultraviolet (XUV)-driven hydrodynamic escape in which atmospheric species flow outward in the form of a Parker wind (Owen & Wu 2013; Jin et al. 2014; Lopez & Fortney 2014; Chen & Rogers 2016; Jin & Mordasini 2018); and core-powered mass loss in which the planetary core's formation energy drives escape over gigayear timescales (Ginzburg et al. 2018; Gupta & Schlichting2019, 2020). An alternative mechanism is the gas-poor (but not gas depleted) formation scenario which suggests a primordial radius valley in which gas accretion onto low-mass cores is limited (≲1–2 M_⊕; Lee & Connors 2021; Lee et al. 2022). These three models predict negative values of β ∈ [−0.15, −0.09] (Lopez & Rice 2018; Gupta & Schlichting 2020; Lee & Connors 2021). Yet a fourth model supposes that enveloped terrestrials form within the first several megayears when the gaseous disk is still present, whereas terrestrial planets form at later times after the dissipation of the gaseous disk in a gas depleted environment. The gas depleted model predicts the opposite sign for the period dependence of the radius valley (β = 0.11; Lopez & Rice 2018).

The slope of the radius valley around Sun-like stars with T_eff > 4700 K has been well characterized and measurements of β take values of [−0.11, −0.09] (Van Eylen et al. 2018; Martinez et al. 2019), consistent with thermally driven mass loss and gas-poor formation model predictions. However, in the lower stellar mass regime of late-K to mid-M dwarfs, tentative evidence suggests that β = 0.06 ± 0.02, which is inconsistent with the values measured around FGK stars and is instead consistent with predictions from gas depleted formation models (Cloutier & Menou 2020).

While occurrence rate measurements around low-mass stars are currently insufficient to resolve a transition in competing radius valley emergence mechanisms between the low stellar mass and high stellar mass regimes, we may gain insight by characterizing keystone planets like TOI-1695 b, which span the model predictions in period–radius space. The distinct slopes of the radius valley's period dependence carve out a wedge in the orbital period–planet radius space, between which the competing models make conflicting predictions (shaded region in Figure 10). At periods less than 23.5 days, thermally driven mass loss and gas-poor formation models predict that planets in the wedge (i.e., keystone planets) are rocky. Conversely, gas depleted formation models predict they should host gaseous envelopes because their sizes exceed the maximum rocky planet mass that can form out of the minimum-mass extrasolar nebula at its observed orbital separation. At P = 3.13 days and R_p = ${1.90}_{-0.14}^{+0.16}$ R_⊕, TOI-1695 b is one such keystone planet, whose composition directly constrains the prominence of the competing physical processes on close-in planets around early-M dwarfs (Figure 10).

Zoom In Zoom Out Reset image size

Figure 10 also features the same population of planets displayed in Figure 7. Intermediate planets may be explained by a variety of compositions including an H/He envelope, a volatile-rich composition, or perhaps an exotic rocky composition that is enhanced in Ca and Al (Dorn et al. 2018). Our analysis reveals that TOI-1695 b is inconsistent with Earth-like and gas-rich compositions and requires an alternative physical interpretation to explain its mass and radius. The bulk composition of TOI-1695 b is therefore inconsistent with predictions from thermally driven mass-loss models.

So what about the prospect that TOI-1695 b formed via a gas depleted formation scenario? If the gas depleted formation mechanism were operating in the TOI-1695 system, then TOI-1695 b must have formed early on before disk dispersal and it subsequently accreted—and at least partially retained—a primordial H/He envelope. However, if TOI-1695 b accreted its primordial envelope at its current location, then such an atmosphere should have been rapidly lost to thermal escape, as demonstrated in Section 5.2. This discrepancy may be reconciled if the planet hosts a high mean molecular weight atmosphere and/or migrated inward to its current location so as to avoid the bulk stellar XUV output during the first few hundreds of megayears. We therefore conclude that the classical picture of the gas depleted formation model, which only produces gas-enveloped terrestrials and terrestrial cores that are born rocky, cannot explain the observed composition of TOI-1695. Instead, we speculate that TOI-1695 b is more likely to be rich in volatiles. However, we emphasize that our conclusion that TOI-1695 b is inconsistent with a thermally driven mass-loss scenario remains.

At the time of this publication, there are seven M dwarf keystone planets for which reliable mass and radius measurements are available (Table 4). Using this sample, we ask the question: what is the probability that the keystone planet population around M dwarfs is sculpted by a thermally driven mass-loss process (P(TDML))? To answer this question, we must first define the probability that the composition of each keystone planet is consistent with a thermally driven mass-loss hypothesis (P_i ). For this we compute the probability that each planet is consistent with having an Earth-like composition based on its mass and radius. We compute each P_i as the fraction of samples from the planet's joint mass–radius posterior that result in a sampled radius that is less than the radius of a pure MgSiO₃ planet at the sampled mass value. This criterion adopted to define the radius upper limit for an Earth-like planet is equivalent to our definition described in Section 5.1 and shown in Figure 7. Our probabilities P_i are included in Table 4.

Table 4. M Dwarf Keystone Planet Parameters

Planet	Stellar Mass	Orbital Period	Planet Radius	Planet Mass	Ref.	Probability of Consistency with
Name	(M⊙)	(days)	(R⊕)	(M⊕)		an Earth-like Composition, Pi
TOI-1235 b	0.640 ± 0.016	3.445	${1.738}_{-0.076}^{+0.087}$	${6.91}_{-0.85}^{+0.75}$	1	0.76
TOI-776 b	0.544 ± 0.028	8.247	1.85 ± 0.13	4.0 ± 0.9	2	0.03
TOI-1695 b	0.513 ± 0.012	3.134	${1.90}_{-0.14}^{+0.16}$	6.36 ± 1.00	3	0.20
TOI-1634 b	0.502 ± 0.014	0.989	${1.790}_{-0.081}^{+0.080}$	${4.91}_{-0.70}^{+0.68}$	4	0.09
TOI-1685 b	0.495 ± 0.019	0.669	1.70±0.07	3.78 ± 0.63	5	0.05
G 9-40 b	0.295 ± 0.014	5.746	1.90 ± 0.07	4.00 ± 0.63	6	7.0 × 10−5
TOI-1452 b	0.249 ± 0.008	11.062	1.67±0.07	4.82 ± 1.30	7	0.42

References. (1) Cloutier et al. (2020); (2) Luque et al. (2021); (3) this work; (4) Cloutier et al. (2021b); (5) Bluhm et al. (2021); (6) Luque et al. (2022); and (7) Cadieux et al. (2022).

Download table as:  ASCIITypeset image

Although the probability that any individual keystone planet is consistent with thermally driven mass loss is often not very illuminating, the statistical statement that we can make from the seven keystone planet sample is meaningful. By treating the measurement of each keystone planet's mass and radius as an independent Bernoulli experiment, we can calculate the probability that the keystone planet population is being sculpted by thermally driven mass loss as the product of the individual probabilities (i.e., $P(\mathrm{TDML})={\prod }_{i=1}^{7}{P}_{i}$). That is, what is the probability that all seven keystone planets have Earth-like compositions and are therefore consistent with a thermally driven mass-loss scenario? We find that P(TDML) = 5.9 × 10⁻¹⁰. Similarly, we evaluate the probability that the keystone planet population is inconsistent with a thermally driven mass-loss scenario and find that $P(\mathrm{not}\ \mathrm{TDML})={\prod }_{i=1}^{7}(1-{P}_{i})=9.4\,\times \,{10}^{-2}$. Comparing these values, we find P(not TDML)/P(TDML) = 1.6 × 10⁸. Thus, the M dwarf keystone planet population strongly disfavors a thermally driven mass-loss scenario. This result is not surprising as only one planet out of seven (i.e., TOI-1235 b) has a greater than 50% chance of having an Earth-like composition and therefore with being consistent with a thermally driven mass-loss process.

Our results are consistent with the emerging picture that the M dwarf radius valley is a by-product of planet formation and is not sculpted by thermally driven mass loss. Early investigations of close-in planet occurrence rates around late-K to mid-M dwarfs suggested that the slope of the radius valley with instellation was inconsistent with predictions from thermally driven mass-loss models (Cloutier & Menou 2020). More recent empirical evidence from Luque & Palle (2022) demonstrated that the sub-Neptune peak in the M dwarf radius valley represents water-rich planets and not gas-enveloped terrestrials. Luque & Palle (2022) refined the masses and radii of small transiting planets around M dwarfs and revealed three distinct planet types: Earth-like, water worlds, and puffy sub-Neptunes. The mass–radius profiles of these subpopulations are consistent with Earth-like, 50% water-dominated ices/50% silicates, and H/He-enveloped compositions respectively, and are interpreted as such. It is concluded that rocky planets must form within the water ice line while water worlds form beyond and migrate inward. In this case, the apparent radius valley around M dwarfs is sculpted by accretion history rather than by atmospheric mass loss. If water worlds are indeed ubiquitous, it is likely that TOI-1695 b belongs to this subpopulation. Our calculated water mass fraction of ${31}_{-22}^{+33}$% supports this interpretation and is consistent with the 50% water mass fraction subpopulation from Luque & Palle (2022). Taken together, the slope of the M dwarf radius valley (Cloutier & Menou 2020), the recovery of a population of likely water-rich planets (Luque & Palle 2022), and our results for the keystone planets suggest that thermally driven mass loss does not explain the origin of the M dwarf radius valley. Instead, it is likely to emerge directly from the planet-formation process.

Exoplanet transit surveys and RV follow-up of transiting systems have shown that multiplanet systems are common around M dwarfs (Dressing & Charbonneau 2015; Gaidos et al. 2016; Cloutier et al. 2021a). It is therefore reasonable to expect additional planets in the TOI-1695 system that evade detection due to their small sizes, long orbital periods, or nontransiting orbital configurations. We assessed the detection sensitivity of our HARPS-N RV data set to place constraints on the presence of additional planets by computing our detection sensitivity as a function of orbital period and planet mass via a set of injection–recovery tests. We took a Monte Carlo approach by injecting synthetic Keplerian signals into the residuals of the HARPS-N RV time series after removing the maximum a posteriori RV solution (i.e., TOI-1695 b plus GP). We inject a single planet in each of the 10⁵ iterations. To generate the Keplerian signal for each iteration, planet masses and orbital periods were sampled uniformly in log space from 1–20 M_⊕ and 1–200 days, respectively. Note that our RV baseline is 428 days. Because compact multiplanet systems are often nearly coplanar, we sampled orbital inclinations from a Gaussian distribution ${ \mathcal N }$(i_b , σ_i ), where i_b = 885 and we adopt a dispersion of mutual inclinations of σ_i = 2° following from studies of multiplanet M dwarf systems (Ballard & Johnson 2016). We sampled the stellar mass from its posteriors and used it to calculate the RV semiamplitude assuming a circular orbit. We injected the resulting Keplerian signals into the RV residuals while preserving the individual measurement uncertainties and timestamps.

Recovery of the injected synthetic planets corresponds to a successful detection and involved a two-step process. To warrant a detection, an injected signal must first produce a significant peak in a GLS periodogram with a false alarm probability (FAP) of ≤1%. The GLS periodogram was constructed for each iteration and the analytical FAP was calculated using the analytical formalism described by Zechmeister & Kurster (2009). Second, the six-parameter Keplerian model must be strongly favored over the null hypothesis (i.e., a flat line with a constant offset). To perform the model comparison, we calculated the Bayesian information criterion (BIC) for each model as BIC $=\,2\mathrm{ln}{ \mathcal L }+x\mathrm{ln}N$, where ${ \mathcal L }$ is the likelihood of the RV data given the assumed model, x is the number of model parameters (i.e., one and six for the null and Keplerian models, respectively), and N = 49 is the number of RV measurements. We claimed the successful recovery of an injected planet signal if both criteria are met: the GLS periodogram power of the largest signal within ± 2% of the injected period has FAP ≤ 1%, and ΔBIC = BIC_Keplerian − BIC_null ≥ 10. The sensitivity of our RV data set is defined as the ratio of recovered planets to injected planets and is presented in Figure 11.

Zoom In Zoom Out Reset image size

Figure 11 shows that the mass and orbital period of TOI-1695 b lay above the 90% detection contour. We also find that we are sensitive to approximately 50% of planets at 3 M_⊕ and 90% of planets at 4 M_⊕ at a 1 day orbital period. Within 10 days, we are sensitive to all planets ≥10 M_⊕. Adopting the empirical recent Venus and early Mars habitable zone (HZ) limits from Kopparapu et al. (2013), i.e., 35–150 days, we find that we are only sensitive to very massive HZ planets (>15 M_⊕). Such planets would likely host massive gaseous envelopes, rendering their surfaces uninhabitable by the traditional definition of the HZ.

Additionally, we performed a blind search over a wide period space with RVSearch, which revealed no significant signals below the 0.1% FAP threshold (Rosenthal et al. 2021). We conclude that additional planetary signals beyond TOI-1695 b are not detectable in our RV data.

Following the announcement of the TOI-1695.01 level-one planet candidate, multiple precise RV instrument teams began following up this target through TFOP. In this study we have presented the efforts from the HARPS-N Collaboration and we acknowledge that another collaboration is also in the process of presenting their own RV time series and analysis using data from the SPIROU near-IR spectropolarimeter (Kiefer et al. 2023). While the submissions of these complementary studies were coordinated between the two groups, their respective data, analyses, and write ups were intentionally conducted independently.