A resonant pair of warm giant planets revealed by TESS

ABSTRACT

We present the discovery of a pair of transiting giant planets using four sectors of TESS photometry. TOI-216 is a 0.87 M_⊙ dwarf orbited by two transiters with radii of 8.2 and 11.3 R_⊕, and periods of 17.01 and 34.57 d, respectively. Anticorrelated TTVs are clearly evident indicating that the transiters orbit the same star and interact via a near 2:1 mean motion resonance. By fitting the TTVs with a dynamical model, we infer masses of |$30_{-14}^{+20}$| and |$200_{-100}^{+170}$| M_⊕, establishing that the objects are planetary in nature and have likely sub-Kronian and Kronian densities. TOI-216 lies close to the southern ecliptic pole and thus will be observed by TESS throughout the first year, providing an opportunity for continuous dynamical monitoring and considerable refinement of the dynamical masses presented here. TOI-216 closely resembles Kepler-9 in architecture, and we hypothesize that in such systems these Saturn analogues failed to fully open a gap and thus migrated far deeper into the system before becoming trapped into resonance, which would imply that future detections of new analogues may also have sub-Jupiter masses.

1 INTRODUCTION

A small fraction of exoplanets in the cosmos have the correct orbital geometry to transit their star as seen from our home. These transiting planets have been a ‘royal road to success’ in planet discovery (Russell 1948; Winn 2010) yielding thousands of discoveries in recent years (see the NASA Exoplanet Archive; Akeson et al. 2013), despite the fact that they represent just a sliver of the total population. An even more rarefied population is that of transiting planets exhibiting transit timing variations (TTVs; Agol et al. 2005; Holman & Murray 2005; Deck & Agol 2015). For these special worlds, not only can one infer the planetary size from the transit depths, but dynamical modelling of the TTVs can often provide planetary masses too – a fact heavily exploited by Kepler (Holman et al. 2010; Lithwick, Xie & Wu 2012; Nesvorný et al. 2012). In such cases, it is therefore possible to confirm the planetary nature of a system almost exclusively from photometric observations (e.g. Ford et al. 2011; Steffen et al. 2013).

Kepler enjoyed many successes with this strategy, largely enabled by its patience to stare at the same stars for over 4 yr continuously. With TESS, the full-sky nature of the survey means that most parts of the sky are observed for much shorter windows,¹ potentially posing a challenge to dynamical confirmation of planetary candidates. However, TESS does maintain a longer vigil on the ecliptic poles, observing these fields for up to a year continuously (Ricker et al. 2016).

In this work, we describe the discovery of two TESS planets near a 2:1 mean motion resonance (MMR) leading to highly significant TTVs. Thanks to the host star’s fortuitous location near the southern ecliptic pole, TESS can observe the target for most of the first year making TOI-216 an excellent target for monitoring planet–planet interactions. We describe the observations by TESS in Section 2 with attention to detrending, contamination, and stellar properties. In Section 3 we regress light-curve models and TTV models, demonstrating that the system is a pair of planet-mass objects gravitationally interacting with one another. Finally, we discuss the possibilities opened up by this exciting new system in Section 4.

2 OBSERVATIONS

2.1 Identification

TIC 55652896 was observed by TESS in the first four sectors of year one, and indeed is scheduled to be observed in every sector of that year except for sector 10 (2019 March 26–April 22). Falling on camera 4, the target is a relatively rare example of a transiting planetary system caught within the TESS continuous viewing zone (CVZ). With an ecliptic latitude of −82.476408°, we highlight that the target also lands close to JWST’s planned CVZ and would be observable for ≳ 260 d per year of the mission.

A TESS alert was issued on 2018 November 30 for two candidate transiting planets associated with TIC 55652896 using sectors 1 and 2, dubbed TOI-216.01 and TOI-216.02. With periods of ∼17.1 (TOI-216.02) and ∼34.5 d (TOI-216.01), the outer candidate was only seen to transit twice during this time (once per sector).

Among the 300 + TOIs identified at the time of writing, this pair stood out as particularly interesting because the planetary candidates lie near a 2:1 period commensurability. If the objects were orbiting the same star, and gravitationally interacting, then it may be possible to confirm the planetary nature of the pair without any ground-based follow-up (see e.g. Steffen et al. 2013). For this reason, we decided to further study this system.

2.2 Stellar properties

TOI-216 has an apparent magnitude of 11.5 in the TESS bandpass. From the TESS Input Catalogue (TIC) version 7 (Stassun et al. 2018), catalogue survey spectroscopy of the star constrains T_eff = (5026 ± 125) K, [M/H] = 0.32 ± 0.10, and log (g) = 4.66 ± 0.20. These properties are used by TIC to infer M_⋆ = (0.879 ± 0.073) M_⊙ and R_⋆ = (0.715 ± 0.166) R_⊙.

