TOI-1749: an M dwarf with a Trio of Planets including a Near-resonant Pair

TOI-1749 (TIC 233602827) was observed by TESS with a 2 minutes cadence in 12 TESS sectors, specifically, Sectors 14–21 (from 2019 July 18 to 2020 February 18) and Sectors 23–26 (from 2020 March 18 to 2020 July 4). The coordinates and magnitudes of the target star are summarized in Table 1. The collected data were processed with a pipeline developed by the TESS Science Processing Operations Center (SPOC) at NASA Ames Research Center (Jenkins et al. 2016), from which two transit signals with orbital periods of 4.49 and 9.05 days were detected (Jenkins 2002; Jenkins et al. 2010; Li et al. 2019) and were labeled as TOI-1749.01 and TOI-1749.02, respectively (Guerrero et al. 2021). In addition, we detected an additional transiting planetary candidate interior to these two candidates as we describe in Section 3.2.1, which we internally labeled as TOI-1749.03. As we describe in Section 3.4, we validate these three planetary candidates as bona fide planets. Hereafter we designate TOI-1749.03, TOI-1749.01, and TOI-1749.02 as TOI-1749b, TOI-1749c, and TOI-1749d, respectively, in order of orbital period.

Table 1. Properties of the Host Star

Parameter	Value	Sources
Stellar name
TOI	1749	(1)
TIC	233602827	(2)
Astrometric and kinematic information
α	18:50:56.93	(2)
δ	64:25:10.08	(2)
Distance (pc)	99.56 ± 0.12	(2)
Parallax (mas)	10.0582 ± 0.0095	(3)
μα cos δ (mas yr−1)	−54.347 ± 0.012	(3)
μδ (mas yr−1)	+61.844 ± 0.012	(3)
U (km s−1)	−21.38 ± 0.17	This work
V (km s−1)	−16.44 ± 2.02	This work
W (km s−1)	27.96 ± 0.92	This work
Magnitudes
TESS	12.2574 ± 0.0073	(2)
G	13.1798 ± 0.0003	(3)
V	13.86 ± 0.072	(2)
J	11.069 ± 0.023	(4)
H	10.446 ± 0.021	(4)
K	10.270 ± 0.019	(4)
Physical parameters
Sp. type	M0V ± 0.5 subtype	This work
Mass ( M⊙)	0.58 ± 0.03	This work
Radius (R⊙)	0.55 ± 0.03	This work
Teff (K)	3985 ± 55	This work
$\mathrm{log}g$ (cgs)	4.70 ± 0.05	This work
[Fe/H] (dex)	−0.26 ± 0.08	This work
Age (Gyr)	>0.8	This work

Note. Sources: (1) Guerrero et al. (2021), (2) TIC v8 (Stassun et al. 2019), (3) Gaia EDR3 (Gaia Collaboration et al. 2021), (4) 2MASS (Skrutskie et al. 2006).

Download table as:  ASCIITypeset image

For further light-curve analyses, we extracted the Presearch Data Conditioning Simple Aperture Photometry (PDC-SAP; Smith et al. 2012; Stumpe et al. 2012, 2014) from the Mikulski Archive for Space Telescopes (MAST) at the Space Telescope Science Institute (STScI). Figure 1 shows an example of the TESS target pixel file (TPF) data around the target star, in which the aperture mask used by the pipeline to extract the light curve is indicated by orange-outlined pixels. There are two other sources of the Gaia DR2 catalog in the aperture masks (labeled as 2 and 3 in Figure 1), the fluxes of which were subtracted from the PDC-SAP fluxes, which are in absolute flux scale. We removed data points that are flagged as suspect in regard to quality from the PDC-SAP light curve, and then normalized it so that the median flux value in each sector is unity.

Zoom In Zoom Out Reset image size

To confirm the transit signals detected by TESS and measure TTVs in the system, we conducted photometric (or spectrophotometric) observations of predicted transits of the TOI-1749 planets using various ground-based telescopes. The telescopes and instruments that we used for these observations are listed in Table 2, and the observed transits are summarized in Table 3. In the following subsections, we briefly describe observations and data reductions for data that are used in the subsequent analyses. Other observations that are not used in the subsequent analyses are described in the Appendix.

Table 2. List of Telescopes and Instruments Used for Ground-based Transit Observations

Telescope/Instrument	Aperture	Observatory	Pixel scale	FoV	# of imaging channels
	(m)		(pixel−1)	(or slit size)	(or spectral resolution)
TCS/MuSCAT2	1.52	Teide	044	74 × 74	4
LCO 1 m/Sinistro	1.0	McDonald	0389	265 × 265	1
LCO 2 m/Spectral	2.0	Haleakala	0304	105 × 105	1
			(2 × 2 binned)
LCO 2 m/MuSCAT3	2.0	Haleakala	0266	91 × 91	4
GMU 0.81 m	0.81	GMU	034	$23^{\prime} \times 23^{\prime}$	1
NOT/ALFOSC	2.56	Roque de los Muchachos	0214	73 × 73	1
OAA 0.4 m	0.4	OAA	144	$36^{\prime} \times 36^{\prime}$	1
GTC/OSIRIS	10	Roque de los Muchachos	0254	($40^{\prime\prime} \times 7\buildrel{\,\prime}\over{.} 4$)	(1122)
			(2 × 2 binned)

Download table as:  ASCIITypeset image

2.2.1. TCS 1.52 m / MuSCAT2

2.2.2. LCO 1 m / Sinistro

2.2.3. LCO 2 m / Spectral

2.2.4. LCO 2 m / MuSCAT3

On 2020 September 5 UT, we obtained the optical low-resolution spectrum of TOI-1749 with the Alhambra Faint Object Spectrograph and Camera (ALFOSC) on the 2.56 m Nordic Optical Telescope (NOT) at Roque de los Muchachos Observatory on La Palma in the Canary Islands under the observing program 59–210. ALFOSC is equipped with a 2048 × 2064 CCD detector with a pixel scale of 02138 pixel⁻¹. We used grism number 5 and a horizontal long slit with a width of 10, which yield a nominal spectral dispersion of 3.53 Å pixel⁻¹ and a usable wavelength space coverage between 5000 and 9100 Å. The ALFOSC spectrum was acquired with an exposure time of 1800 s at the parallactic angle and an airmass of 1.26. On the same night, we acquired an ALFOSC spectrum of the spectrophotometric standard star G191–B2B (white dwarf) with the same instrumental setup as TOI-1749, at an airmass of 1.18 and exposure time of 180 s. Raw ALFOSC images were reduced following standard procedures at optical wavelengths: bias subtraction, flat-fielding using spectral halogen flat fields, and optimal extraction using appropriate packages within the IRAF ⁴⁰ environment. Wavelength calibration was performed with a precision of 0.65 Å using He i and Ne i arc lines observed on the same night. The instrumental response was corrected using observations of the spectrophotometric standard star G191–B2B. Because the primary target and the standard star were observed close in time and at a similar airmass, we corrected for telluric lines absorption by dividing the target spectrum by the spectrum of the standard normalized to the continuum. The final low-resolution spectrum of TOI-1749 depicted in Figure 2 has a spectral resolution of 16 Å (R ≈ 450 at 7100 Å).

Zoom In Zoom Out Reset image size

