Near-resonance in a system of sub-Neptunes from TESS

TL;DR

This paper reports the discovery and validation of a multi-planet system around TOI 125 from TESS data, including sub-Neptunes and candidates, with analysis indicating the system is likely not in resonance despite near 2:1 period ratios.

Contribution

First detailed characterization of a multi-planet system with sub-Neptunes around TOI 125, including validation and dynamical analysis of period ratios.

Findings

Three validated sub-Neptune planets with radii around 2.8 R⊕.

Detection of two low signal-to-noise planet candidates, including a super-Earth.

The system's period ratios are close to, but not in, exact resonance, and are dynamically stable.

Abstract

We report the Transiting Exoplanet Survey Satellite ($TESS$) detection of a multi-planet system orbiting the $V=10.9$ K0 dwarf TOI 125. We find evidence for up to five planets, with varying confidence. Three high signal-to-noise transit signals correspond to sub-Neptune-sized planets ($2.76$, $2.79$, and $2.94\ R_{\oplus}$), and we statistically validate the planetary nature of the two inner planets ($P_b = 4.65$ days, $P_c = 9.15$ days). With only two transits observed, we report the outer object ($P_{.03} = 19.98$ days) as a high signal-to-noise ratio planet candidate. We also detect a candidate transiting super-Earth ($1.4\ R_{\oplus}$) with an orbital period of only $12.7$ hours and a candidate Neptune-sized planet ($4.2\ R_{\oplus}$) with a period of $13.28$ days, both at low signal-to-noise. This system is amenable to mass determination via radial velocities and transit timing variations, and provides an opportunity to study planets of similar size while controlling for age and environment. The ratio of orbital periods between TOI 125 b and c ($P_c/P_b = 1.97$) is slightly smaller than an exact 2:1 commensurability and is atypical of multiple planet systems from $Kepler$, which show a preference for period ratios just $wide$ of first-order period ratios. A dynamical analysis refines the allowed parameter space through stability arguments and suggests that, despite the nearly commensurate periods, the system is unlikely to be in resonance.

Tables

Table 1. Table 1: Ground-based Photometry of TOI 125

Planet	Facility	Filter	Typeaa F: full transit (covering ingress and egress); I: ingress only; E: egress only	$T_c$
				BTJDTDBbb Times of conjunction are given in the standard TESS-reported format, which is BJD${}_{\mathrm{TDB}}{}^{}-2457000$.
b	TRAPPIST-S $0.6$-mcc TRAPPIST (Jehin et al. 2011); SPECULOOS (Delrez et al. 2018)	z′	F	1378.6245
b	LCO-SS $0.4$-m	i′	I	1383.2783
b	SLR2 $0.5$-m	V	E	1383.2783
b	SLR2 $0.5$-m	V	I	1392.5860
b	LCO-SAAO $1.0$-m	zs	I	1392.5860
b	LCO-SAAO $1.0$-m	i′	I	1406.5475
c	MKO CDK700 $0.7$-m	r′	I	1371.0589
c	LCO-SAAO $1.0$-m	zs	F	1389.3608
c	IRSF $1.4$-m	J, H, Ks	F	1398.5117
c	SSO/Europa $1.0$-mcc TRAPPIST (Jehin et al. 2011); SPECULOOS (Delrez et al. 2018)	z′	F	1407.6626
.03	LCO-CTIO $1.0$-m	i′	I	1442.7549

Table 2. Table 2: CORALIE radial velocities of TOI 125

${\mathrm{BJD}}_{\mathrm{TDB}}$	RV	${\sigma }_{\mathrm{RV}}$	BIS	${\sigma }_{\mathrm{BIS}}$
	km s-1	km s-1	km s-1	km s-1
$2458368.687418$	$11.071$	$0.011$	$-0.082$	$0.011$
$2458379.908910$	$11.047$	$0.018$	$-0.077$	$0.018$
$2458393.713423$	$11.064$	$0.013$	$-0.073$	$0.013$

Table 3. Table 3: Literature Properties for TOI 125

	TIC 52368076
Other	TYC 88956-00192-1
identifiers	2MASS J01342273-6640328
	Gaia DR2 4698692744355471616
Parameter	Description	Value	Source
${\alpha }_{J2000}$	Right Ascension (RA)	01:34:22.735	1
${\delta }_{J2000}$	Declination (Dec)	-66:40:32.95	1
$T$	TESS $T$ mag.	$10.138\pm 0.017$	6
$B_T$	Tycho BT mag.	11.882 $\pm$ 0.077	2
$V_T$	Tycho VT mag.	11.102 $\pm$ 0.065	2
$B$aa The uncertainties of the photometry have a systematic error floor applied. Even still, the global fit requires a significant scaling of the uncertainties quoted here to be consistent with our model, suggesting they are still significantly underestimated for one or more of the broad band magnitudes.	APASS Johnson $B$ mag.	11.701 $\pm$ 0.025	3
$V$	APASS Johnson $V$ mag.	10.892 $\pm$ 0.016	3
$G$	Gaia $G$ mag.	10.7180 $\pm$ 0.0004	1
$g^{\prime }$	APASS Sloan $g^{\prime }$ mag.	11.268$\pm$0.019	3
$r^{\prime }$	APASS Sloan $r^{\prime }$ mag.	10.458$\pm$0.041	3
$i^{\prime }$	APASS Sloan $i^{\prime }$ mag.	10.662$\pm$0.017	3
$J$	2MASS $J$ mag.	9.466 $\pm$ 0.02	4
$H$	2MASS $H$ mag.	9.112 $\pm$ 0.03	4
$K_S$	2MASS $K_S$ mag.	8.995 $\pm$ 0.02	4
WISE1	WISE1 mag.	8.945 $\pm$ 0.03	5
WISE2	WISE2 mag.	9.006 $\pm$ 0.03	5
WISE3	WISE3 mag.	8.944 $\pm$ 0.03	5
WISE4	WISE4 mag.	8.613 $\pm$0.262	5
${\mu }_{\alpha }$	PM in RA (mas yr-1)	-119.800 $\pm$ 0.066	1
${\mu }_{\delta }$	PM in DEC (mas yr-1)	-122.953 $\pm$ 0.080	1
$\pi$	Parallax (mas)	8.976 $\pm$ 0.036	1
$RV$	Systemic RV ( km s-1)	11.062 $\pm$ 0.012	7

Table 4. Table 4: TOI 125 stellar parameters: median values and 68% CI

Parameter	Units	Values
Stellar Parameters:
$M_{\ast }$	Mass ($M_{☉}$)	${0.871}_{-0.040}^{+0.046}$
$R_{\ast }$	Radius ($R_{☉}$)	${0.852}_{-0.016}^{+0.017}$
$L_{\ast }$	Luminosity ($L_{☉}$)	${0.509}_{-0.025}^{+0.024}$
${\rho }_{\ast }$	Density (cgs)	$1.99\pm 0.15$
$\log g$	Surface gravity (cgs)	$4.518\pm 0.027$
$T_{\mathrm{eff}}$	Effective Temperature (K)	${5282}_{-75}^{+67}$
$\left[\mathrm{Fe}/H\right]$	Metallicity (dex)	${0.069}_{-0.081}^{+0.083}$
$Age$	Age (Gyr)	${6.6}_{-4.2}^{+4.6}$
$EEP$	Equal Evolutionary Point	${348}_{-27}^{+28}$
$A_V$	V-band extinction (mag)	${0.024}_{-0.017}^{+0.019}$
${\sigma }_{SED}$	SED photometry error scaling	${2.45}_{-0.62}^{+1.1}$
$\pi$	Parallax (mas)	${8.975}_{-0.100}^{+0.099}$
$d$	Distance (pc)	${111.4}_{-1.2}^{+1.3}$

Table 5. Table 5: TOI 125 planetary and transit parameters: median values and 68% confidence interval

Parameter	Description (Units)	Values
		Low SNR	Validated	Validated	Marginal	High SNR
		.04	b	c	.05aaThe values reported for TOI 125.05 are those from EXOFASTv2 without any radial velocity constraint. We note that when we exclude those solutions inconsistent with the CORALIE RVs, the most probable size is $R_P\sim 4.2$ $R_{\oplus }$ (see Fig. 4). Moreover, planets like TOI 125.05 are a priori likely to be small; if real, its true size is probably even smaller.	.03