We also queried the star within Gaia DR2 (Gaia Collaboration 2018) and find a parallax measurement of 5.591 ± 0.028 mas (GAIA DR2 4664811297844004352). Rather than use the TIC-7 summary statistics for stellar mass and radius, we would prefer to work with posterior samples – as well as include the Gaia parallax – and so we elected to perform our own Bayesian isochrone fitting. To this end, we used the isochrones package by T. Morton using the previously listed constraints on V, T_eff, [M/H], log g, and parallax with the Dartmouth stellar evolution models.

The inputs to our fits (the ‘star.ini’ file) are given in the top panel of Table 1, and the derived parameters of interest in the lower panel. As expected, our results closely agree with those listed in the TIC, although our inference is more precise as a result of using the GAIA parallax which was not available when TIC-7 was compiled.

2.3 Contamination

With a pixel size of 21 arcsec, there is a greater chance of crowding with TESS than Kepler. The aperture used in each sector varies slightly but is approximately 4 by 3 pixels and thus sources out to 84 arcsec can contaminate the aperture.

Fortunately, there are no comparably bright stars that lie within this region. The nearest star listed in the TIC to our target is TIC 55652894, separated by 48.2 arcsec but far fainter with an apparent TESS magnitude of 16.3 (1.2 per cent the brightness level). Gaia DR2 (Gaia Collaboration 2018) reports 46 stars within 84 arcsec, with G-band magnitudes from 17.1 to 21.1 (TOI-216 is 12.2). Together, these sources could maximally dilute the target by 5.6 per cent in G, although the true value will be less due to colour correction to the redder TESS bandpass, location of the sources, and the finite PSF widths.

These contaminating source cumulatively lead to a small amount of dilution of TOI-216, which is estimated within the TESS aperture to be 0.33 per cent, 0.46 per cent, 0.29 per cent, 0.26 per cent, 0.58 per cent, and 0.31 per cent for sectors 1 to 6, respectively, (values taken form the TESS light-curve files) and these are included in our later light-curve fits using the prescription of Kipping & Tinetti (2010).

We also highlight that an unresolved companion may reveal itself through centroid shifts during the moments of transit, but all centroid shifts for this star are below 1 σ as reported by the TESS sector 1–3 cumulative Data Validation (DV) report.

2.4 Light-curve detrending

We downloaded sectors 1–6 short-cadence (2 min) data for TIC 55652896 and work with the PDC (Pre-search Data Conditioning) product in what follows (Jenkins et al. 2017). Any bad data flags were removed, and outliers filtered with 5 σ clipping against a 20-point moving median. Using the ephemeris and duration for TOI-216.01 and TOI-216.02, we de-weight the in-transit points for the purposes of detrending. The PDC light curve for all six sectors is shown in Fig. 1.

We next detrend each transit epoch of each planet independently, using four different algorithms following Teachey & Kipping (2018) – CoFiAM, PolyAM, local polynomials, and a Gaussian Process using a squared exponential kernel. The four light curves are then combined into a single time series – a method marginalized light curve – by taking the median at each time stamp and inflating the formal uncertainties by adding the inter-method standard deviation in quadrature. We direct the reader to Teachey & Kipping (2018) for a detailed description of the four algorithms as well as the method marginalization process. The resulting light curves from all four methods, as well as the method marginalized light curves, are made available at this URL.

We find that the inter-method standard deviation is many times smaller than the formal uncertainties, indicating a highly stable detrending. The median formal uncertainty is 2430 ppm but the median inter-method standard deviation is 16 times smaller at 150 ppm. After adding this extra component in quadrature to the formal uncertainties, the error increases by just 0.2 per cent.

3 ANALYSIS

3.1 Light-curve model

We initially built light-curve models for the system which treated each planetary candidate as orbiting an independent star. It became immediately clear that the transiters displayed strong TTVs, as can be seen by simple inspection of Fig. 2. More detailed analysis of these light curves presented in Section 3.3 reveals strong evidence for anticorrelation – the hallmark of dynamically interacting planets (Steffen et al. 2013). This is only possible if the two transiters are orbiting the same primary, and thus in our final light-curve modelling we decided to treat the objects as sharing a common host star.

This is particularly useful for modelling the inner candidate, TOI-216.02, whose light curve displays a V-shaped morphology consistent with a grazing geometry. Treated as an independent body, V-shaped transits display strong degeneracies between size, impact parameter, limb darkening coefficients, and host star density (Carter et al. 2008). Since the anticorrelated TTVs imply a common host star, and the other transiter is non-grazing, the conditional relationship greatly aids in the inference of a unique light-curve solution for TOI-216.02.

Our light-curve model is that of the classic Mandel & Agol (2002) quadratic limb darkening code, which is oversampled by a factor of 5 to correct for the slight distortion of finite integration time via the prescription of Kipping (2010b). The quadratic limb darkening coefficients are reparametrized to q₁ and q₂ following Kipping (2010a).