As part of our standard process for validating transiting exoplanets to assess the possible contamination of bound or unbound companions on the derived planetary radii (Ciardi et al. 2015), TOI 1749 was observed with infrared high-resolution adaptive optics (AO) imaging at the Keck Observatory. The observations were made with the NIRC2 instrument on Keck II behind the natural guide star AO system (Wizinowich et al. 2000) on 2020 May 28 UT in the standard 3-point dither pattern that is used with NIRC2 to avoid the left lower quadrant of the detector, which is typically noisier than the other three quadrants. The dither pattern step size was 3'' and was repeated four times.

The observations were made in the narrowband Br − γ filter (λ_o = 2.1686; Δλ = 0.0326 μ m) with an integration time of 1 s with one coadd per frame for a total of 3 s on target. The camera was in the narrow-angle mode with a full FOV of ∼ 10'' and a pixel scale of 0099442 per pixel. The AO data were processed and analyzed with a custom set of IDL tools. The science frames were flat fielded and sky subtracted. The flat fields were generated from a median average of dark-subtracted flats taken on-sky. The flats were normalized such that the median value of the flats is unity. The sky frames were generated from the median average of the 15 dithered science frames; each science frame was then sky subtracted and flat fielded. The reduced science frames were combined into a single combined image using an intra-pixel interpolation that conserves flux, shifts the individual dithered frames by the appropriate fractional pixels, and median-coadds the frames. The final resolution of the combined dithers was 0050, which was determined from the FWHM of the point-spread function (PSF).

The Keck AO observations revealed no additional stellar companions within a resolution of ∼0051 FWHM (Figure 3). The sensitivities of the final combined AO image were determined by injecting simulated sources azimuthally around the primary target every 20° at separations of integer multiples of the central source's FWHM (Furlan et al. 2017). The brightness of each injected source was scaled until standard aperture photometry detected it with 5σ significance. The resulting brightness of the injected sources relative to the target set the contrast limits at that injection location. The final 5σ limit at each separation was determined from the average of all of the determined limits at that separation and the uncertainty on the limit was set by the rms dispersion of the azimuthal slices at a given radial distance. The sensitivity curve is shown in Figure 3 along with a "zoomed" inset image centered on the primary target showing no other companion stars.

Zoom In Zoom Out Reset image size

3.1.1. Spectroscopic Properties and Kinematics

3.1.2. Photometric Properties

The spectral type and age of TOI-1749 estimated from the ALFOSC spectrum are also supported by the location of TOI-1749 in the Gaia color–magnitude diagram (CMD) displayed in Figure 4; this star has colors and absolute magnitudes compatible with its being a main-sequence M0 star. All stellar sequences shown in Figure 4 were built using Gaia photometry and parallaxes (see Luhman 2018 for the young isochrones and Cifuentes et al. 2020 for the main-sequence track). Note that TOI-1749's absolute G-band magnitude is slightly fainter than that of the main sequence for its color, which is consistent with a slightly sub-solar metallicity.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

We estimated the physical parameters of the host star empirically via the spectral energy distribution (SED) constructed from catalog broadband photometric magnitudes, according to the methodology described in Stassun & Torres (2016) and Stassun et al. (2017). The available photometry includes the near-UV magnitude from GALEX, the GG_BP G_RP magnitudes from Gaia, the ugri magnitudes from SDSS and Pan-STARRS, the JHK_S magnitudes from 2MASS, and the W1–4 magnitudes from WISE, as shown in Figure 5. In this fitting, we imposed a Gaussian prior of [Fe/H] = −0.05 ± 0.20 dex on the metallicity, which represents the metallicity distribution of nearby M dwarfs (Gaidos & Mann 2014), and we also fit for T_eff and $\mathrm{log}g$. We also assumed no interstellar extinction due to the proximity of the star (100 pc). Note that according to the 3D extinction map of Green et al. (2019), ⁴¹ the reddening along the line of sight toward the minimum reliable distance of 274 pc is estimated to be E(g − r) = 0.009 ± 0.009, which justifies the above assumption. From this analysis, we obtained T_eff = 3985 ± 55 K, $\mathrm{log}g=4.70\pm 0.05$, and [Fe/H] = −0.26 ± 0.08 dex.

According to the TESS Input Catalog ver. 8 (TICv8, Stassun et al. 2019), the mass and radius of TOI-1749 are M_s = 0.555 ± 0.020M_⊙ and R_s = 0.561 ± 0.017R_⊙, respectively, which were calculated using empirical mass−M_K (Mann et al. 2019) and radius−M_K (Mann et al. 2015) relations for M dwarfs, respectively, where M_K is the K-band absolute magnitude. Here, we updated these estimates using the revised distance from Gaia DR2 and the [Fe/H] value estimated from the SED fitting, obtaining M_s = 0.58 ± 0.03M_⊙ and R_s = 0.55 ±0.03R_⊙. The derived physical parameters of TOI-1749 are summarized in Table 1.

3.1.3. Stellar Variability

3.2.1. Searching for a Third Planetary Candidate

The candidate transit signals of TOI-1749c and TOI-1749d were first identified from the first six sectors (Sectors 14–19) by the SPOC pipeline, and were later confirmed with additional six sectors (Sectors 20–21 and 23–26). However, no additional threshold crossing event (TCE) was detected by the pipeline, not even from all the 12 sectors. We confirmed the transit signals of TOI-1749c and TOI-1749d in the PDC-SAP light curve using the transit least squares (TLS) method (Hippke & Heller 2019), yet failed again to find an additional planetary signal in the residual light curve from which the transit signals of TOI-1749c and d were subtracted.

To further search for additional planetary candidates, we removed systematic trends in the residual PDC-SAP light curve that are apparent in several sectors. To model the systematic trends, we applied a Gaussian process (GP) implemented in celerite (Foreman-Mackey et al. 2017) with a stochastically driven, damped simple harmonic oscillator (SHO) as a kernel function. The power spectral density of SHO is written as

$$\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}/{{ \mathcal Q }}^{2}},\end{eqnarray} \tag{ 1 }$$

where ω₀ is the frequency of undamped oscillation, S₀ is a scale factor to the amplitude of the kernel function, and ${ \mathcal Q }$ is a quality factor (Anderson et al. 1990; Foreman-Mackey et al. 2017). This model can capture damping oscillatory behaviors with a characteristic frequency of ${\omega }_{0}\sqrt{1-1/4{{ \mathcal Q }}^{2}}$ if ${ \mathcal Q }\gt 1/2$ (Foreman-Mackey et al. 2017). Because the amplitudes and timescales of the systematics can be different from sector to sector, we independently modeled S₀ and ω₀ for each sector, while fixing ${ \mathcal Q }$ at unity for all sectors to avoid over fitting.

After removing the systematic trends using the Gaussian process model with the derived median values, we performed the TLS analysis again, finding an additional transit signal with an orbital period of 2.389 days, as shown in the left panel of Figure 6. This signal has a signal detection efficiency (SDE) of 18.2, which corresponds to a false alarm probability of 8 × 10⁻⁵.

Zoom In Zoom Out Reset image size

The same transit signal was also detected with an independent pipeline of Open Transit Search (OpenTS ⁴² ) with a signal-to-noise (S/N) ratio of the transit signal of 4.8.