$P$	Period (days)	${0.528474}_{-0.000030}^{+0.000040}$	${4.65382}_{-0.00031}^{+0.00032}$	${9.15067}_{-0.00069}^{+0.00062}$	${13.2781}_{-0.0019}^{+0.0020}$	${19.9807}_{-0.0049}^{+0.0045}$
$R_P$	Radius ($R_{\oplus }$)	${1.36}_{-0.16}^{+0.14}$	${2.755}_{-0.079}^{+0.091}$	$2.79\pm 0.10$	${8.8}_{-4.4}^{+4.7}$	$2.94\pm 0.16$
$T_C$	Time of conjunction (${\mathrm{BJD}}_{\mathrm{TDB}}$)	$2458350.8394\pm 0.0011$	${2458355.35520}_{-0.00087}^{+0.00093}$	${2458352.7582}_{-0.0013}^{+0.0014}$	${2458365.0560}_{-0.0020}^{+0.0019}$	${2458342.8514}_{-0.0033}^{+0.0034}$
$a$	Semi-major axis (AU)	${0.01222}_{-0.00019}^{+0.00021}$	${0.05210}_{-0.00082}^{+0.00090}$	${0.0818}_{-0.0013}^{+0.0014}$	${0.1048}_{-0.0016}^{+0.0018}$	${0.1376}_{-0.0022}^{+0.0024}$
$i$	Inclination (Degrees)	${72.80}_{-0.70}^{+0.72}$	${88.99}_{-0.81}^{+0.70}$	${88.52}_{-0.19}^{+0.32}$	${87.70}_{-0.14}^{+0.15}$	${88.753}_{-0.081}^{+0.080}$
$e$	Eccentricity	–	${0.183}_{-0.098}^{+0.14}$	${0.065}_{-0.046}^{+0.067}$	${0.037}_{-0.027}^{+0.046}$	${0.075}_{-0.051}^{+0.056}$
${\omega }_{\ast }$	Argument of Periastron (Degrees)	–	$-{91}_{-56}^{+57}$	${90}_{-98}^{+97}$	$50\pm 120$	${90}_{-110}^{+100}$
$ecos{\omega }_{\ast }$		–	$-0.00\pm 0.17$	$-{0.000}_{-0.054}^{+0.055}$	${0.000}_{-0.034}^{+0.036}$	$0.000\pm 0.063$
$esin{\omega }_{\ast }$		–	$-{0.114}_{-0.098}^{+0.057}$	${0.014}_{-0.042}^{+0.075}$	${0.001}_{-0.031}^{+0.037}$	${0.011}_{-0.053}^{+0.072}$
$⟨F⟩$	Incident Flux (109 erg s-1 cm-2)	$4.63\pm 0.25$	${0.243}_{-0.019}^{+0.017}$	${0.1026}_{-0.0055}^{+0.0056}$	${0.0627}_{-0.0033}^{+0.0034}$	${0.0362}_{-0.0019}^{+0.0020}$
$T_{eq}$	Equilibrium temperature (K)	${2126}_{-29}^{+28}$	$1029\pm 14$	$821\pm 11$	${725.8}_{-9.8}^{+9.5}$	${633.5}_{-8.5}^{+8.3}$
$R_P/R_{\ast }$	Radius of planet in stellar radii	${0.0146}_{-0.0016}^{+0.0014}$	${0.02962}_{-0.00063}^{+0.00077}$	${0.02998}_{-0.00086}^{+0.00087}$	${0.095}_{-0.047}^{+0.051}$	${0.0317}_{-0.0016}^{+0.0015}$
$a/R_{\ast }$	Semi-major axis in stellar radii	$3.085\pm 0.077$	$13.15\pm 0.33$	$20.65\pm 0.52$	$26.47\pm 0.66$	$34.75\pm 0.87$
$d/R_{\ast }$	Separation at mid transit	$3.085\pm 0.077$	${14.29}_{-0.97}^{+1.3}$	${20.2}_{-1.5}^{+1.1}$	${26.4}_{-1.2}^{+1.1}$	${34.1}_{-2.5}^{+2.2}$
$\delta$	Transit depth (fraction)	${0.000214}_{-0.000045}^{+0.000043}$	${0.000877}_{-0.000037}^{+0.000046}$	${0.000899}_{-0.000051}^{+0.000053}$	${0.0090}_{-0.0067}^{+0.012}$	${0.001004}_{-0.000096}^{+0.00010}$
$\tau$	Ingress/egress transit duration (days)	${0.00206}_{-0.00039}^{+0.00049}$	${0.00376}_{-0.00024}^{+0.00095}$	${0.00488}_{-0.00078}^{+0.00070}$	${0.0232}_{-0.0015}^{+0.0018}$	${0.0086}_{-0.0014}^{+0.0016}$
$T_{14}$	Total transit duration (days)	${0.0255}_{-0.0027}^{+0.0022}$	${0.1233}_{-0.0026}^{+0.0025}$	${0.1227}_{-0.0028}^{+0.0025}$	${0.0464}_{-0.0030}^{+0.0037}$	${0.1284}_{-0.0057}^{+0.0055}$
$T_{FWHM}$	FWHM transit duration (days)	${0.0235}_{-0.0031}^{+0.0025}$	${0.1192}_{-0.0024}^{+0.0023}$	${0.1179}_{-0.0030}^{+0.0026}$	${0.0232}_{-0.0015}^{+0.0019}$	${0.1198}_{-0.0063}^{+0.0060}$
$b$	Transit Impact parameter	${0.912}_{-0.022}^{+0.023}$	${0.25}_{-0.17}^{+0.24}$	${0.524}_{-0.14}^{+0.078}$	${1.056}_{-0.057}^{+0.055}$	${0.745}_{-0.060}^{+0.047}$
$M_P$	Predicted Mass ($M_{\oplus }$)	${2.65}_{-0.56}^{+0.94}$	${8.5}_{-1.8}^{+2.8}$	${8.6}_{-1.9}^{+2.8}$	${61}_{-42}^{+200}$	${9.5}_{-2.1}^{+3.2}$
$K$	Predicted RV semi-amplitude (m/s)	${2.19}_{-0.47}^{+0.78}$	${3.65}_{-0.80}^{+1.2}$	${2.88}_{-0.64}^{+0.95}$	${18}_{-13}^{+59}$	${2.45}_{-0.56}^{+0.84}$
$logK$	Log of RV semi-amplitude	${0.34}_{-0.10}^{+0.13}$	${0.56}_{-0.11}^{+0.12}$	${0.46}_{-0.11}^{+0.12}$	${1.26}_{-0.52}^{+0.63}$	${0.39}_{-0.11}^{+0.13}$
${\rho }_P$	Predicted Density (cgs)	${5.9}_{-1.3}^{+2.0}$	${2.21}_{-0.47}^{+0.72}$	${2.17}_{-0.47}^{+0.70}$	${0.69}_{-0.36}^{+0.92}$	${2.04}_{-0.44}^{+0.68}$
$log{g}_P$	Predicted Surface gravity	${3.147}_{-0.092}^{+0.12}$	${3.04}_{-0.10}^{+0.12}$	${3.03}_{-0.10}^{+0.12}$	${2.94}_{-0.15}^{+0.26}$	${3.03}_{-0.10}^{+0.12}$
$M_P/M_{\ast }$	Predicted Mass ratio	${0.0000091}_{-0.0000020}^{+0.0000033}$	${0.0000291}_{-0.0000064}^{+0.0000098}$	${0.0000295}_{-0.0000066}^{+0.0000098}$	${0.00021}_{-0.00015}^{+0.00069}$	${0.0000325}_{-0.0000075}^{+0.000011}$
Wavelength Parameters:		TESS
$u_1$	linear limb-darkening coeff	$0.391\pm 0.037$
$u_2$	quadratic limb-darkening coeff	$0.231\pm 0.036$
Transit Parameters:		TESS Sector 1	TESS Sector 2
${\sigma }^2$	Added Variance	${4.0}_{-1.7}^{+1.8}\times {10}^{-8}$	$1.1\pm 1.7\times {10}^{-8}$
$F_0$	Baseline flux	$1.000033\pm 0.000013$	$1.000074\pm 0.000012$

Equations