Our later fits (see Section 3.2) reveal that our light-curve-derived density is measured to a precision of 6 per cent and thus from equation (1) one may show that e < 0.14 should be expected to lead to a less than 1 σ systematic error in the inferred density. As an apparently fairly compact, multiple planet system, we consider this assumption is reasonable on the grounds of dynamical stability, and indeed our later TTV fits favour low eccentricities (see Section 3.3). Nevertheless, we choose not to use the light-curve-derived density in any attempt to refine the isochrone modelling from Section 2.2.

3.2 Light-curve fits

Our light-curve model has a total of 17 free parameters: two ratio-of-radii (p_inner & p_outer), two impact parameters (b_inner & b_outer), two limb darkening coefficients (q₁ & q₂), a mean stellar density (ρ_⋆), six times of transit minimum for the inner transiter (τ_{inner, i}) and four times of transit minimum for the outer transiter (τ_{outer, i}). We assume uniform priors for all parameters except for ρ_⋆ for which we adopt a broad loguniform prior. All of the priors are listed in Table 2

Fits were conducted using multinest (Feroz & Hobson 2008; Feroz, Hobson & Bridges 2009) with 4000 live points,² an evidence tolerance of 1.0 and in non-constant efficiency mode. The maximum a posteriori light-curve solution is plotted in Fig. 2. We make the full posterior samples available at this URL but list the credible intervals on the 10 transit times in Table 3 and the other seven global parameters in Table 5.

Considering the seven global parameters, there are two noteworthy conclusions that can be drawn from the results. The first is that the impact parameter of the inner planet is unusually high at |$b_{\mathrm{inner}}=0.957_{-0.022}^{+0.047}$|⁠. Given that the ratio-of-radii is measured to be |$p_{\mathrm{inner}}=0.089_{-0.012}^{+0.032}$|⁠, then we have b_inner > 1 − p_inner and thus this is a definitively grazing transit. Such transits are rare and have been hypothesized to be powerful probes of nodal variations (Kipping 2009), and thus TOI-216.02 should be carefully monitored in the future for such changes.

Second, the a posteriori mean stellar density is found to be |$2380_{-140}^{+100}$| kg m⁻³. We remind the reader that this was using a loguniform prior and thus was inferred agnostically. The only assumption in the model is that both transiters orbit the same star, which is established from the anticorrelated TTVs, and that the eccentricities are small (e ≲ 0.14), which is reasonable given the system’s compactness and multiplicity for orbital stability. This density is consistent with the independently derived value from our earlier isochrone analysis (Section 2.2), which yielded |$2090_{-300}^{+270}$| kg m⁻³ and thus adds further credence to the hypothesis that both transiters are orbiting the target star TIC 55652896, rather than some unresolved companion.

3.3 Transit timing variations

We plot the TTVs in Fig. 3, where one can clearly see the strong case for anticorrelation mentioned earlier in this work. This establishes that the transiters orbit the same star (Steffen et al. 2013), although this point alone does not establish the planetary nature of the two transiters. Critically, their masses could potentially be consistent with a brown dwarf or a late-type star, especially for TOI-216.01 whose radius is similar to Jupiter.

The TTVs may be modelled by considering two masses orbiting a primary with an N-body integrator. Dynamical analysis of the observed transit times was performed with a symplectic N-body integrator code described in Nesvorný et al. (2012). The code was instructed to simultaneously fit all transit times of both transiters, using multinest to perform the regression. The integration time-step was initially set to 0.6 d, but we also repeated for the final fits moving to a higher resolution of 0.3 d to ensure we recover consistent results – which indeed we do. We also note that the results are consistent to that derived when considering just sectors 1–4, which was done in a previous draft of this paper.

Our dynamical model has 14 parameters: mass ratios M_inner/M_⋆ and M_outer/M_⋆, orbital periods P_inner and P_outer, eccentricities e_inner and e_outer, longitudes of periapsis ϖ_inner and ϖ_outer, impact parameters b_inner and b_outer, difference in nodal longitudes Ω_outer − Ω_inner, stellar density ρ_⋆, and reference epochs τ_inner,ref and τ_outer,ref between a reference time and the first observed transit of each planet. All orbital elements are given at the reference time 2457 000 BJD_UTC (TBJD).

We used uniform priors for all parameters except for b_inner, b_outer, and ρ_⋆, since our earlier light-curve fits provide strong constraints which can be leveraged here. Since multinest requires simple parametric forms of the priors for the purposes of inverse transform sampling, we approximated the marginal posteriors from our earlier fits such that b_outer is uniform between 0 and 0.4, b_inner is a Gaussian centred on 0.95 with 0.025 standard deviation, and ρ_⋆ is a Weibull prior with shape parameters 22.7 and 2425.1. These priors are listed in Table 4.