Motivated by the detection of this marginal signal of a third planet, we revisited the transit signal search by the SPOCs Transiting Planet Search (TPS) module. This module iteratively searches a light curve for transiting planet signatures until no TCE can be returned (Twicken et al. 2018; Li et al. 2019), recording the strongest signals even when none meet the detection criteria. Based on all available data in sectors 14–26, we found that the strongest signal for a third transiting planet was consistent with the signal detected by the TLS and OpenTS analyses, with a multiple event detection statistic (MES) of 7.02. This maximum MES value was just over 1% below the transit detection threshold for the TESS transit search of 7.1, and was not initially sent to the SPOC Data Validation (DV) process. Note that we estimate the SDE of the signal detected by TPS as 18.2 (as shown in the right panel of Figure 6), which is identical to the SDE value of the signal detected by the GP+TLS analysis. ⁴³ This coincidence can be understood as follows. Both TPS and GP+TLS explicitly account for the correlation structure of the observation noise (residual systematic errors plus stellar variability). GP+TLS explicitly models the correlation structure of the observation noise by construction. That said, TPS models the Power Spectral Density (PSD) of the observation noise process using a nonparametric model via an adaptive, wavelet-based matched filter (Jenkins 2002). TPS thus performs a joint noise-characterization/signal detection task that explicitly accounts for time-varying, non-white observation noise in the formulated detector. Since the PSD is defined as the Fourier Transform of the correlation function of a wide-sense stationary process (Kay 1999), the approach implemented in TPS can be viewed as the dual of the use of a GP to model the light curve prior to running a search with TLS on the residual light curve.

We then ran the DV module in a supplemental mode for the candidate planets in TOI-1749 including the third planetary candidate, TOI-1749b. Based on a limb-darkened transiting planet model fit (Li et al. 2019) in the DV module, the signal-to-noise (S/N) ratio of the transit signal of TOI-1749b was calculated as 7.5. Furthermore, the transit signal was fully consistent with that of a transiting planet for all supplemental DV diagnostics: odd/even transit depth consistency test, difference image centroid offset analysis, and optical ghost diagnostic test (Twicken et al. 2018).

Given that this candidate signal was detected by three independent analyses and passed the above DV tests, we consider TOI-1749b as a robust planetary candidate and include it in the subsequent analyses. It should be noted, however, that although the difference image centroid offsets indicated that the location of the transit source was consistent with that of the target star, they were unable to eliminate all nearby stars as the source of the transit signal. We statistically validate that the target star is the source of the transit signal in Section 3.4. This is also supported by a marginal detection of a transit signal of TOI-1749b on target from the MuSCAT2 observation as described in Section 3.3.2.

3.2.2. Transit Model Fit

To estimate the transit parameters of the three planetary candidates, we first modeled the PDC-SAP light curve with transit models assuming that all three planetary candidates have constant orbital periods. We modeled the transit light curves with a Mandel & Agol model implemented by PyTransit (Parviainen 2015) with the following parameters: scaled semimajor axis a/R_s , impact parameter b, planet-to-star radius ratio R_p /R_s , eccentricity e, longitude of periastron ϖ, orbital period P, reference transit time t₀, and two coefficients u₁ and u₂ for the stellar limb-darkening effect for which we assume a quadratic limb-darkening law. We assumed a circular orbit for all planets by fixing e and ϖ to zero. In addition, we placed informative Gaussian priors on a/R_s , u₁, and u₂ as follows. For the mean value and standard deviation of a/R_s , we utilized the stellar density estimated from M_s and R_s derived in Section 3.1, ρ_s = 4.91 ± 0.84 g cm⁻³, which was converted to a/R_s for the circular orbits of TOI-1749b, TOI-1749c, and TOI-1749d as a/R_s = 11.40 ± 0.65, 17.35 ± 0.99, and 27.7 ± 1.6, respectively. For u₁ and u₂, we adopted u₁ = 0.320 ± 0.012 and u₂ = 0.255 ± 0.017, which were calculated by LDTk (Parviainen & Aigrain 2015) based on PHOENIX stellar models (Husser et al. 2013) for the range of stellar parameters estimated in Section 3.1. Note that we enlarged the uncertainties of u₁ and u₂ provided by LDTk by a factor of 3 considering the systematic differences between stellar models and calculation methods for a given set of stellar parameters found in the tabulated values of Claret (2017). We assumed uniform priors for the other parameters. Simultaneously with the transit models, we also modeled systematic trends in the PDC-SAP light curve using the Gaussian process model in the same way as in Section 3.2.1.

To estimate the posterior distributions of the total of 53 free parameters, we ran a Markov Chain Monte Carlo (MCMC) analysis using emcee (Foreman-Mackey et al. 2013). From an MCMC run with 106 walkers and 2 × 10⁴ steps after convergence, we derived the median values and 1σ uncertainties of the parameters as listed in Table 4.

Next, to investigate the effect of possible TTVs, we fit the same light curve with the same transit+systematic models, but set aside the assumption of linear ephemerides and, instead, let individual mid-transit times (T_c ) be free for TOI-1749c and TOI-1749d. Note that because the S/N ratios of individual transits of TOI-1749b are too low to measure T_c , we held the assumption of a linear ephemeris for this planetary candidate. Given the number of transits that the TESS data cover of 58 and 28 for TOI-1749c and TOI-1749d, respectively, the number of free parameters in this model is 135. We ran an MCMC analysis in the same way as before, but this time with 270 walkers and 4 × 10⁴ steps. The measured individual transit times are listed in Table 5, and the other transit parameters are reported in Table 4. The best-fit light curve model is shown in red in Figure 7, and phase-folded light curves are shown in Figure 8. The individual-T_c model improves the best χ² value by 438.0 with 82 additional free parameters over the linear ephemeris model, for a total of 188,918 data points. Among these data points, most of the χ² improvement comes from only the data around individual transits, which corresponds to ∼20,000 data points. Therefore, the difference of Bayesian information criterion (BIC $\equiv {\chi }^{2}+k\mathrm{ln}N$, where k is the number of free parameters and N is the number of data points; Schwarz 1978) between the two models is estimated to be ∼−370. This indicates that the individual-T_c model does not better describe the data over the linear ephemeris model; i.e., we do not detect significant TTV signals in the TESS data from this analysis. Note that the same conclusion was obtained by fitting a straight line to the set of derived T_c values and transit epochs, where the uncertainty of T_c was approximated to be a Gaussian; the fit gives a reduced χ² of 0.69 and 1.3 for TOI-1749c and d, respectively.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Table 5. Mid-times of Individual Transits

Epoch	Tc	1σ error	1σ error	Instrument	Planet
	(BJDTDB )	(upper)	(lower)
0	2458685.9405	0.0398	0.0128	TESS	c
1	2458690.4394	0.0359	0.0735	TESS	c
2	2458694.9077	0.0089	0.0090	TESS	c
3	2458699.4013	0.0112	0.0117	TESS	c
4	2458703.9008	0.0087	0.0263	TESS	c
⋯

Note. This table is available in its entirety in machine-readable form. A portion is shown here for guidance regarding its form and content.

Only a portion of this table is shown here to demonstrate its form and content. A machine-readable version of the full table is available.

Download table as:  DataTypeset image

Table 4. Results of Transit Model Fitting

Parameter	Unit	TESS	TESS	Ground-based
		(linear ephemeris fit)	(individual Tc fit)