$${\left(\frac{1}{e_b}\frac{d{e}_b}{dt}\right)}^{-1}=83\text{Myr}\times \left(\frac{Q/k_2}{10^3}\right)\left(\frac{m_b}{10M_{\text{earth}}}\right).$$

$$\delta =\frac{P_2}{P_1}-\frac{(q+1)}{q}>-\mathrm{0.04.}$$

Full text

planetary systems, planets and satellites: detection, stars: individual (TOI 125, TIC 52368076), planets and satellites: dynamical evolution and stability

††facilities: TESS, Euler1.2m (CORALIE), Gemini:South (DSSI), SOAR (HRCam), VLT:Antu (NACO), IRSF, LCOGT, SSO:1m, TRAPPIST††software: EXOFASTv2 (Eastman et al. 2013; Eastman 2017), vespa (Morton 2015), Mercury6 (Chambers 1999), AstroImageJ (Collins et al. 2017)

\AuthorCallLimit

=15

\noaffiliation

1 Introduction

NASA’s Transiting Exoplanet Survey Satellite (TESS) (Ricker et al. 2015) is a (nearly) all-sky survey, the primary objective of which is to discover and characterize transiting planets smaller than Neptune orbiting the nearest and brightest stars in the sky. While the space-based transit survey carried out by Kepler (Borucki et al. 2010) led to breakthroughs in our understanding of the occurrence rates of planetary systems (e.g., Fressin et al. 2013) and the various dynamical configurations of multi-planet systems (e.g., Lissauer et al. 2011b; Fabrycky et al. 2014), TESS is designed to discover the planetary systems most amenable to detailed characterization through follow-up observations. Typical TESS planet hosts will be several magnitudes brighter than those from Kepler (Sullivan et al. 2015; Barclay et al. 2018; Huang et al. 2018a), and these statistically rare systems are amenable to the most precise radial-velocity (RV) mass measurements and more efficient atmospheric characterization. This expectaction is borne out by the experience of the K2 mission (Howell et al. 2014), which surveyed a larger area of sky than Kepler and discovered a number of bright planetary systems like the ones expected from TESS (see, e.g., Rodriguez et al. 2018a, and the overview presented therein). TESS should also detect objects that are intrinsically rare, such as events occurring on astronomically short timescales, or the unlikely outcomes of dynamical interactions. Further study of these benchmark objects may lead to breakthroughs in our understanding of the fundamental processes that govern the formation and evolution of planetary systems.

Data from just the first TESS observing sector ($27.4$ days, or two spacecraft orbits) have already begun to fulfill this promise. A $2.1$-$\,R_{\rm\earth}$ planet transiting the 5th magnitude star $\pi$ Mensae already has a mass measurement ($4.8\ \,M_{\rm\earth}$) because the star was previously known to host a long-period giant planet and there exist extensive archival RV measurements (Huang et al. 2018b; Gandolfi et al. 2018). The star is one of the very brightest known to host a transiting planet, which will enable further detailed characterization. A second system (LHS 3844; Vanderspek et al. 2018) is an Ultra-Short Period (USP) Earth-sized planet ($R_{p}=1.32\ \,R_{\rm\earth}$) in an $11$-hour orbit around a late M dwarf $15$ pc away. It is one of the nearest planetary systems, and in many respects the USP planet that is most amenable to follow-up studies. The remaining $25$ observing sectors in the two-year prime TESS mission will survey additional bright stars, some for longer periods of time, and will lead to the discovery of many more benchmark planetary systems.

Among the myriad discoveries from the Kepler mission, the detection of systems of multiple transiting planets and the subsequent study of their ensemble properties remain among the most impactful results (e.g., Steffen et al. 2010; Lissauer et al. 2011a; Latham et al. 2011). Multi-planet transiting systems allow investigations of formation and evolution processes through measurements of mutual inclinations, adjacent planet sizes, planet spacings, stellar obliquities, mass measurements via transit timing variations, and more. Given the prevalence of multi-transiting systems, TESS will build upon the Kepler legacy by discovering the nearest and brightest such systems, as well as the rare examples.

In this paper, we present the discovery and validation of a system of multiple transiting planets orbiting the star TIC 52368076, which has been assigned TESS Object of Interest (TOI) number TOI 125. The proposed architecture of the system is illustrated in Figure 1. We identify three candidates with high signal-to-noise (SNR) transits (filled circles), as well as two low-SNR candidates (open circles), one of which is an ultra-short-period terrestrial candidate.

We describe the analysis of data from TESS, ground-based follow-up, and archival imaging in Section 2, fit a global model to all available data in Section 3, present a statistical validation of the planets in Section 4, investigate the dynamics in Section 5, and discuss the properties of the system and prospects for future characterization in Section 6.

2 Data & Analysis

2.1 TESS Photometry

TIC 52368076 was observed by TESS in the first two sectors of the prime mission (2018 Jul 25 through 2018 Sep 21), on CCD 1 of Camera 3 in Sector 1 and CCD 2 of Camera 3 in Sector 2. The CCDs obtain images at a two-second cadence, which are summed on board the spacecraft to produce images with the appropriate effective exposure time. All stars within the TESS field of view are observed with an effective exposure time of $30$ minutes, but a subset of stars (including TIC 52368076) were pre-selected, primarily on the basis of planet detectability (Stassun et al. 2018), for data to also be returned to Earth at a two-minute cadence.

The two-minute data were reduced with the Science Processing Operations Center (SPOC) pipeline (Jenkins et al. 2015, 2016), adapted from the pipeline for the Kepler mission at the NASA Ames Research Center (Jenkins et al. 2010). Two transit signals were strongly detected with periods of $4.65$ and $9.15$ days and signal-to-noise ratios (SNR) of $20.1$ and $16.4$  respectively. These candidates were assigned identifiers TOI 125.01 and TOI 125.02 by the TESS team. An additional signal, TOI 125.03, was detected with only two transit-like events at a period of $19.98$ days and SNR of $9.8$. In the analyses that follow, we use the Pre-search Data Conditioning (PDC) light curve from SPOC (see Stumpe et al. 2012). Figure 2 shows the PDC light curve after flattening (we note that the raw PDC light curve looks nearly identical, as the star is photometrically very quiet). Interruptions in data acquisition occur at the perigee of each TESS orbit (once every 13.7 days) and last approximately 1 day, during which time the spacecraft reorients to downlink data. The second orbit of Sector 1 included a two-day period during which the data were of lesser quality due to a one-time occurrence of abnormally unstable spacecraft pointing. The worst of these data were flagged by SPOC and removed, which can be seen as an under-density of points in Figure 2 just before BTJD 1350. Fortuitously, none of the transits of these three candidates occurred at this time. During Sectors 1 and 2, the spacecraft thrusters were fired periodically (approximately every 2.5 days) to reduce the speed of the reaction wheels, allowing them to operate at frequencies that introduced less pointing jitter. In $10$- to $15$-minute intervals around these "momentum dumps", we removed data from our analysis.

2.2 Additional Planet Candidates

Following the convention from Kepler, we adopt a formal significance threshold of $7.1\sigma$, and the three candidate planets described above are the only formally significant periodic signals in the data. However, we do detect lower SNR transit-like signals at two other periods. The first, TOI 125.04, has a period of $0.52854$ days and a depth of $180$ ppm ($5.2$-$\sigma$), which corresponds to a planet radius of $1.36$ $\,R_{\rm\earth}$ (see Figure 3). The second, TOI 125.05, was detected at a period of $13.2780$ days and a depth of $675$ ppm ($5.1$-$\sigma$). If this were a central transit, it would correspond to a planet radius of $\mathord{\sim}2.2$ $\,R_{\rm\earth}$, but it is best modeled as a grazing transit (see Figure 4), so when unconstrained by other data even giant planets are allowed. However, our CORALIE spectroscopy (Section 2.6) rules out the largest companions, and we ultimately derive a radius of $\mathord{\sim}4.2$ $\,R_{\rm\earth}$, still with large error bars. Given the low SNR of these candidates and the non-Gaussian noise of the first TESS Sectors, we do not consider the signals strong enough to be validated as planets, particularly TOI 125.05, which only shows four transits of varying quality. Nevertheless, the presence of three other strong planet candidates makes these signals more intriguing, and we note them here so that they can be taken into account during follow-up observations and subsequent analysis of the first three candidates.

2.3 Ground-based Photometry