The fits converged to a unique solution and the joint posteriors are depicted in Fig. 4 for reference. The maximum a posteriori solution is plotted in Fig. 3, which shows how the model is able to fully describe the observed deviations.

Combining the derived mass ratios with the stellar mass derived earlier (see Section 2.2) allows us to measure that |$M_{\mathrm{inner}} = 30_{-14}^{+20}$| M_⊕ and |$M_{\mathrm{outer}} = 200_{-100}^{+170}$| M_⊕, which establishes that the masses are far below the deuterium burning limit and these objects may be classified as ‘planets’. Accordingly, in what follows, we refer to TOI-216.02 as TOI-216b (the inner planet), and TOI-216.01 as TOI-216c (the outer planet).

3.4 Final parameters

To complete our analysis, we combine the fundamental stellar parameters derived earlier (see Section 2.2) with the relative radii (from Section 3.2) and relative masses (from Section 3.3) to calculate physical properties for both planets. Our final planet properties are listed in Table 5.

4 DISCUSSION

We have demonstrated that the TESS planetary candidates TOI-216.01 and .02 must orbit the same primary star given their anticorrelated TTVs. The light-curve-derived stellar density found by fitting both signals yields a value almost precisely equal to the target star’s density from an isochrone analysis, establishing that the objects indeed orbit the target rather than a contaminant. Finally, we have regressed an N-body dynamical model to the observed TTVs to demonstrate that the masses of each body are far below the deuterium burning limit making these bona fide ‘planets’.

The TOI-216 planetary system displays some close similarities to the Kepler-9 system (Holman et al. 2010), but is 1.6 magnitudes brighter in V. In both cases, one finds low-density gas giants in a 2:1 MMR orbiting a Sun-like star at similar periods (∼20 and 40 d). To a lesser degree, the system also resembles KOI-872 (Nesvorný et al. 2012). In both of these cases, the MMR pair of planets are accompanied by a short-period super-Earth and thus it is natural to wonder if perhaps TOI-216 may also be accompanied by a small and currently unresolved terrestrial planet. We ran a box-least-squares search (Kovács, Zucker & Mazeh 2002) for such a signal but find no significant peaks with the available TESS data.

The TTVs of TOI-216 are characterized by a superperiod as the longitude of conjunctions circulates with a time-scale of |$\mathcal {O}[10^3]$| d. Although TOI-216 will be monitored throughout the first year of TESS observations, the superperiod looks likely to exceed this baseline and thus continuous monitoring from the ground in 2020 would greatly benefit the determination of precise orbital elements. With the limited phase coverage available at the time of writing, the masses quoted in this work will surely be refined considerably in the future.

The resonance between the gas giants is consistent with dissipative processes in disc–planet interaction during their presumably inward migration from beyond the snow line (Crida, Sándor & Kley 2008; Havel et al. 2011; Cimerman, Kley & Kuiper 2018). In the Grand Tack hypothesis of the Solar system (Hansen 2009; Walsh et al. 2011), Jupiter is thought to have opened up a gap, migrating slower (type II) than Saturn (type I), which likely failed to fully open a gap. This enabled Saturn to catch up to Jupiter, trapping the pair in resonance when Jupiter was at ∼1.5 au, which reversed subsequent migration. In the case of both Kepler-9 and TOI-216, the gas giants have maximum a posteriori masses similar to Saturn, and thus may have failed to have opened up full gaps, causing them to type I migrate far deeper. As a larger sample of such systems is found in the future, it will be interesting to test if these giants tend to have sub-Jupiter masses in a statistically significant manner, as would be expected under this hypothesis.

Disc migration is likely a more favourable scenario for the survival of satellite systems around these giants than planet–planet scattering (Gong et al. 2013). Further, the fact that TOI-216 is relatively bright means that follow-up with larger facilities would be well suited to make a search for a satellite system. Similarly, the low-density (high scale height), deep transits, and reasonably bright target star would make TOI-216 a potential target for atmospheric characterization of gas giants in a cooler regime to their hot-Jupiter counterparts. Further, like Kepler-9, the two planets provide an opportunity for differential transit spectroscopy alleviating systematic effects.

ACKNOWLEDGEMENTS

DMK is supported by the Alfred P. Sloan Foundation Fellowship. DN’s work was supported by the NASA Exoplanet Research Program (XRP). AT and TJ are supported through the NSF Graduate Research Fellowship (DGE 16-44869).

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 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 SIMBAD data base, 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.

Facilities: TESS

Software:multinest (Feroz & Hobson 2008; Feroz et al. 2009), isochrones (Morton 2015), mercury6 (Chambers 1999), corner.py (Foreman-Mackey 2016).

Footnotes

REFERENCES