TOI-1749b
P	days	$2.388821{\ }_{-0.000073}^{+0.000057}$	$2.388843{\ }_{-0.000076}^{+0.000044}$	2.388843 (fixed)
t0	BJDTDB − 2458680	$4.3613{\ }_{-0.0049}^{+0.0070}$	$4.3597{\ }_{-0.0040}^{+0.0066}$	⋯
a/Rs		$11.53{\ }_{-0.65}^{+0.61}$	11.60 ± 0.55	11.60 (fixed)
b		$0.76{\ }_{-0.15}^{+0.08}$	$0.74{\ }_{-0.16}^{+0.10}$	0.74 (fixed)
Rp /Rs (TESS)		0.0239 ± 0.0025	$0.0241{\ }_{-0.0022}^{+0.0025}$	⋯
Rp /Rs (ground)		⋯	⋯	$0.0244{\ }_{-0.0071}^{+0.0056}$
TOI-1749c
P	days	$4.489010{\ }_{-0.000061}^{+0.000040}$	4.489010 (fixed)	4.489010 (fixed)
t0	BJDTDB − 2458680	$5.9386{\ }_{-0.0019}^{+0.0028}$	⋯	−–
a/Rs		17.7 ± 0.8	17.2 ± 0.8	17.3 ± 0.6
b		$0.29{\ }_{-0.18}^{+0.16}$	$0.21{\ }_{-0.13}^{+0.15}$	0.22 ± 0.12
Rp /Rs (TESS)		0.0344 ± 0.0011	0.0337 ± 0.0012	−–
Rp /Rs (g)		⋯	⋯	$0.0336{\ }_{-0.0045}^{+0.0039}$
Rp /Rs (r)		⋯	⋯	0.0331 ± 0.0020
Rp /Rs (i)		⋯	⋯	0.0385 ± 0.0018
Rp /Rs (zs )		⋯	⋯	0.0336 ± 0.0019
TOI-1749d
P	days	9.04455 ± 0.00012	9.04455 (fixed)	9.04455 (fixed)
t0	BJDTDB − 2458680	8.7805 ± 0.0028	⋯	−–
a/Rs		28.0 ± 1.5	28.2 ± 1.3	27.6 ± 1.2
b		$0.719{\ }_{-0.049}^{+0.042}$	$0.733{\ }_{-0.058}^{+0.049}$	$0.718{\ }_{-0.034}^{+0.031}$
Rp /Rs (TESS)		$0.0387{\ }_{-0.0018}^{+0.0016}$	$0.0399{\ }_{-0.0022}^{+0.0017}$	−–
Rp /Rs (g)		⋯	⋯	$0.0416{\ }_{-0.0038}^{+0.0030}$
Rp /Rs (r)		⋯	⋯	$0.0403{\ }_{-0.0032}^{+0.0030}$
Rp /Rs (i)		⋯	⋯	$0.0435{\ }_{-0.0022}^{+0.0019}$
Rp /Rs (zs )		⋯	⋯	$0.0447{\ }_{-0.0032}^{+0.0030}$

Download table as:  ASCIITypeset image

3.3.1. TOI-1749c and TOI-1749d

We modeled the light curves of TOI-1749c and TOI-1749d obtained from the ground-based follow-up observations as follows. All of the photometric data sets of each planetary candidate that show significant transit signals (check-marked in Table 3) were simultaneously fit with transit+systematic models. For the transit model, we used the same one described in Section 3.2.2, but allowed R_p /R_s to be free for each band to check for a chromaticity dependence in this parameter, which can be a sign of flux contamination in the photometric aperture. We let T_c be free for each transit epoch, while treating b and a/R_s as common parameters for all light curves. The limb-darkening parameters of u₁ and u₂ were fixed for each band to the theoretical values given by LDTk, specifically, (u₁, u₂) = (0.656, 0.129), (0.541, 0.183), (0.374, 0.248), and (0.302, 0.253) for the g, r, i, and z_s bands, respectively.

We modeled the systematics in each light curve by a combination of a linear function of ΔX and ΔY, which are the stellar displacements on the detector in the X and Y directions, respectively, and a Gaussian process model as a function of time with an approximated Matérn 3/2 kernel implemented in celerite. The former function takes account of systematics originating from stellar movements on the detector, which are typically within a few pixels during each transit observation. When stars move on the detector, they are subject to photometric systematics caused by imperfections in the flat-field corrections. The Gaussian process models other nonoscillating, time-correlated noise, which presumably mostly originates in effects of observing through the Earth's atmosphere. The approximated Matérn 3/2 kernel is written by

$$\begin{eqnarray}\begin{array}{rcl}k(\tau ) & = & {\sigma }^{2}\left[(1+\displaystyle \frac{1}{\epsilon }){e}^{-(1-\epsilon )\sqrt{3}\tau /\rho }\right.\\ & & \left.\times \,(1-\displaystyle \frac{1}{\epsilon }){e}^{-(1+\epsilon )\sqrt{3}\tau /\rho }\right],\end{array}\end{eqnarray} \tag{ 2 }$$

where τ is the distance of two data points in time, and σ, ρ, and are coefficients. We fixed , which is a quality factor for approximation, to a default value of the code of 0.01. Each light curve was modeled by the GP model with a mean function of

$$\begin{eqnarray}&&\mu ={{ \mathcal F }}_{\mathrm{transit}}\times ({c}_{0}+{c}_{x}{\rm{\Delta }}X+{c}_{y}{\rm{\Delta }}Y),\end{eqnarray} \tag{ 3 }$$

where ${{ \mathcal F }}_{\mathrm{transit}}$ is the transit model and c₀, c₁, and c₂ are the coefficients for the linear systematic model. Among the parameters for the systematic model, we forced ρ to be a common parameter for each transit assuming that the timescale of the time-correlated noise in each night is shared among the different bands, while letting σ, c_x , and c_y be free for individual light curves. Note that c₀ was obtained for a given set of parameters by taking a median value of $F/{{ \mathcal F }}_{\mathrm{transit}}-({c}_{x}{\rm{\Delta }}X+{c}_{y}{\rm{\Delta }}Y)$, where F is the observed flux. We also modeled a white jitter noise for each light curve, σ_jitter, in the form of ${\sigma }_{\mathrm{flux}}=\sqrt{{\sigma }_{\mathrm{calc}}^{2}+{\sigma }_{\mathrm{jitter}}^{2}}$, where σ_flux is the uncertainty of the flux of each data point and σ_calc is the theoretical uncertainty calculated by the photometric pipeline.