Given the $21\arcsec$ TESS pixels and a PSF a few pixels wide, the light from an individual star on the detector can extend well beyond $1\arcmin$. In order to capture most of the light from the target star, the TESS photometric apertures must also be large. Therefore, even apparently isolated stars may be contaminated by relatively distant neighbors, with the exact contamination fraction depending upon the aperture choice and the magnitudes of the stars (see, for example, the size of the PSF in the TESS image of TOI 125 shown in Figure 5). If a neighboring star is an eclipsing binary (EB), deep eclipses can be diluted to resemble shallow planetary transits. While previous experience with Kepler candidate multi-planet systems shows the vast majority to be real planets (e.g., Lissauer et al. 2012), the larger TESS pixels and aperture create more opportunity for unassociated EBs to contaminate real planetary systems (producing candidate multi-planet systems containing both real planets and false positives). Centroid analysis of the TESS difference images (comparing the in-transit to out-of-transit flux) is often effective at identifying nearby EBs, but transit signals with a small number of events or contaminants within about a pixel might not be robustly detected. We therefore observed TOI 125 with ground-based facilities at predicted times of transit to search for deep eclipses in nearby stars. We enumerate these observations in Table 1. In order to produce the $1$ mmag events of TOI 125 b, c, or .03 with even a $50$% eclipse, a nearby star must be no more than $6.9$ magnitudes fainter (and fully blended in the TESS aperture). Among Gaia DR2 sources, only one nearby star is bright enough, but at a distance of $75\arcsec$ (to the SSW), it is not fully blended within the TESS aperture (see Figure 5). Nonetheless, we search this and other nearby stars for evidence of deep eclipses and we find no indication of contamination from nearby EBs in our timeseries observations.

2.4 Archival Imaging

Ground-based photometric follow-up can rule out the presence of EBs at modest separations, but a physically unassociated background star within a few arcseconds of the location of TOI 125 could plausibly produce a transit-like signal that would not be resolved as a separate source in the few-arcsecond PSFs of follow-up images. To address this possibility, we examine archival images, in which the proper motion of TOI 125 has carried it away from its current location. Figure 5 shows the TESS image from Sector 1 along with images from the ESO/SERC Southern Sky Atlas (SERC-J; taken in 1975) and the Anglo-Australian Observatory Second Epoch Survey (AAO-SES; 1993). The most constraining of these is the SERC-J image, enlarged in the left panel. The proper motion of $\mu_{\alpha}=-120$ ${\rm mas\,yr}^{-1}$ and $\mu_{\delta}=-123$ ${\rm mas\,yr}^{-1}$ leads to motion of $1\farcs 7$ per decade; in the $43$ years since the SERC-J image was obtained, TOI 125 has moved $7\farcs 4$. A background source at the current location of TOI 125 should be seen as elongation of the PSF or a nearly resolved source. There is no indication of such features in either the blue-sensitive SERC-J or the red-sensitive AAO-SES images, so we conclude that there is no background source coincident with the present-day location of TOI 125.

2.5 High Angular Resolution Imaging

Ground-based photometry rules out EBs at modest separations and archival imaging rules out background sources, but there may still be a bound stellar companion at small angular separation. An unresolved companion may itself be an EB responsible for one of the transit-like signals, but even if it is not, the dilution must be taken into account in the light curve fit in order to derive accurate radii (e.g., Buchhave et al. 2011), and the presence (or absence) of a binary companion can help us understand the formation of compact planetary systems (e.g., Kraus et al. 2016). Fortunately, bound companions to TESS planet hosts will be more easily revealed by high-resolution imaging than the typical Kepler system because they are, on average, more nearby (e.g., Ciardi et al. 2015; Matson et al. 2018).

We searched for close companions to TOI 125 in $I$-band using the HRCam speckle imager on the $4.1$-m Southern Astrophysical Research (SOAR) telescope (Tokovinin 2018; Ziegler et al. 2018b) on 2018 September 25 UT, in narrow-band $Br\gamma$ using the NaCo adaptive optics imager (Rousset et al. 1998; Lenzen et al. 1998) on the 8-m UT1 of the VLT on 2018 October 23, and simultaneously in $R$- and $I$-band using the DSSI speckle imager (Horch et al. 2009, 2012) on the 8-m Gemini South Telescope on 2018 October 31 UT. We detected no companions in any of these images down to contrast ratios of more than $5$ magnitudes outside of $0\farcs 2$ of TOI 125. Outside of $1\farcs 5$, Gaia DR2 can exclude the presence of stellar sources bright enough to cause the $\mathord{\sim}1$ mmag transit signals when blended with TOI 125 (Ziegler et al. 2018a). The $5$-$\sigma$ contrast curves are shown in Figure 6.

2.6 Reconnaissance Spectroscopy

We obtained three spectra with the CORALIE spectrograph (Queloz et al. 2000; Pepe et al. 2017) on the Swiss Euler 1.2-m telescope of the ESO-La Silla Observatory (Chile) between UT 2018 Sept 07 and 2018 Oct 02. We used simultaneous Fabry-Pérot calibration for intrinsic drift measurement. The SNR per pixel of the individual spectra was $\mathord{\sim}20$. Data were reduced using an adapted version of the HARPS pipeline: the average stellar line profiles, or cross-correlation functions (CCFs), were computed by cross-correlating the CORALIE spectra with a weighted binary G2 mask from which various tellurics and ISM lines were removed (Pepe et al. 2002). We see no evidence for multiple peaks in the CCF, suggesting that TOI 125 does not have a bright, unresolved stellar companion. The RVs, reported in Table 2, show no significant velocity variation. We derive spectroscopic parameters using SpecMatch (Petigura et al. 2015; Yee et al. 2017), and find $T_{\rm eff}=5187\pm 110$ K, $\log g=4.52\pm 0.12$, $\left[{\rm Fe}/{\rm H}\right]=0.06\pm 0.09$, and ${v\sin{i_{\star}}}<2$  km s*-1*. We use these values as starting guesses for our global model, and apply the derived $\left[{\rm Fe}/{\rm H}\right]$ as a prior (see Section 3).

3 EXOFASTv2 Global Fit

To gain a full understanding of the system parameters, we globally fit the available photometric and spectroscopic data using the publicly available exoplanet modeling suite, EXOFASTv2 (Eastman et al. 2013; Eastman 2017). Specifically, we fit the TESS light curves from observing Sectors 1 and 2 for planets b, c, and candidates .03, .04, and .05 (See Figures 2, 3, and 4), while constraining the host star parameters using the spectral energy distribution (SED) and the MESA Isochrones and Stellar Tracks (MIST) stellar isochrones (Dotter 2016; Choi et al. 2016; Paxton et al. 2011, 2013, 2015). The broadband photometry is given in Table 3 and shown along with the final model in Figure 7. We enforce Gaussian priors $T_{\rm eff}=5187\pm 110$ K and $\left[{\rm Fe}/{\rm H}\right]=0.06\pm 0.09$ dex from the analysis of the CORALIE spectra. We also place a conservative Gaussian prior on the parallax from Gaia DR2 of $8.976\pm 0.1$ mas, since all possible uncertainties should total to less than 0.1 mas (Gaia Collaboration et al. 2016, 2018). Lastly, we enforce an upper limit on the extinction of AV = 0.0521 from the Schlegel Galactic dust reddening and extinction maps (Schlegel et al. 1998). All other parameters were allowed to vary without prior constraints. We allowed an error scaling term for the SED photometry (reported in Table 4) and a variance term for each sector of the TESS photometry (Table 5). Limb darkening parameters are interpolated using the current $\log g$, $T_{\rm eff}$, and $\left[{\rm Fe}/{\rm H}\right]$ at each step and the limb-darkening tables of Claret (2017). We adopt the strict convergence criteria recommended by Ford (2006) in order to assure the global minimum has been identified and covariances are well characterized: the Gelman-Rubin statistic for all parameters must be less than $1.01$, and the number of independent draws (chain length divided by correlation length) must exceed $1000$. We ran a fit that allowed transit timing variations (TTVs) but find no significant TTVs and no changes to the derived parameters, so for simplicity we adopt the solution that assumes periodic ephemerides. We also ran fits including only the two, three, or four strongest signals, and we find stellar and planetary results fully consistent with the five planet solution. Because inclusion of the two more marginal candidates does not affect our conclusions about TOI 125 b, c, and .03, we present herein the 5-planet fit. The final system parameters determined by the EXOFASTv2 fit, including predicted masses using the relations of (Chen & Kipping 2017), are shown in Tables 4 and 5. We again refer the reader to the top view of the system architecture in Figure 1.

4 Statistical Validation with VESPA

We used the vespa package (Morton 2015) to assess the statistical likelihood that the transits of TOI 125 b and c are caused by planets, rather than false positives. vespa simulates stellar and planetary systems to generate transits (and eclipses) to compare against the observed data of TOI 125. Rejecting systems that are inconsistent with the observations, vespa then calculates the false positive probability (FPP) for each candidate. The FPP depends not only on the transit shape, but also the position of the star on the sky (to assess the likelihood of background blends), the stellar parameters (which inform not only transit and eclipse shapes but also the likelihood of stellar companions), the extent to which nearby EBs can be excluded by high resolution imaging, and the presence or absence of features in the light curve that might indicate the presence of a binary (such as depth differences in alternating transits or the presence of a secondary eclipse). We therefore provide to vespa the sky coordinates, stellar parameters, literature photometry, high resolution imaging contrast curves, and the flattened TESS light curve. Our RVs rule out an EB (as opposed to a blended background or hierarchical EB).

After running vespa, we adjusted the FPP by excluding the scenario in which a direct EB companion to TOI 125 causes one of the transit signals. The resulting FPPs are $6\times 10^{-5}$ and $9\times 10^{-5}$ for TOI 125.01 and .02, respectively. Note that we did not remove the contribution from background EBs even though our inspection of the archival imaging suggests there are no background stars at the current location of TOI 125; the FPPs would be even lower with this adjustment. We therefore conclude that these are statistically validated planets, and now refer to them as TOI 125 b and c. We do not attempt to validate TOI 125.03 despite its high SNR, because we observed only two transits. Similarly, we do not attempt to validate TOI 125.04 and .05 because vespa does not assess the likelihood that a signal is an instrumental false alarm, and we cannot fully exclude this possibility given the low SNR of these events.

5 Dynamics

This section considers dynamical aspects of the TOI 125 system, considering the three planets with the highest S/N detections. We first note that TOI 125 b has non-zero eccentricity, but a relatively short time scale for the damping of its eccentricity. As a result, the system is dynamically interesting, suggesting some type of planet-planet interactions. These types of interactions could lead to TTVs (Section 5.1), although they are not observed in the present data set. We also need to consider possible instability of the system (Section 5.2), which puts additional constraints on the allowed current-day values of the orbital eccentricities. Finally, we consider the question of mean motion resonances, and find that the system is highly unlikely to currently be in a true resonant configuration (Section 5.3).

5.1 Transit Timing Variations (TTVs)

The proximity of TOI 125 b and c to a 2:1 MMR means that their mutual gravitational perturbations add in a nearly-coherent manner that can lead to significant and potentially-measurable TTVs (e.g., Agol et al. 2005; Holman & Murray 2005). An attempt to model the planets’ timing variations in the first two sectors of TESS did not yield significant dynamical constraints, as the uncertainty on the transit times from TESS are larger than the expected TTVs. Here we briefly consider the prospects for extracting dynamical information by combining the first two sectors of TESS observations with future follow-up observations.

To assess the possibility for extracting dynamical information from TOI 125 b and c’s TTVs, we need to know how the expected TTV signal amplitude compares to the precision of any future transit time measurements. We employ the analytic TTV model detailed in Hadden et al. (2018) in order to predict the planets’ TTV signals and the accompanying dynamical constraints derived from them. Figure 8 shows the TTV signals predicted for b and c, assuming no free eccentricity for either planet111 Due to its short orbital period, the eccentricity damping timescale of planet b is short. Adopting the best-fit stellar mass, planet radius, and orbital period from Table 5, the tidal eccentricity damping timescale for TOI 125 b is given by

$\displaystyle\left(\frac{1}{e_{b}}\frac{de_{b}}{dt}\right)^{-1}=83~{}\text{Myr}\times\left(\frac{Q/k_{2}}{10^{3}}\right)\left(\frac{m_{b}}{10\,M_{\rm\earth}}\right).$

(1)

where $Q$ is the planet’s tidal quality factor and $k_{2}$ its tidal Love number (Goldreich & Soter 1966). The parenthetical terms on the right hand side are of order unity. With the same assumptions about $Q/k_{2}$ for planet c, its nominal eccentricity damping timescale is $1.9~{}\text{Gyr}$. Dynamical coupling between b and c should enhance the efficiency with which planet c’s eccentricity is damped.

and fiducial masses of $3\times 10^{-5}\,M_{*}$ for both planets. The resulting TTVs are approximately sinusoidal with amplitudes of $\mathord{\sim}3$ minutes. The planets’ TTV signals will, to excellent approximation, simply scale linearly with the planet masses. To estimate the timing precision one could achieve using Spitzer, we scale SNR from existing transit observations and employ the analytic formulae of Carter et al. (2008) and Price & Rogers (2014), which lead to a per-transit uncertainty of $1$ min. The bottom panel of Figure 8 shows the estimated precision with which planet c’s mass could be measured using planet b’s TTV by obtaining a series of follow-up observations centered on the peaks of the approximately sinusoidal TTV and assuming transit mid-times are measured with $1\sigma$ uncertainties of $1$ minute. Precisions approaching $\sim 1M_{\earth}$ are achievable with a series of transit observations consisting of $\sim 3-5$ transit timing measurements made at 2 or 3 successive peaks of the TTV signal. Note, though, that these mass measurement precisions are based on a TTV model that assumes zero free eccentricity for TOI 125 b and c; relaxing this assumption would significantly increase the mass uncertainty due to the mass-eccentricity degeneracy inherent to TTVs of planets near MMRs (e.g., Lithwick et al. 2012). On the other hand, the combination of TTVs and RVs would provide the strongest possible mass and eccentricity constraints, and in Section 6.1 we argue that the planets should be amenable to RV follow-up.

5.2 Dynamical Stability

In the absence of external perturbations, a planet with an eccentric orbit residing close to its host star would generally become tidally circularized on astronomically short timescales (the exact timescale depending on the tidal quality factor of the planet and properties of its host star; see, e.g., Equation 1). Therefore, significant eccentricities for compact systems often require that the planets be located in regions of resonance (Beaugé et al. 2006) and can result in significant transfer of angular momentum between planets (Kane & Raymond 2014; Antoniadou & Voyatzis 2016). Thus, a dynamical analysis of a proposed orbital solution can be used to validate or potentially revise the allowed architecture for the system. The EXOFASTv2 global model of the TOI 125 multi-planet system cannot exclude a moderately large eccentricity of TOI 125 b ($\sim 0.18$), which would be unusual for a planet in a tightly-packed, multi-planet system (Kane et al. 2012; Hadden & Lithwick 2014; Van Eylen & Albrecht 2015; Xie et al. 2016). We therefore set out to investigate whether such high eccentricities can be ruled out through dynamical simulations, or if there exists evidence that TOI 125 b and TOI 125 c may be in resonance.

In the analysis that follows, we consider only the three high-SNR transit signals (TOI 125 b, c, and .03). The candidate USP planet (TOI 125.04), given its small size and orbital period, is effectively dynamically decoupled from the other planets and if it exists should not affect our conclusions. TOI 125.05 would have to be incorporated in a dynamical investigation if it is shown to be real, but given its low SNR, we do not include it in our analysis for now. If follow-up observations show that TOI 125.05 is real, it would reside slightly interior to the 3:2 mean motion commensurability (with a period ratio of $\sim 1.45$), close to the regime in which resonant interactions are most relevant.

To evaluate the dynamical stability and orbital evolution of the three planets in the TOI 125 system, we performed 3000 numerical simulations using N-body code Mercury6 (Chambers 1999), altered to include the effect of general relativistic precession. The simulations were performed using a hybrid symplectic and Bulirsch-Stoer (B-S) integrator with a time-step of 30 minutes for a total integration time of slightly more than 1 million years per integration (which is roughly 80 million orbits of the innermost planet), and energy was conserved to better than one part in $10^{9}$ (for energy changes due to the integrator). For each integration, we drew one link from the EXOFASTv2 transit fit posterior, and assigned planet and star properties equal to those in the chosen posterior link.