With this model, we ran MCMC for a total of 23 light curves (6 transit epochs) and 19 light curves (7 transit epochs) for TOI-1749c and TOI-1749d, with numbers of free parameters of 110 and 96, respectively. For a/R_s and b, we imposed informative Gaussian priors using the results of the MCMC analysis for the TESS light curve with the individual T_c model; that is, (a/R_s , b) = (17.2 ± 0.8, 0.21 ± 0.14) and (28.2 ± 0.13, 0.733 ± 0.053) for TOI-1749c and TOI-1749d, respectively. We applied uniform priors for the other parameters, with natural logarithmic form for ρ, σ, and σ_jitter. Using emcee, we calculated posterior probability distributions of the parameters from 4 × 10⁴ MCMC steps with 220 and 192 walkers for TOI-1749c and TOI-1749d, respectively. The derived median and 1σ confidence intervals of the mid-transit times and the other transit parameters are appended to Tables 4 and 5, respectively. We show 100 light-curve models randomly selected from the posterior distributions along with the individual ground-based light curves in Figures 9 and 10, and also show the best-fit transit model along with phase-folded, systematic-corrected light curves for each planet and each band in Figure 11. We find that the R_p /R_s values of respective planets measured in five different bands (one from TESS and four from the ground) are consistent with each other within ∼2σ as shown in Figure 12, without any particular trend as a function of wavelength. We thus find no evidence of flux contamination within the photometric aperture of ground-based observations (26–52). We note that the atmospheric scale heights of TOI-1749c and TOI-1749d are expected to be ∼5 × 10⁻⁴ and ∼4 × 10⁻⁴ in units of R_s , assuming hydrogen-rich atmospheres with planetary masses of 5 and 7 M_⊕, respectively (the masses predicted using an empirical relation of Chen & Kipping 2017). Thus, the effect of atmospheric opacity on the measured R_p /R_s , which can be ∼5 times the scale height, is comparable to or smaller than the measurement uncertainties, unless the planetary masses are unusually small.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

3.3.2. TOI-1749b

The light curves of TOI-1749b obtained with MuSCAT2 show no apparent transit signal, as shown in the top panel of Figure 13 (a). To see if there is the expected transit signal of this candidate, which has a depth of ∼0.058% (∼580 ppm) and an expected mid-transit time of ${T}_{c}={2459133.462}_{-0.015}^{+0.011}$, we fit a transit+systematics model to the r-, i-, and z_s -band light curves in the same way as in Section 3.3.1, but fixing b and a/R_s to the values obtained from the TESS light curves, i.e., 0.74 and 11.60, respectively. The parameters of R_p /R_s and T_c were shared between the three bands. To avoid sampling unrealistic parameter space, we limited the allowed range of T_c to ±3σ from the above expected value. We note that the expected TTV amplitude for the transits of this candidate is at most a few minutes, as shown in Section 3.5, which is negligible compared to the timing uncertainty from the linear ephemeris. As a result, we find a marginal transit signal at ${T}_{c}=2459133.4612{\ }_{-0.0027}^{+0.0079}$ at a significance of ∼3σ. The radius ratio is measured to be ${R}_{p}/{R}_{s}=0.0244{\ }_{-0.0071}^{+0.0056}$. The posterior transit model and systematic-corrected light curves are shown in the bottom panel of Figure 13 (a), and the posterior distributions of T_c and R_p /R_s are shown in Figure 13 (b). Although the significance of the signal is marginal, the measured T_c and R_p /R_s are consistent with the predicted values from TESS, which are shown in blue and red lines in Figure 13 (b), supporting the detection of this near-threshold planetary candidate in the TESS data (Section 3.2.1).

Zoom In Zoom Out Reset image size

We have confirmed the transit signals of all three planetary candidates identified from the TESS data by the ground-based photometric observations (although marginal for TOI-1749b) on target within apertures of 26–52 radius. The observed transit depths from the ground are all consistent with those from TESS, with no apparent chromaticity for TOI-1749c and TOI-1749d, as shown in Figure 12. We have also confirmed from the Keck AO observation that there is no companion or background star down to a magnitude difference of 5 within 015 from TOI-1749 (Figure 3). Therefore the detected transit signals most likely originate from TOI-1749 itself. However, there is still some chance that the each signal comes from an eclipsing binary that contaminates the 015 aperture.

We computed false positive probabilities (FPPs) for each planetary candidate using the Python package vespa (Morton 2015). The vespa package was developed for the statistical validation of planets in bulk, e.g., from the Kepler mission (Morton et al. 2016), which were too numerous or faint to permit detailed follow-up observations sufficient to robustly ascertain their dispositions. The vespa package employs a robust statistical framework to compare the likelihood of the planetary scenario to likelihoods of several astrophysical false positive scenarios involving eclipsing binaries, relying on simulated eclipsing populations based on the TRILEGAL Galaxy model (Girardi et al. 2005). As input data to vespa, we used the phase-folded TESS light curve for each planet candidate (width ∼ 4 × T₁₄), the Keck K-band contrast curve (see Figure 3), the 3σ upper limit on the secondary eclipse depth (0.34, 0.28, and 0.46 ppt for TOI-1749b, c, and d, respectively), and a radius of exclusion for the transit signal of 26, which corresponds to the minimum aperture radius for ground-based transit observations (see Section 2.2). The FPPs from vespa for planets TOI-1749b, TOI-1749c, and TOI-1749d are 0.092, 1.6 × 10⁻¹², and 6.8 × 10⁻⁵, respectively. Moreover, since vespa does not account for multiplicity, these FPPs are overestimated: Lissauer et al. (2012) demonstrated that a candidate in a system with one (two) or more additional transiting planet candidates is 25 (50) times more likely to be a planet based on multiplicity alone. Application of the 25 (50) multi-boost factor to TOI-1749b yields a probability of 99.6% (99.8%) that the candidate is a bona fide planet. Given these considerations, we find that the FPPs of all three candidates fall below the standard validation threshold of 1%.

Because the three planets in this system orbit in a compact region (≤0.07 au), and also the orbits of TOI-1749c and TOI-1749d are close to the 2:1 commensurability, they may exhibit measurable TTVs due to gravitational interactions between the planets. If one assumes that TOI-1749c and TOI-1749d each have a mass of 7 M_⊕ (a mass predicted for a 2.5 R_⊕ planet by Chen & Kipping 2017) and zero free eccentricity, and they are not locked in a MMR, then one can predict that the amplitudes of TTVs are 7.4 and 4.1 minutes for TOI-1749c and d, respectively, with a timescale of the super period of ∼610 days (Lithwick et al. 2012).

The TESS observation of this target spanned almost one year, which is more than half of the super period. However, as discussed in Section 3.2.2, no evidence of TTVs was found in the TESS data. However, we calculate the orbital periods of TOI-1749c and TOI-1749d based on the individual T_c values from the ground-based observations assuming linear ephemerides to be 4.48934 ± 0.00011 days and 9.04535 ± 0.00032 days, respectively. These values are both slightly inconsistent with those derived from the TESS data (see Table 5) by ∼ 3σ, which might be a sign of TTV signals. Nevertheless, TTV signals are not apparent even in the ground-based data when we compare them with the linear ephemerides derived using both TESS and ground-based observations (see Figure 14), which are

$$\begin{eqnarray*}\begin{array}{rcl}{T}_{c,{\rm{c}}}({\mathrm{BJD}}_{\mathrm{TDB}}) & = & 2458685.9350\ (26)+4.489093\ (31)\times {E}_{{\rm{c}}},\\ {T}_{c,{\rm{d}}}({\mathrm{BJD}}_{\mathrm{TDB}}) & = & 2458688.7790\ (28)+9.044669\ (76)\times {E}_{{\rm{d}}},\end{array}\end{eqnarray*}$$

where E_c and E_d are transit epochs of TOI-1749c and d, respectively.

Zoom In Zoom Out Reset image size