These numerical simulations allow us to impose an additional level of constraints beyond those derived from the transit shapes: some planetary eccentricities will lead to dynamical instabilities in the system (which occur when scattering events or orbit crossing leads to physical collisions between planets, the ejection of a planet from the system, or collision of a planet with the central star). The composite eccentricity distributions (stacked stable + unstable) in Figure 9 show the eccentricity draws for our 3000 simulations. The stable subset of each distribution contains those which allow the planets to remain dynamically stable for the entire 1 Myr integration. Not shown is the variation in other orbital elements (also drawn from the EXOFASTv2 posteriors), but the overall trend in stability fraction is shown in the top panels.

Of these 3000 simulations, 32% remained dynamically stable for the entire 1 Myr integration, and Figure 9 show a higher stability fraction at lower values of eccentricity. Dynamical stability considerations thus prefer the lower values of eccentricity, with eccentricities above 0.25–0.3 disallowed for each planet. The conclusion from this analysis is that although the EXOFASTv2 posteriors allow an unusually large range of eccentricities, the true values are likely on the lower ends of these ranges.

5.3 Resonant Behavior and Formation History

The orbital periods of the two planets b and c (4.65437 and 9.1536 days) lie close to the mutual 2:1 commensurability. As a result, it is natural to wonder whether the two planets are trapped in mutual 2:1 mean motion resonance (with a period ratio of 1.967). As orbital elements (including semi-major axis) librate while planets reside in resonance, it is possible to reside in resonance even without a perfect 2:1 period ratio (Batygin & Morbidelli 2013, i.e.,).

However, the majority of planets near, but not in, orbital resonance reside slightly outside of a resonant configuration. The results of Terquem & Papaloizou (2018, see Figure 9) show that for planets migrating in a disk, it is very easy for the 2:1 resonance to be disrupted when the inner planet enters a disk cavity, as might occur at small orbital radii. The subsequent evolution of the system is more difficult to model, but it is expected that departures from resonance will move towards the outside of resonance in cases where

[TABLE]

when the resonance considered is the $q+1:q$ first order resonance. Since the inner two planets in the TOI 125 system have $\delta$ = -0.033, they would be expected to fit this trend and reside slightly outside of the 2:1 MMR period ratio; however, these planets instead appear to have a period ratio slightly less than this value. There are three potential explanations for this:

1) These planets are currently in true orbital resonance. A system with a similar architecture to TOI 125 is the Gliese-876 system (Marcy et al. 2001; Rivera et al. 2010), which also has two planets close to the 2:1 resonance (Gliese-876 c and b, at $\sim 30$ and $\sim 60$ days), the inner of which exhibits significant (0.2) eccentricity. In the Gliese-876 system, those two planets form a Laplace resonance with a third planet. The non-zero planetary eccentricities are maintained through the resonance. The TOI 125 system also has three high-SNR planets with orbital periods moderately close to 2:1 resonances. An orbital resonance which persisted after the disk dissipated could explain the large eccentricity of the inner planet. While the disk is present, the orbital eccentricity of both planets in the system will be damped. Once the disk has dissipated, the eccentricities of the resonant pair is free to grow (and experience secular cycles) due to interactions with other extra planets in the system (in the TOI 125 system, there appear to be several additional planets; additionally, an increased eccentricity for either planet involved in the resonance may change the resonance width and lead eventually to the disruption of the resonance; Wittenmyer et al. 2012; Malhotra et al. 2018).

Beaugé et al. (2003) presents solutions, inspired by Gliese-876, for stable aligned pericenters in the 2:1 resonance. Notably, stable solutions must have non-null eccentricities for the inner planet, and for an inner eccentricity of 0.3 or so, the outer planet must have a non-null eccentricity as well (for equal mass planets). The EXOFASTv2 posteriors indicate that TOI 125 c could have a non-zero eccentricity; however, the pericenters are not well enough constrained to determine whether this aligned scenario occurs.

Using our suite of N-body simulations, we can evaluate the fraction of fitted posteriors which are consistent with a resonant configuration for planets b and c. Of the 3000 trials considered, 32% of them are dynamically stable. Of this stable subset, only one began in a true resonance (as defined by a librating resonant angle) and remained so for the entire integration. However, some of the integrations show that the planets can attain and lose a resonant configuration through their natural orbital evolution: in 7% of the dynamically stable integrations, TOI 125 b and c attain a true 2:1 mean motion resonance for at least some of the integration (generally for periods around $10^{5}$ years at a time) but are subsequently disrupted from that resonance. A further 1% attain resonance during the integration and remain stable in that resonance for the entire rest of the 1 Myr integration. However, 86% of the dynamically stable integrations never attain a resonant configuration. Barring one single integration that, after a scattering interaction, attained the 5:3 true resonance, the remaining $\sim$5% of the stable simulations exhibit (for at least some of the integration, but not a majority) a ‘nodding’ behavior (i.e., Ketchum et al. 2013) in and out of the 2:1 resonance.

From this suite of simulations, it appears that the vast majority of the dynamically stable posteriors are fully non-resonant. However, a true orbital resonance could explain both the eccentricity of the inner planet and the continued stability of the system. The simulations show that this system, if in resonance, is likely characterized by non-consistent attainment of true resonance.

2) These planets formed in-situ or via scattering inwards, and do not have resonance in their history. Terquem & Papaloizou (2018) note that only roughly 15% of systems are consistent with smooth, disk-driven migration, which results in systems with $0<\delta<$ 0.04 (the ’outside of resonance’ population common in the observational sample). If TOI 125 was in resonance while the disk was still present (required if it assembled via disk-driven migration), it should have moved toward positive $\delta$ while in resonance and ended with orbits consistent with this population. Its small, negative value of $\delta$ can be explained if the system did not assemble via smooth migration, is not in resonance, and reached its current proximity to resonance by chance.

3) These planets formed via disk-driven migration and were in resonance, but are no longer in resonance. As discussed in Adams et al. (2008) and Batygin & Adams (2017), turbulent fluctuations in the disk can destabilize resonances for small planets. These planets are both 2.7 Earth radii, slightly larger than should be possible to turbulently force out of resonance according to Batygin & Adams (2017). However, the 2:1 resonance is rather weak. As such, the resonant angles for this resonance generally have a large libration amplitude, potentially permitting either liberation from true resonance with minor perturbations, or large excursions in orbital element libration. As demonstrated by the numerical simulation, a sizable fraction of the posteriors attain and subsequently lose the 2:1 resonance during the integrations (sometimes multiple times).

6 Discussion

Of the $4723$ Kepler planets and planet candidates discovered to date222Kepler Objects of Interest reported as Confirmed or Candidate planets on the NASA Exoplanet Archive on 2018 Oct 22., the majority are smaller than Neptune and larger than the Earth, and orbit within a few tenths of an AU (e.g., Thompson et al. 2018), a class of planet not seen in the Solar System. The occurrence rates of these short-period super-Earths and mini-Neptunes indicate that they are common byproducts of star formation (e.g., Fressin et al. 2013; Petigura et al. 2013; Dressing & Charbonneau 2015). As such, their physical and orbital properties hold a wealth of information about the processes governing planet formation and evolution that were previously unconstrained by observation. The ensemble properties of these planets have revealed some of their fundamental characteristics (e.g., Fulton et al. 2017; Fulton & Petigura 2018; Berger et al. 2018), and detailed investigations of individual systems can complement the information gained from population studies.

6.1 Radial Velocity Characterization