Although there are no apparent TTV signals in the measured transit times, TTV signals could be enhanced by the so-called photodynamical modeling (Carter et al. 2011, 2012; Almenara et al. 2015), which directly models the light curves with both transit light-curve and TTV models simultaneously. This method has several advantages. First, this method does not require measurements of individual T_c 's from light curves but instead models the light curves with only the physical parameters of the planetary system. This effectively reduces the total number of free parameters and thus improves the determination of the physical parameters. Next, with this method the eccentricity (e) and longitude of periastron (ϖ) can be shared between the TTV and transit light-curve models. Because these parameters are difficult to constrain from the transit light-curve model alone, they are often simply assumed to be zero. However, because a/R_s can be converted from e, ϖ, P, and stellar density ρ_s , the constraints on e and ϖ from the TTV model allow one to use ρ_s as a common parameter within multiple planets in the same system instead of a/R_s for each planet in the transit light-curve model; this also helps to reduce the number of free parameters. Finally, this technique allows one to properly estimate the uncertainties of the parameters even if individual T_c 's have asymmetric (non-Gaussian) posterior distributions due to, e.g., partial transit coverage or systematics in the light curve.

To constrain the mass and eccentricity of the planets and obtain better estimates of other transit parameters, we performed a joint analysis of the TESS and ground-based (those check-marked in Table 3) light curves by photodynamical modeling. In this analysis, we used the detrended TESS light curves, in which the systematics in the original light curves are corrected by the systematic model derived from the individual T_c fit in Section 3.2.2, while we used the un-detrended light curves (same as the ones used in Section 3.3.1) for the ground-based data, in which the systematic and photodynamical models could be non-negligibly correlated due to the similar timescales between the systematics and transits. We calculated TTV models using TTVFast (Deck et al. 2014), where we used as the parameters the mean orbital period P, planet-to-star mass ratio M_p /M_s , eccentricity and argument of periastron in the forms of $\sqrt{e}\sin \varpi$ and $\sqrt{e}\cos \varpi$, and orbital phase $\varpi +{ \mathcal M }$, where ${ \mathcal M }$ is the mean anomaly at a reference time (BJD = 2458683.354296). We assumed coplanar orbits with an orbital inclination of 90° and argument of ascending node of 180° for all planets, and used a time step of 0.05 days for the dynamical calculation. For the transit model, we applied the same model as described in Section 3.2.2 and Section 3.3.1, but dropped the assumption of circular orbits and allowed nonzero values for e and ϖ, which were shared with the TTV model. The T_c of individual transits were calculated by the TTV model. We assumed an achromatic transit depth and used a single R_p /R_s for all bands for each planet. For the systematic models of ground-based data, we fixed c_x and c_y to the values determined in Section 3.3.1 to reduce the number of free parameters, while letting ρ and σ be free in the same way as in Section 3.3.1. We also fixed the white-noise jitter values to those determined in Sections 3.2.2 and 3.3.1 for TESS and ground-based data, respectively. The total number of free parameters for three planets was 23.

We then estimated the posterior distributions of the parameters in an MCMC analysis using emcee. Because the likelihood function with respect to the parameters of a TTV model tends to be multimodal, it is not efficient to start MCMC chains from one single location in the parameter space. Alternatively, we ran MCMC with 400 walkers in which M_p , e, and ϖ were initialized with uniformly random values in the ranges of [1, 20] M_⊕, [0.01, 0.2], and [−π, π], respectively, while the other parameters were initialized around the median values of the posteriors calculated in Section 3.2.2. In addition, we sampled MCMC chains with the differential evolution algorithm, which can efficiently find global maxima of the likelihood function even if it is multimodal. After running MCMC for 4 × 10⁴ steps, we sorted the walkers by likelihood of the last step and discarded the lower half of the walkers, because many of them were stuck at local minima of the likelihood function. We then continued MCMC runs with the remaining walkers for a sufficient number of steps until the chains converged, and sampled the last 8 × 10⁴ steps, which span intervals roughly 20 times the auto-correlation lengths, to construct the posterior probability distributions of the parameters.

We summarize the median values and 1σ confidence intervals of the parameters in Table 6. In Table 7 we also show the derived parameters of e, ϖ, R_p , and M_p , where we apply R_s and M_s derived in Section 3.1 to calculate R_p , and M_p , respectively. In Figure 14, we display randomly selected 100 posterior TTV models for each planet along with the individually measured transit times in the previous sections for comparison.

Table 6. Results of Photodynamical Analysis

Parameter	Unit	Host star	TOI-1749b	TOI-1749c	TOI-1749d
ρs	g cm−3	$5.02{\ }_{-0.46}^{+0.59}$	⋯	⋯	⋯
u1 (TESS)		$0.329{\ }_{-0.013}^{+0.015}$	⋯	⋯	⋯
u2 (TESS)		0.249 ± 0.019	⋯	⋯	⋯
P	days	⋯	$2.38839{\ }_{-0.00066}^{+0.00031}$	$4.4929{\ }_{-0.0027}^{+0.0038}$	$9.0497{\ }_{-0.0032}^{+0.0049}$
b		⋯	$0.75{\ }_{-0.20}^{+0.11}$	$0.36{\ }_{-0.31}^{+0.09}$	$0.717{\ }_{-0.032}^{+0.028}$
Rp /Rs		⋯	$0.0232{\ }_{-0.0029}^{+0.0032}$	0.0353 ± 0.0008	0.0421 ± 0.0010
Mp /Ms	10−5	⋯	$10{\ }_{-7}^{+11}$	$1.1{\ }_{-0.8}^{+2.9}$	$2.2{\ }_{-1.8}^{+3.2}$
$\sqrt{e}\cos \varpi$		⋯	−0.05 ± 0.12	$0.025{\,}_{-0.043}^{+0.059}$	−0.02 ± 0.08
$\sqrt{e}\sin \varpi$		⋯	0.08 ± 0.13	$0.03{\,}_{-0.06}^{+0.09}$	−0.03 ± 0.11
$\varpi +{ \mathcal M }$	radians	⋯	$2.051{\,}_{-0.069}^{+0.035}$	1.102 ± 0.010	$0.943{\,}_{-0.022}^{+0.019}$

Download table as:  ASCIITypeset image

Table 7. Derived Parameters from the Results of Photodynamical Analysis

Parameter	Unit	TOI-1749b	TOI-1749c	TOI-1749d
Planetary radius (Rp )	R⊕	$1.39{\ }_{-0.19}^{+0.21}$	2.12 ± 0.12	2.52 ± 0.15
Planetary mass (Mp )	M⊕	$19{\ }_{-13}^{+22}$ (57)	$2.1{\ }_{-1.6}^{+5.7}$ (14)	$4.3{\ }_{-3.5}^{+6.2}$ (15)
Eccentricity (e)		$0.02{\ }_{-0.02}^{+0.06}$ (0.14)	$0.007{\ }_{-0.005}^{+0.007}$ (0.019)	$0.015{\ }_{-0.011}^{+0.017}$ (0.062)
Longitude of periastron (ϖ)	radians	$2.3{\,}_{-1.0}^{+1.7}$	$1.8{\,}_{-1.0}^{+3.2}$	$3.8{\,}_{-2.2}^{+1.3}$
Orbital inclination (i)	degrees	$86.4{\ }_{-0.6}^{+0.9}$	$88.8{\ }_{-0.3}^{+1.0}$	$88.53{\ }_{-0.09}^{+0.11}$
Semimajor axis (a)	au	0.0291 ± 0.0005	0.0443 ± 0.0008	0.0707 ± 0.0012
Equilibrium temperature (Teq) a	K	831 ± 18	673 ± 15	533 ± 12

Note. The value in parentheses represents the 95% confidence upper limit.

^a Calculated for a planet with zero albedo and a constant surface temperature.

Download table as:  ASCIITypeset image

From this analysis, we find that a nonzero mass is preferred for all planets at the ∼1σ level (${19}_{-13}^{+22}$ M_⊕, ${2.1}_{-1.6}^{+5.7}$ M_⊕, and ${4.3}_{-3.5}^{+6.2}$ M_⊕ for TOI-1749b, c, and d, respectively), which is indicative of TTV signals in the data. To further check for the significance of the TTV signals, we also fit a constant period model to the same data sets assuming all planets have circular orbits, and compared the minimum χ² values with that of the photodynamical model. As a result, we find that the χ² value of the photodynamical model improves over the null-TTV model by 762. Given the total numbers of data points relevant to individual transits of ∼45,000 (∼25,000 from the ground-based data) and the number of additional free parameters of 9, this χ² improvement corresponds to a BIC improvement of ∼666, which indicates that the TTV model better describes the observed data over the null-TTV one. Note that this χ² improvement entirely comes from the ground-based data. The portion of χ² from the TESS data increased by 108, which may indicate the presence of uncorrected systematics in the TESS data, uncorrected systematics in the ground-based data that generates a bias in the TTV model in the region of the TESS data, or an additional body in the system. Further TTV observations are required to investigate these possibilities.

While the planetary masses are not detected at high significance, we place strong upper limits on the mass of 57 M_⊕, 14 M_⊕, and 15 M_⊕ (95% confidence level) for TOI-1749b, c, and d, respectively, confirming that all three planets have a mass well within the planetary range (<13 M_Jup = 4132 M_⊕).

We carried out a set of dynamical simulations to study the long-term stability of the system. We took the parameters reported in Tables 6 and 7 and drew 60,000 samples from the parameters posteriors as initial parameters for the dynamical simulation. Each parameter set was integrated for 10⁹ orbits of the inner planet orbital period using the Stability of Planetary Orbital Configurations Klassifier (SPOCK; Tamayo et al. 2020). The system was found to be dynamically stable for the whole parameter posterior space with a median value of the spock stability probability of 0.84.

Since the photodynamical analysis in Section 3.5 gave only upper limits (95% confidence level) on the masses for TOI-1749b, c, and d, we employed a set of dynamical simulations as described above combined with an MCMC sampling using emcee to explore the overall stability of the system. We sampled the SPOCK stability probability using 50 walkers, running the MCMC over 15,000 iterations of which we used the last 10⁷ and used the first 5 × 10⁶ as burn-in. We allowed for eccentricities up to 0.5 and planet masses up to 100 M_⊕. The posteriors for these two parameters are show in Figure 15 where we report the 95% confidence level as the upper limits. Although the eccentricities have upper limits of 0.19, 0.17, and 0.25 for TOI-1749b, c, and d, respectively, the planet masses for the innermost planets are not constrained at all and only marginally constrained for the outer planet. Those simulations show that the adopted solution (Table 6) lies within the stable regions.

Zoom In Zoom Out Reset image size

In the period−radius plane, the innermost planet (TOI-1749b) and the outer two planets (TOI-1749c and d) are located below and above, respectively, the proposed location of the radius valley for planets around M dwarfs (Van Eylen et al. 2021) as shown in Figure 16. Several studies have shown that this radius valley can be a consequence of photoevaporation and/or core-powered loss of hydrogen envelopes on top of rocky cores; the hydrogen envelopes of planets smaller than a certain size for a given orbital period are selectively dissipated (e.g., Owen & Wu 2017; Ginzburg et al. 2018). Considering the fact that the three planets have similar orbits in the same system, they could have similar initial compositions. If so, then their current radii can be explained by a scenario wherein they all initially consisted of a rocky core surrounded by a thin hydrogen envelope, and then only the innermost planet had lost its envelope due to the photoevaporation and/or core-powered mass-loss mechanisms. For comparison, we also show the locations of TOI-270b, c, and d in Figure 16, which are a benchmark triple orbiting the nearby M dwarf TOI-270 (also known as L231-32; Günther et al. 2019). This trio also has the same trend that the innermost and the outer two planets are located below and above the predicted radius valley, respectively.

Zoom In Zoom Out Reset image size

In the top panel of Figure 17, we show the densities of the posterior probabilities of the mass and radius of TOI-1749b (green), TOI-1749c (blue), and TOI-1749d (green) derived from the photodynamical analysis, along with the known planets with mass and radius both measured taken from the NASA Exoplanet Archive. ⁴⁴ We also show theoretical mass–radius relations from Zeng et al. (2019) for planets with an Earth-like rocky composition (32.5wt% Fe + 67.5wt% MgSiO₃), a rocky+icy composition (50wt% Earth-like rocky core + 50wt% H₂O), an Earth-like rocky core + 0.1wt% hydrogen envelope, and an Earth-like rocky core + 1wt% hydrogen envelope (light-blue solid, dashed, dashed–dotted, and dotted lines, respectively). Although the uncertainties in the masses of the TOI-1749 planets are large enough for the masses and radii to be consistent with a range of compositions, the posterior probabilities of the innermost (TOI-1749b) and the outer two (TOI-1749c and d) planets are the most dense below and above the theoretical line for an Earth-like composition, respectively. This result is consistent with the above scenario, in which TOI-1749b currently consists of a bare rocky core, while TOI-1749c and d still have a hydrogen envelope on top of a rocky core. We note that the median value of the mass of TOI-1749b (19 M_⊕) is unusually large for its size (∼1.4 R_⊕), which, however, could be biased by any systematics in the data, given the low significance of the TTV signals. In particular, there is a well-known mass−eccentricity degeneracy in TTV models, and a smaller mass and larger eccentricity of TOI-1749b could also explain the data with only a slightly lower likelihood. Additional observations are required to test for this possibility.

Zoom In Zoom Out Reset image size

We note that, thanks to its bright host star, the masses and radii of the benchmark trio of the TOI-270 system have been precisely measured through radial velocity and TESS photometry, respectively (Van Eylen et al. 2021), from which the innermost and the outer two planets of this system were confirmed to be consistent with a rocky and rocky+hydrogen envelope compositions, respectively (see Figure 17).

One interesting feature of this system is that the planetary pair of TOI-1749c and d has a period ratio of 2.015, which is very close to the exact 2:1 commensurability with an outward departure of only 0.7%. This feature is shared with two other pairs of planets around M dwarfs, TOI-270 c and d and TOI-175 (also known as L98-59) c and d (Kostov et al. 2019; Cloutier et al. 2019), which have period ratios of 2.011 and 2.019, respectively. These specific period ratios are a priori less likely to be of primordial origin than the consequence of resonant capture by convergent migration in protoplanetary disks, which might be followed by some resonant repulsion effects.

A period ratio just outside of a commensurability has also been commonly observed in planetary pairs around FGK stars, which were mainly discovered by Kepler (Lissauer et al. 2011). However, in the case of the planetary pairs around FGK stars, it is relatively rare that the period ratio of a pair departs from the exact 2:1 commensurability by only 1% or less, as shown in Figure 18. In this figure, the period ratio of planetary pairs with radii not larger than 4 R_⊕ are shown as a function of orbital period of the inner planet for M dwarf and FGK dwarf systems in the top and bottom panels, respectively. This possible difference between M dwarf and FGK dwarf systems implies that the dominant mechanisms that repel planetary pairs around FGK dwarfs from exact commensurabilities would not work in the same way for those around M dwarfs, or at least for the above three pairs.