The observed bimodality in the radius distribution of Kepler planets (e.g., Owen & Wu 2013; Zeng et al. 2017; Fulton et al. 2017), with peaks at $\mathord{\sim}1.3$ and $2.4$ $\,R_{\rm\earth}$, can be reproduced theoretically from the photoevaporation of close-in, low-mass planets, which are stripped to their bare ($\mathord{\sim}1.3$ $\,R_{\rm\earth}$) cores, while more massive planets hold onto their H/He envelopes (e.g., Owen & Wu 2017; Jin & Mordasini 2017). If this were universally true, then the larger of these planets should be more or less the same mass as their smaller counterparts, as a H/He envelope will contribute a significant fraction of a planet’s radius but very little mass. However, some planets with radii between $2$ and $3$ $\,R_{\rm\earth}$ appear too dense for this scenario (see, e.g., the recent TESS discovery of HD 21749 b; Dragomir et al. 2019). One explanation is that these denser sub-Neptune planets correspond to those in the large-core tail of the distribution. Another is that planet formation proceeds hierarchically, first accreting a rocky core, followed by CNO (e.g., water), and finally H/He, suggesting that planets of $\mathord{\sim}2.4$ $\,R_{\rm\earth}$ correspond to “water worlds”—planets with a high mean-molecular-weight envelope (Zeng et al. 2017, 2018). However, this alternative would not explain the low-mass, large planets more consistent with an envelope of H/He. If both modes of planet formation operate, observation can constrain their relative occurrence. Systems like TOI 125, which contains three sub-Neptune planets of similar size and with a range of orbital periods, provide a good opportunity to measure the densities of these planets under controlled conditions. Having been subjected to the same stellar environment, conclusions based on the relative properties of the planets orbiting TOI 125 are less affected by assumptions about stellar evolution and the history of stellar irradiation. As such, spectroscopic follow-up of TOI 125 may provide insight into planet formation and evolution.

TOI 125 is well suited to precise radial velocity measurements to determine the mass of its planets. The star is a bright ($V=11.0$; $T=10.1$), slowly rotating (${v\sin{i_{\star}}}<2$  km s*-1*), late-G dwarf with very little photometric variation ($\sigma_{\rm phot}$). The Chen & Kipping (2017) planetary mass-radius relationship predicts masses of $8.5$, $8.6$, and $9.5\,\,M_{\rm\earth}$, corresponding to RV semi-amplitudes $3.7$, $2.9$, and $2.5$  m s*-1*. Given an instrumental precision of $\mathord{\sim}1$  m s*-1* for facilities such as HARPS and PFS, we expect that all three planets should have detectable RV signals.

We can estimate the time requirements to characterise the TOI 125 system with the HARPS spectrograph using the RVFC tool developed by Cloutier et al. (2018). RV noise sources are estimated as a combination of the instrument noise floor ($0.$5  m s*-1*), the photon noise ($2.51$  m s*-1* for $30$-minute exposures), stellar activity ($0.5$  m s*-1* for a worst case $v\sin{i_{\star}}$ of $2$  km s*-1*) and the RV rms caused by additional unseen planets (typically $0.4$  m s*-1* in this case). Details of how these noise sources are generated from the known stellar parameters, including Gaussian Process (GP) trials to simulate the stellar activity, are given in Cloutier et al. (2018). A complication for this system is its known multiplicity. Additional planets which are perfectly modelled do not impact the characterization of a planet, but given no model is perfect some additional rms will be present. We use the unseen planet RV rms estimate from RVFC as a zero-order guess at this contribution to the noise budget. We take the longest period planet as the driver of the necessary observations, implicitly assuming that observations sample well the orbital phase curve of that and each interior planet, and that all planets are on circular orbits. Although photon noise dominates for this apparently low activity star, the effect of stellar activity can be large depending on the rotation period of the star, and in particular whether it is near a harmonic of the planet orbital periods. We present several representative cases, each calculated using the given stellar parameters, estimated planetary masses and $10$ GP trials to estimate stellar activity. To characterize the system with a $5$-$\sigma$ detection of the semi-amplitude for the outer planet (TOI 125.03) RVFC predicts $68\pm 9$ RV observations for a stellar rotation period of $25$ days, rising to $141\pm 37$ observations for a difficult case rotation period of $40$ days, double the orbital period of candidate .03. For the case of a $20$-day stellar rotation period characterizing candidate .03 becomes untenable, with RVFC predicting $45\pm 5$ observations to characterize only planets b and c.

If real, the low SNR events TOI 125.04 and TOI 125.05 may complicate the RV detection of the other three planets, but they are also predicted to have detectable RV signals. With a timescale very different from the other planets and a predicted semi-amplitude of $2.2$  m s*-1*, TOI 125.04 could be detected with a dedicated high-cadence RV campaign, which would ultimately benefit the detection of the outer signals, as it would otherwise enter as an additional source of noise. TOI 125.05 is the lowest SNR candidate and least likely to prove real, but may also produce a detectable RV signal. Its predicted semi-amplitude given its derived size ($4.2$ $\,R_{\rm\earth}$) is $\mathord{\sim}4.8$  m s*-1*. On the other hand, such a large planet is not typically seen in tightly packed systems. Therefore, we expect that if it is real, the planet is more likely to reside at the smaller end of the fit posteriors (see Figure 4), with an RV signal on the order of $1$–$2$  m s*-1*.

6.2 TOI 125 among the Kepler multis

Nearly $2000$ of the $4723$ Kepler Objects of Interest (KOIs) reside in multi-planet systems, and their orbital architectures provide clues to their formation and evolution: the typical mutual inclination of short-period systems can be derived from the number of planets per system (e.g., Lissauer et al. 2011b; Fang & Margot 2012; Ballard & Johnson 2016), and from the ratio of transit durations within each system (e.g., Fang & Margot 2012; Fabrycky et al. 2014), and inform the dynamical histories of planetary systems; the sizes and orbital spacing of neighboring planets relay information about formation and physical evolution (e.g., Weiss & Marcy 2014); the assembly of planets from planetesimals in the inner region of the protoplanetary disk can be examined through the lens of the present-day properties of short-period planets (e.g., Lee & Chiang 2017).

One striking feature of the population of Kepler multi-planet systems is the distribution of period ratios near first-order mean motion resonances. As discussed in Section 5, there is an underdensity of planet pairs just interior to first-order resonances—particularly the 2:1 resonance—and an excess of systems just exterior to resonance. We present an updated histogram of period ratios for Kepler systems in Figure 10 (see Steffen et al. 2013; Fabrycky et al. 2014, for a broader discussion of these data), and we note that TOI 125 b and c have a period ratio that falls right in the gap interior to the 2:1 resonance. We have explored possible causes for this in Section 5, and here we compare TOI 125 to the small handful of other systems in or near this gap, seen in the inset of Figure 10.

None of the five other systems interior to, but within $2\%$ of 2:1, is quite like TOI 125, which is larger and/or shorter period than the others. Kepler-176 d (KOI-520.03) is the most similar in size to TOI 125 b and c, but it, along with Kepler-334 (KOI-1909), is longer period (weeks, rather than days; Rowe et al. 2014). Kepler-271 (KOI-1151) and KOI-4504 have similar periods to TOI 125, but the planets are much smaller. One other system (KOI-1681) has an interesting architecture, with one small planet and one hot Jupiter (similar to WASP-47; Becker et al. 2015), but it also has a nearby stellar companion, a third signal that is likely a false positive, and a fourth signal corresponding to another giant planet candidate. If real, it would be unlike any other system that we know. TOI 125 is thus a unique system even within its sparsely populated region of parameter space, and worthy of additional study.

6.3 The Candidate USP Planet TOI 125.04

While TOI 125.04 was only detected with ${\rm SNR}\mathord{\sim}5.2$, the presence of three high-SNR transit signals in the system makes it more likely that the signal is real compared to an isolated signal of similar strength. Moreover, the architecture of the TOI 125 planets would match that of other known USP systems, both in semi-major axis and in mutual inclination. Dai et al. (2018) find that Kepler and K2 planets in multi-planet systems tend to exhibit large mutual inclinations when the innermost planet has a small semi-major axis ($a/R_{*}<5)$, and that among these USP planets, the systems with the largest mutual inclinations also have large period ratios (see, e.g., K2-266, with an extreme mutual inclination of $\mathord{\sim}12$ degrees Rodriguez et al. 2018b). TOI 125.04 orbits at ($a/R_{*}\mathord{\sim}3.1$) and with a projected mutual inclination of $\mathord{\sim}16$ degrees. The period ratio between TOI 125 b and TOI 125.04 is $8.8$, similar to the other such misaligned USP planets. TOI 125.04 would be the USP planet with the largest known mutual inclination, but as discussed above, TOI 125 b is larger than most inner planets in packed systems. We speculate that in the framework suggested by Dai et al. (2018), one would expect the interaction between the two planets that leads to the inclined USP planet to produce a more extreme outcome when the adjacent planet is more massive than usual.

6.4 Dynamical Results