Zoom In Zoom Out Reset image size

So far, various mechanisms have been proposed to explain the paucity and pileup of planetary pairs at and just outside the exact commensurabilities, respectively, observed in the Kepler-discovered multiplanet systems. These include tidal interactions between the planets and central star (Lithwick & Wu 2012; Batygin & Morbidelli 2013; Delisle & Laskar 2014), planet−disk interactions (Baruteau & Papaloizou 2013), planet-planetesimal interactions (Chatterjee & Ford 2015), and dynamical instabilities of planets in resonant chains (Izidoro et al. 2017).

One of the most frequently discussed mechanisms among them has been tidal interaction. In this scenario, the forced eccentricities of two near-resonant planets raised by the resonance are damped by the tidal effect, which works more efficiently for the inner planet, leading to a repulsion from the resonance. Using Equation (26) of Lithwick & Wu (2012), one can express the fractional distance of the period ratio from the exact 2:1 commensurability that a planetary pair with the initial period ratio of 2 obtains by the tidal effect after the time t as

$$\begin{eqnarray}\begin{array}{rcl}{{\rm{\Delta }}}_{2:1} & \sim & 0.01{\left(\displaystyle \frac{Q{{\prime} }_{1}}{150}\right)}^{-1/3}{\left(\displaystyle \frac{{M}_{p,1}}{10{M}_{\oplus }}\right)}^{1/3}{\left(\displaystyle \frac{{R}_{p,1}}{2{R}_{\oplus }}\right)}^{5/3}\\ & & \times \,{\left(\displaystyle \frac{{M}_{* }}{{M}_{\odot }}\right)}^{-8/3}{\left(\displaystyle \frac{{P}_{1}}{5\mathrm{days}}\right)}^{-13/9}{\left(\displaystyle \frac{t}{5\mathrm{Gyr}}\right)}^{1/3},\end{array}\end{eqnarray} \tag{ 4 }$$

where $Q^{\prime} \equiv 3Q/2{k}_{2}$, Q is tidal quality factor, k₂ is tidal Love number, M_p is planetary mass, R_p is planetary radius, P is orbital period, M_* is stellar mass, and the subscript "1" denotes the inner planet. This equation indicates that the resonant repulsion effect has a strong, negative dependence on the stellar mass, making a factor of ∼6 difference in Δ_2:1 between M (M_* = 0.5M_⊙) and G (M_* = 1M_⊙) dwarf systems for a given set of the other parameters, as shown by dashed lines in Figure 18. Note that the other parameters in Equation (4),i.e., M_p , R_p , $Q^{\prime}$, and t, on which Δ_2:1 has weaker dependence than on M_*, could also be somewhat different between the M and FGK samples. The important point here is that despite the strong and negative dependence of Δ_2:1 on M_*, there do exist planetary pairs closer to the exact 2:1 commensurability around M dwarfs than any of those around FGK dwarfs for a given orbital period of the inner planet. This fact implies that there are some mechanisms that counteract the tidal repulsion effect selectively in the M dwarf systems, or alternative resonant repulsion mechanisms that do not depend on (or have a positive dependence on) the stellar mass play a major role in both types of systems.

Another common feature among the near 2:1 planetary pairs in the TOI-175, TOI-270, and TOI-1749 systems is that they all have an additional inner planet. The period ratios of the inner adjacent pairs (i.e., the period ratio of the middle to the innermost planets) are all between 1.5 and 2, as shown in the top left panel of Figure 19. However, near 2:1 period-ratio planetary pairs in FGK systems tend to have an additional planet at an outer orbit rather than inner one, as shown in the bottom panels of Figure 19. Note that the TOI-178 system is an exception to this trend; i.e., the near 2:1 period ratio pair of TOI-178 c and d has an additional inner planet (TOI-178b) despite the fact that the host star is a K dwarf; however, the host star has an effective temperature of T_eff = 4316 K (Leleu et al. 2021), which is close to the temperature range of M dwarfs, <4000 K.

Zoom In Zoom Out Reset image size

Several previous works have pointed out that planetary pairs in high multiplicity (>2 planets) systems tend to avoid resonant repulsion (Lithwick & Wu 2012; Steffen & Hwang 2015), probably due to extra dynamical forces caused by the additional planets. The different tendency in the orbital architecture of planetary trios between M and FGK dwarf systems, seen in Figure 19, might suggest that the specific orbital configuration of the TOI-175, TOI-270, and TOI-1749 systems has played a role in competing with the resonant repulsion effects. Given that TOI-1749 is the third such example, planetary trios in which the outer pair has a period ratio just beyond 2 might be relatively common around M dwarfs compared to FGK dwarfs, while a rigorous statistical test is beyond the scope of this paper. If confirmed, this fact might reflect the planetary formation and migration processes, which should depend on the stellar mass, or more specifically, the properties of the protoplanetary disks.

Given that TOI-1749 is one of a few M dwarf systems that host multiple transiting planets including a pair close to the 2:1 MMR, further follow-up observations are encouraged.

Continuous monitoring of TTVs is of particular importance in confirming and narrowing down the masses and eccentricities of the planets. As demonstrated in this work, although the TTV amplitudes in the outer two planets are too small for TESS to detect, ground-based 1–2 m class telescopes can measure the times of individual transits with a high enough precision. Multiband cameras like the MuSCATs are particularly efficient for this purpose. As seen in Figure 14, the timing precisions achieved for TOI-1749d from MuSCAT2 (red points) are better than any other data from single-band imagers (blue ones), regardless of whether the transit coverage is full or not. In addition, detecting transits of TOI-1749b would be challenging for single-band imagers mounted on 1–2 m or smaller class telescopes.

Radial velocity (RV) observations can also yield independent confirmation of the planetary masses. The expected RV amplitudes from the individual planets are ∼2–3 m s⁻¹ each if they have a mass of ∼3–4 M_⊕, which is within reach of the current world-best instruments. One drawback is the faintness of the host star (V = 13.9, J = 11.1), which would require instruments on large aperture (∼8 m) telescopes such as Subaru/IRD and Gemini/Maroon-X.

This system may not be a prime target for planetary atmospheric study due to the relative faintness of the host star. Using Equation (1) of Kempton et al. (2018), we estimate the transmission spectroscopy metric (TSM) of TOI-1749c and TOI-1749d both to be 32, where we apply the masses estimated using the empirical relation of ${M}_{p}=1.436{R}_{p}^{1.70}{M}_{\oplus }$ from Chen & Kipping (2017; 5.1 and 6.9 M_⊕ for TOI-1749c and d, respectively). These values are a factor of 3–4 smaller than those for their analogs of TOI-270c and TOI-270d (Chouqar et al. 2020; Van Eylen et al. 2021). However, because TSM is inversely proportional to the planetary mass, they would become good atmospheric targets if the true mass is unusually small, such as 1–2 M_⊕, which is still allowed by the current observations (see Figure 17). We note that TOI-1749 is located very close to the northern continuous viewing zone of the James Webb Space Telescope, allowing flexible scheduling with the promising new space telescope.