The EXOFASTv2 transit fit allowed a relatively large range in eccentricities for TOI 125 b, c, and d, but the largest of the allowed values can be excluded due to dynamical stability constraints. Of the entire EXOFASTv2 posteriors, approximately 32% of draws result in integrations that remain dynamically stable for 1 Myr. Preferentially, the stable subset are those with lower eccentricities: eccentricities above 0.25–0.3 are disallowed for each planet (see Figure 5.2).

Using the results of the numerical simulations to inform the resonant behavior paints a largely non-resonant picture of the posteriors. Of the stable subset of integrations, roughly 86% exhibited nearly exclusively non-resonant behavior for planets b and c. In these cases, the median period ratio libration has a (min-to-max) amplitude of $\Delta\delta\approx 0.01$, surrounding a median $\delta$ value of $\delta\approx-0.032$ (compared to the measured current-day value of $\delta\approx-0.033$). From these simulations, it is reasonable that the observed three planet system can remain dynamically stable in the observed orbits for a relatively large fraction of the EXOFASTv2 posteriors. In comparison, among the $\sim 12\%$ of integrations that reside in or nod in/out of the 2:1 resonance for some or all of of the 1 Myr, while in resonance the size of $\Delta\delta$ depends on the libration width of the resonant angle. For nodding, this value is as large as $\approx 0.08$, around a median $\delta$ value of $\delta\approx 0.001$. In true resonance, the typical $\Delta\delta$ values range between 0.02 and 0.04, as the median $\delta$ value becomes closer to $\delta\approx 0$. It is possible that TOI 125 b and c previously resided in the 2:1 resonance and naturally lost the resonance at some point, becoming trapped in a dynamically stable but non-resonant orbit Although the errors on the currently measured orbit are too large to conclusively determine the current-day resonance behavior of TOI 125 b and c, the simulations suggest they are most likely not currently in resonance, but have eccentricities near the lower end of the measured posteriors.

7 Summary

In this paper, we have presented the TESS discovery of the TOI 125 multi-planet system, and we fit a global model to the TESS data, spectroscopic stellar parameters, literature photometry, and the Gaia parallax in conjunction with stellar models to characterize the planet candidates. We then statistically validated the planetary nature of TOI 125 b and c using vespa with the aid of archival imaging and our photometric, spectroscopic, and high resolution imaging observations. We demonstrated that the system is likely amenable to mass determination via both TTV and precise RV follow-up, and that the planets are worthy of such additional study. The three strongest transit signals are caused by planets with radii $2.8$–$2.9$ $\,R_{\rm\earth}$, a class of planet not seen in the Solar System but abundant in the galaxy. These planets have been proposed as the progenitors (via photoevaporation) of the terrestrial planets commonly found in short periods around nearby stars, and studying three of them in the controlled environment of the same host star can help illuminate the formation and evolution processes at play. The candidate terrestrial USP planet, with an orbital period less than $13$ hours and a mutual inclination of $16$ degrees with the other planets, is an extreme example of the trend toward such architectures among other known USP planets in multiple systems, and may be the end result of dynamical interaction with its much larger sub-Neptune neighbors. Finally, the period ratio between planets b and c is very near, but just interior to, a 2:1 commensurability, which is quite unusual compared to known Kepler systems. While one possible explanation is that the system is in—and librating about—the 2:1 resonance, our dynamical analysis suggests that it is unlikely that the system is currently in true resonance. The discovery of the TOI 125 system demonstrates that TESS continues in its early days to deliver on its promise to identify rare systems of small planets amenable to follow-up observations and detailed characterization.

JCB is supported by the NSF Graduate Research Fellowship Grant No. DGE 1256260 and by the Leinweber Center for Theoretical Physics. Work performed by JER was supported by the Harvard Future Faculty Leaders Postdoctoral fellowship. DJA gratefully acknowledges support from the STFC via an Ernest Rutherford Fellowship (ST/R00384X/1). AV’s contribution to this work was performed under contract with the California Institute of Technology (Caltech)/Jet Propulsion Laboratory (JPL) funded by NASA through the Sagan Fellowship Program executed by the NASA Exoplanet Science Institute. Work by JNW was supported by the Heising-Simons Foundation. DD acknowledges support for this work provided by NASA through Hubble Fellowship grant HST-HF2-51372.001-A awarded by the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., for NASA, under contract NAS5-26555. AG is supported by the Ida M. Green Fellowship. M.T. is supported by MEXT/JSPS KAKENHI grant Nos. 18H05442, 15H02063, and 22000005. This work is partly supported by JSPS KAKENHI Grant Numbers JP15H02063, JP18H01265, JP18H05439, JP18H05442, and JST PRESTO Grant Number JPMJPR1775. We acknowledge the support provided by the Polish National Science Center through grants 2016/21/B/ST9/01613 and 2017/27/B/ST9/02727. We thank Zach Hartman, Dan Nusdeo, and Jen Winters for help with the Gemini-South observations. We thank Akihiko Fukui, Nobuhiko Kusakabe, Kumiko Morihana, Tetsuya Nagata, Takahiro Nagayama, Taku Nishiumi, and the staff of SAAO for their kind support for IRSF SIRIUS observations and analyses. Funding for the TESS mission is provided by NASA’s Science Mission directorate. We acknowledge the use of public TESS Alert data from pipelines at the TESS Science Office and at the TESS Science Processing Operations Center. This paper includes data collected by the TESS mission, which are publicly available from the Mikulski Archive for Space Telescopes (MAST). This work used the Extreme Science and Engineering Discovery Environment (XSEDE), which is supported by National Science Foundation grant number ACI-1053575. This research was done using resources provided by the Open Science Grid, which is supported by the National Science Foundation and the U.S. Department of Energy’s Office of Science, through allocation TG-AST150033. The Swiss Euler telescope is supported by the Swiss National Science Foundation. This paper includes observations obtained under ESO/VLT program 0102.C-0503(A), and Gemini program GS-2018B-LP-101. Gemini Observatory is operated by the Association of Universities for Research in Astronomy, Inc., under a cooperative agreement with the NSF on behalf of the Gemini partnership: the National Science Foundation (United States), National Research Council (Canada), CONICYT (Chile), Ministerio de Ciencia, Tecnología e Innovación Productiva (Argentina), Ministério da Ciência, Tecnologia e Inovação (Brazil), and Korea Astronomy and Space Science Institute (Republic of Korea). This work makes use of observations from the LCOGT network. The IRSF project is a collaboration between Nagoya University and the South African Astronomical Observatory (SAAO) supported by the Grants-in-Aid for Scientific Research on Priority Areas (A) (Nos. 10147207 and 10147214) and Optical & Near-Infrared Astronomy Inter-University Cooperation Program, from the Ministry of Education, Culture, Sports, Science and Technology (MEXT) of Japan and the National Research Foundation (NRF) of South Africa. Based in part on data collected with Solaris network of telescopes of the Nicolaus Copernicus Astronomical Center of the Polish Academy of Sciences. The research leading to these results has received funding from the European Research Council under the European Union’s Seventh Framework Programme (FP/2007-2013) ERC Grant Agreement 336480, from the ARC grant for Concerted Research Actions, financed by the Wallonia-Brussels Federation, and from a research grant from the Balzan Prize Foundation. TRAPPIST is funded by the Belgian Fund for Scientific Research (Fond National de la Recherche Scientifique, FNRS) under the grant FRFC 2.5.594.09.F, with the participation of the Swiss National Science Fundation (SNF). M.G. and E.J. are FNRS Senior Research Associates. This work has made use of data from the European Space Agency (ESA) mission Gaia (https://www.cosmos.esa.int/gaia), processed by the Gaia Data Processing and Analysis Consortium (DPAC, https://www.cosmos.esa.int/web/gaia/dpac/consortium). Funding for the DPAC has been provided by national institutions, in particular the institutions participating in the Gaia Multilateral Agreement. This research has made use of NASA’s Astrophysics Data System. This research has made use of the VizieR catalogue access tool, CDS, Strasbourg, France. The original description of the VizieR service was published in A&AS 143, 23. This research has made use of the SIMBAD database, operated at CDS, Strasbourg, France. This research has made use of the NASA Exoplanet Archive and the Exoplanet Follow-up Observation Program website, which are operated by the California Institute of Technology, under contract with the National Aeronautics and Space Administration under the Exoplanet Exploration Program.