EDP Sciences logo
Subscriber Authentication Point
Search
Menu
Home All issues Volume 669 (January 2023) A&A, 669 (2023) A109 Full HTML
Open Access
Issue
A&A
Volume 669, January 2023
Article Number A109
Number of page(s) 23
Section Planets and planetary systems
DOI https://doi.org/10.1051/0004-6361/202243879
Published online 24 January 2023

© The Authors 2023

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

This article is published in open access under the Subscribe-to-Open model. Subscribe to A&A to support open access publication.

1 Introduction

So far, 860 multi-planetary systems (with two or more planets) around stars other than the Sun have been found, hosting a total of 2050 planets according to the NASA Exoplanet archive (Akeson et al. 2013). In contrast to our Solar System, the components of most of these known multi-planet systems reside in close-in compact configurations (e.g., Gillon et al. 2016), mainly biased by the observational techniques. This is especially relevant in the case of late-type stars, with effective temperatures below 4600 K, for which 90% of multi-planetary systems have all their detected planets interior to the snow line (the distance from the parent star within which volatiles can condense into ice grains), and 91% of the known long-period planets are massive and isolated (Lillo-Box et al. 2022). Only a small number of systems in this stellar regime show components on both sides of the snow line (e.g., WASP-107 by Anderson et al. 2017; Piaulet et al. 2021, or GJ 1148 by Trifonov et al. 2018). This is a clear fingerprint of the planet migration and formation history. Therefore, finding new systems located at both sides of the snow line is an important step towards a comprehensive understanding of planetary histories.

Exoplanet exploration has also revealed different exoplanet populations, including zones in the parameter space where there is a clear dearth of extrasolar planets. One of the most intriguing is the hot-Neptune desert (Szabó & Kiss 2011; Owen & Lai 2018), a deficit of icy to sub-Jovian planets (Rp ~ 2–10 R⊕) at close-in orbits (P < 5 days). Several theories have been proposed to explain the emptiness of this region, such as atmospheric mass loss of low-mass gaseous giants unable to retain their atmospheres due to the strong stellar irradiation suffered along their inward migration (Lopez & Fortney 2014; Owen & Wu 2013); the core accretion scenario forming gas giants more preferentially than Neptune-like planets (e.g., Ida & Lin 2008); or a faster migration of the lowest mass end of gas giants, preventing them from accreting large amounts of gas (Flock et al. 2019).

Interestingly, both the Kepler (Borucki et al. 2010) and TESS (Transiting Exoplanet Survey Satellite, Ricker et al. 2014) missions, assisted by ground-based follow-up observations, have sparsely populated this desert with planets displaying unusual properties. Two clear examples are TOI-849 b (Armstrong et al. 2020) and LTT 9779 b (Jenkins et al. 2020), two ultra-short-period planets in the middle of the hot-Neptune desert. The advantage of the planets detected by TESS in this regime is that they usually orbit stars much brighter than those detected by Kepler. This, combined with the transiting nature of these detections and the precise masses achievable due to their proximity to their host star, opens the possibility to study their interior properties (through internal structure analysis, e.g., Delrez et al. 2021) and composition with the James Webb Space Telescope (JWST, Gardner et al. 2006).

In this paper, we present both confirmation and a characterization of the planetary system TOI-969 (also known as TIC 280437559, ‘TYC 183-755-1, Gaia DR3 3087206553745290624), which is composed of at least two very distinct companions. The architecture of this system and the properties of the individual components could be a suitable test bench for studying theories of planet formation and migration in the cool regime of solar-like stars. The inner component, a short-period mini-Neptune, was detected by the TESS mission. The outer component is an eccentric massive body in the boundary between the planet and brown dwarf domains, and here we report, for the first time, its orbit beyond the snow line, finding a long orbital period. In Sect. 2, we describe the observations from different instruments and the techniques we used to characterize the system. In Sects. 3 and 4, we use these observations to constrain the stellar parameters and to unveil the properties of the planetary system, respectively. We discuss the results of this analysis in Sect. 5, and finish in Sect. 6 by providing our final conclusions.

2 Observations

2.1 High-precision photometry

2.1.1 TESS

Observations of TOI-969 are available in sectors 7 (from Cycle 1) and 34 (Cycle 3), in both cases with camera #1. In the case of Sector 7, the cadence of the time series corresponds to 1800s while in the case of Sector 34 a cadence of 600s is available. We downloaded the light-curve files from the MAST archive1 corresponding to the data extraction by the TESS-SPOC pipeline (Caldwell et al. 2020), which uses the same codebase as the Science Processing Operations Center (SPOC; Jenkins et al. 2016). We use the Presearch Data Conditioned Simple Aperture Photometry (PDCSAP; Stumpe et al. 2012, 2014; Smith et al. 2012) flux and remove any data points marked as outliers or having non-numeric values.

We first inspect the target pixel files (TPFs) for each sector using the tpfplotter algorithm2 (Aller et al. 2020) to check for possible contaminant sources within the pipeline aperture. In Fig. 1 we show the average TPF for each sector. As shown, both sectors include three additional sources in the aperture. However, their magnitude contrasts in the Gaia passband are 6.6 mag (source labeled as #2 in Fig. 1), 6.1 mag (source #3), and 7.6 mag (source #4). Overall, the contamination coming from these three sources is well below 0.1%. Consequently, we assume the flux as originating only from TOI-969.

The photometric time series for both sectors are presented in Fig. 2 and Table A.1, together with their corresponding Lomb-Scargle periodograms. They show clear large photometric variations that can be attributed to the stellar activity, with the main periodicities around 8 and 12 days. These might be the harmonics of the possible rotation period at 24 days, as we see in the spectroscopic activity indicators (see Sect. 2.3). In addition, the light curve clearly shows the dips corresponding to the transits of the 1.8-day period planet candidate alerted by the TESS Science Office (see vertical line marks in Fig. 2).

Based on this alert, we observed the transit of the planet candidate from different ground-based facilities.

thumbnail Fig. 1

Target pixel file corresponding to TOI-969 observations by the TESS mission in sectors 7 (upper panel) and 34 (lower panel). In both panels, the sources from Gaia DR3 are shown as red circles, with their size corresponding to the magnitude contrast against TOI-969 (marked with the label “1” and a white cross). The aperture used by the TESS-SPOC pipeline is shown as a shaded red region in each panel. The figures are computed using the tpfplotter algorithm by Aller et al. (2020).
thumbnail Fig. 2

TESS photometric time series from sectors #7 (left panels) and #34 (right panels). Upper panels display the PDCSAP flux and a simple median filter computed trend (blue solid line). The detrended light curves are shown in the lower panels and the transits from TOI-969 b are marked by vertical red lines at the bottom of each panel.

2.1.2 LCOGT-MuSCAT3

A full transit of TOI-969 b was observed simultaneously in Sloan g′, r′, i′, and Pan-STARRS z-short bands on UTC 2021 January 8 using the Las Cumbres Observatory Global Telescope (LCOGT; Brown et al. 2013) 2 m Faulkes Telescope North at Haleakala Observatory on Maui, Hawai’i. The telescope is equipped with the MuSCAT3 multi-band imager (Narita et al. 2020). We used the TESS Transit Finder – a customized version of the Tapir software package (Jensen 2013) – to schedule our transit observations. The images were calibrated using the standard LCOGT BANZAI pipeline, and photometric data were extracted using AstroImageJ. The images were mildly defo-cused and have typical stellar point spread functions with a full width at half maximum (FWHM) of ~ 1′′.6, and circular apertures with radius 4′′ were used to extract the differential photometry. The apertures exclude all of the flux from the nearest Gaia Data Release 3 (DR3) neighbor (TIC 280437561) 9′′ northwest of the target (labeled #2 in Fig. 1). The transit was detected on target in all four filter bands. The final light curves are presented in Table A.2.

2.1.3 MEarth

One transit of TOI-969 b was observed using the MEarth-South telescope array (Irwin et al. 2015) at Cerro Tololo Inter-American Observatory (CTIO), Chile, on 4 February 2020. Six telescopes were used to obtain photometry of the target star, defocused to half flux diameter (HFD) of 6 pixels, equivalent to 5 arcsec given the pixel scale of 0.84 arcsec pix−1. Exposure times were 60 s. One telescope was used to obtain photometry of any contaminating fainter sources at higher angular resolution, with the target saturated, using the same exposure time. Observations were gathered starting at evening twilight and continuing until the target set below an airmass of 2.

Data were reduced using the standard MEarth processing pipeline (e.g., Irwin et al. 2007; Berta et al. 2012) with a photometric extraction aperture of r = 10 pixels (8.4 arcsec). The light curves contain meridian flips during transit. These were taken into account in the data reduction by using separate magnitude zero points for each combination of telescope and meridian side to remove residual flat fielding error. The final extracted photometric time series is presented in Table A.2.

2.1.4 WASP

The field of TOI-969 was observed between 2009 and 2012 by the WASP (Wide Angle Search for Planets) transit-search survey (Pollacco et al. 2006). A total of 14 600 photometric data points were obtained, observing on clear nights when the field was visible with a typical cadence of 15 min. The 200mm f/1.8 Canon lenses were backed by 2k×2k CCDs. TOI-969 is five magnitudes brighter than any other star in the 48arcsec extraction aperture. Table A.3 presents the photometric time series. Following the methods described in Maxted et al. (2011), the Lomb-Scargle periodogram of this long-term dataset reveals a signal at a period of 24.6 ± 0.6 d, together with power at the 12.3-day first harmonic. The modulation has an amplitude of 8 mmag and has a false-alarm probability (FAP) of below 10−3. This signal is similar to that described above in the TESS data, and its origin clearly points to the rotational periodicity of the star (see Sect. 3 for a detailed discussion about the stellar rotation period).

2.2 High-spatial-resolution imaging

The presence of close stellar companions to a star hosting a planet candidate (either physically bound or chance aligned) can create false-positive exoplanet detections under specific configurations (e.g., eclipsing binaries) or provide additional flux leading to an underestimated planetary radius and incorrect exo-planet parameters (see, e.g., Lillo-Box et al. 2012, 2014; Ciardi et al. 2015; Furlan et al. 2017). Besides the large-separation sources described in the previous section through the analysis of the Gaia Early Data Release (EDR3) catalog (Gaia Collaboration 2021), we explored the close-in region (<3 arcsec) around the candidate host star to unveil previously unresolved companions through high-spatial-resolution images at different wavelengths.

TOI-969 was observed on 14 January 2020 using the Zorro speckle instrument on Gemini South (Scott et al. 2021). Zorro provides simultaneous speckle imaging in two bands (562 nm and 832 nm) with output data products including a reconstructed image and robust contrast limits on companion detections (Howell et al. 2011, 2016). The night was clear, had a slight breeze, and good seeing (0.6–1.0 arcsec) during the observations. Figure 3 (left panel) shows our limiting magnitude contrast curves and our 832 nm reconstructed speckle image. We find that TOI-969 is a single star with no companion brighter than 5 to 8 magnitudes from the diffraction limit out to 1.2 arcsec. At the distance to TOI-969 (77 pc, see Sect. 3), these angular limits correspond to spatial limits of 1.5 to 92 AU.

We also observed TOI-969 on 15 October 2019 using the ShARCS camera on the Shane 3m telescope at Lick Observatory (Kupke et al. 2012; Gavel et al. 2014; McGurk et al. 2014). Observations were taken with the Shane adaptive optics system in natural guide star mode in order to search for nearby, unresolved stellar companions. We collected two sequences of observations, one with a KS filter (λ0 = 2.150 μm, Δλ = 0.320 μm) and one with a J filter (λ0 = 1.238 μm, Δλ = 0.271 μm), and reduced the data using the publicly available SImMER pipeline (Savel et al. 2020)3. Our reduced images and corresponding contrast curves are shown in Fig. 3 (central and right panels). We find no nearby stellar companions within our detection limits.

thumbnail Fig. 3

High-spatial-resolution contrast curves for TOI-969 from the Zorro speckle instrument on Gemini South (blue curve for the 562 nm band and green curve for the 832 nm band), and the 3 m Shane telescope in the J-band (orange) and Ks-band (red). The two dotted lines correspond to 10% (upper) and 5% (lower) contamination levels in the light curve. The dashed line corresponds to the maximum contrast that a blended binary could have to mimic the transit depth of TOI-969 b.

2.3 High-resolution spectroscopy

2.3.1 HARPS

TOI-969 was observed intensively with the HARPS spectrograph (Mayor et al. 2003) on the ESO 3.6 m telescope at La Silla Observatory, Chile, from 2020 November 10 to 2021 March 24. In total, 66 spectra were obtained under the programmes 1102.C-0249 (PI: Armstrong), and 106.21TJ.001 (PI: Gandolfi). HARPS is a stabilized high-resolution spectrograph with a resolving power of R ~ 115 000, capable of sub-meter-per-second RV precision. We used the instrument in high-accuracy mode with a 1 arcsec science fibre on the star and a second fibre on sky to monitor the sky-background during exposure.

Radial velocities were determined with the standard (online) HARPS data reduction pipeline using a K5 binary mask for the cross-correlation (Pepe et al. 2002), and a K5 template for flux correction to match the slope of the spectra across Echelle orders. With a typical signal-to-noise ratio (S/N) per pixel of 30, we achieved an RV precision of 2.5 m s−1. For each epoch, the bisector span (BIS), contrast, and FWHM of the CCF were calculated, as well as the chromospheric activity indicators CaII H&K, Hα, and Na. These RVs and the activity indicators are presented in Table A.4. We clearly detect an RV signal for TOI-969 b, as well as an additional long-term trend (see Fig. 4). When removing the latter signal from a simple Keplerian fit (a linear or quadratic fit is not enough to remove the long-term variability), the Lomb-Scargle periodogram of the HARPS measurements (see Fig. 5) shows three peaks at around 8, 12, and 24 days, which are compatible with those seen in the WASP photometric data (see Sect. 2.1.4 and the periodogram in the WASP panel of Fig. 5). Indeed, these periodicities also appear in the activity indicators corresponding to the FWHM of the CCF and the bisector span (third and fourth panels of Fig. 5, respectively). The signal of TOI-969 b is not statistically significant in the periodogram of this dataset (FAP > 10%) alone but a clear peak at the 1.8-day periodicity clearly stands above other signals.

2.3.2 Planet Finder Spectrograph

We observed TOI-969 with the Planet Finder Spectrograph (PFS; Crane et al. 2006, 2008, 2010), which is mounted on the 6.5 m Magellan II (Clay) Telescope at Las Campanas Observatory in Chile. PFS is a slit-fed echelle spectrograph with a wavelength coverage of 3910-7340 Å. We used a 0.3′′ slit and 3 × 3 binning, which yields a resolving power of R ≈ 110 000. Wavelength calibration is achieved via an iodine gas cell, which also allows characterization of the instrumental profile. We obtained ten spectra, observed through iodine, between 27 October 2020 and 2 January 2021. Exposure times ranged from 20 to 30 min. We also obtained an iodine-free template observation with a 90-min exposure time. The RVs were extracted using a custom IDL pipeline following the prescriptions of Marcy & Butler (1992) and Butler & Marcy (1996), and achieved a mean internal precision of 1.2 ms−1. The velocities are presented in Table A.4 and Fig. 4.

2.3.3 CORALIE

TOI-969 was monitored by the CORALIE high-resolution echelle spectrograph mounted on the 1.2 m Euler telescope at La Silla Observatory (Queloz et al. 2001). We obtained 22 spectra between 13 February 2021 and 15 April 2022 each with S/N of 10–20 depending on sky-conditions and exposure time, which was set between 1800 and 2700 seconds. The spectrum corresponding to the night of 16 January 2022 was discarded due to high contamination from the Moon. Radial velocity measurements were extracted by cross-correlating each spectra with a binary G2 mask (Baranne et al. 1996) using the standard CORALIE data-reduction pipeline. Given the relatively low S/N, a modest RV precision of 30–50 m s−1 was achieved. Line-shape diagnostics such as bisector-span and FWHM were derived for the cross-correlation function (CCF). The extracted RVs and associated uncertainties as well as the corresponding activity indicators obtained from the CCF analysis are presented in Table A.4.

thumbnail Fig. 4

Radial velocity time series and model inference from the joint radial velocity and light curve analysis in Sect. 4.3, including the GP noise (blue dashed line) and full RV model (black solid line for the median posterior model and gray-shaded regions for 68.7 (dark gray) and 99.7% (light gray) confidence intervals).

3 Stellar parameters

TOI-969 is a moderately bright (V=11.6, Høg et al. 2000) star in the late-K dwarf regime. Table 1 summarises the main properties of this star. According to the Gaia (Gaia Collaboration 2016) EDR3 data release (Gaia Collaboration 2022), this star is located at a distance of 77.82 ± 0.14 pc (corresponding to a parallax of 12.849 ± 0.023 mas, Lindegren et al. 2021). Based on the Gaia proper motions, we determine the projected galactic velocities of TOI-969 (see Table 1) and, using the relations from Bensby et al. (2003), we can conclude that this star likely belongs to the Galactic thin disk (with a probability 103 times higher than belonging to the thick disk and more than 104 times higher than belonging to the halo).

We derived the spectroscopic stellar parameters (Teff, log g, microturbulence, [Fe/H]) using the HARPS spectra and ARES+MOOG code, following the same methodology described in Sousa (2014) and Santos et al. (2013). The equivalent widths (EWs) of iron lines were measured on the combined HARPS spectrum of TOI-969 using the ARES code4 (Sousa et al. 2007, 2015). The best set of spectroscopic parameters was found when we reach the ionization and excitation equilibrium. This process makes use of a grid of Kurucz model atmospheres (Kurucz 1993) and the radiative transfer code MOOG (Sneden 1973). As the star is cooler than 5200 K, we used the appropriate iron linelist for our method presented in Tsantaki et al. (2013). Following the same methodology as described in Sousa et al. (2021), we used the Gaia EDR3 parallax and estimated a trigonometric surface gravity of log ɡ = 4.55 ± 0.06 dex. Stellar abundances of the elements were derived using the classical curve-of-growth analysis method assuming local ther-modynamic equilibrium (e.g., Adibekyan et al. 2012, 2015). The same codes and models were used for the abundance determinations. Finally, the mass, radius, and age of the star were derived using the PARAM 1.3 web interface5 (da Silva et al. 2006) using our spectroscopic parameters as input values. The results from this analysis are shown in Table 2 (second column).

We use the S-index provided by the HARPS DRS and calibrated to the Mount Wilson scale (Vaughan et al. 1978) to determine the log R′HK${R'_{{\rm{HK}}}}$ for each of the individual spectra. We obtain an average value of –4.381 ± 0.002 dex with a dispersion of 0.026 dex. This dispersion is significantly larger than the uncertainty of the individual measurements which points to a certain level of activity. Indeed, in the periodogram of the log R′HK${R'_{{\rm{HK}}}}$ time series (see corresponding panel of Fig. 5), a clear peak at 24.6 days stands out, again pointing to this periodicity as the rotaton period of the star. We can use the average value of the log R′HK${R'_{{\rm{HK}}}}$ to estimate the rotation period of the star. Given the B – V color of this star being 1.09 ± 0.33 mag, we can apply the Noyes et al. (1984) activity-rotation relationships through the pyrhk code (da Silva et al. 2006) and using the Middelkoop (1982) bolometric corrections. Following this, we obtain a rotation period of 9.9 ± 3.2 days. If we instead use the empirical relations from Suárez Mascareño et al. (2016) for K-dwarfs, we obtain a rotation period of 9.5 ± 1.0 days, which is compatible with the previous value. Hence, the activity-rotation empirical relations point to a period shorter than the one measured from the periodogram of the activity indicators (see Sect. 2.3 and panels 3-5 in Fig. 5) and the photometric time series from WASP (see Sect. 2.1.4 and panel 6 in Fig. 5). As an additional test, we computed the υ sin i⋆ of TOI-969 from the HARPS spectrum and obtained a value of 2.87 +/–0.37 km s−1. Assuming an inclination of 90°, this corresponds to a rotation period for this star of 12.1 ± 1.7 days, again compatible within the uncertainties with the previous results from empirical relations. For a rotation period of around 24 days, and assuming the same spin axis inclination of 90°, we obtain a projected rotation velocity of around 1.3 km s−1. We then warn that for very slow-rotating stars (υ sin i⋆< 3 km s−1), measuring accurate υ sin i⋆ values is a difficult task and our determination should be taken with care. Consequently, the rotation period of the star is unclear from the present data. In our analysis, we assume the rotation period is Prot = 24.6 ± 0.6 days but we account for the power at 12.3 days in the different indicators and empirical relations by using a rotation kernel that also accounts for Prot/2.

We can now use gyrochronology to estimate the age of the system. By applying the Angus et al. (2019) relations through the stardate6 code, we obtain an age of 2.03 ± 0.17 Gyr when assuming Prot = 24.6 ± 0.6 and the Gaia Bp – Rp color of 1.42 mag Gaia Collaboration (2021). If we use the 12.3-day rotation period obtained from the u sin i⋆, we obtain an age of 670 ± 150 Myr. Given this young age estimation in the second case, we checked for the presence of lithium in the HARPS spectra and found none within the sensitivity limits. However, for this effective temperature and age it is expected that the Li has already been depleted, and so we cannot discard this young estimation based on the absence of lithium.

As an independent determination of the basic stellar parameters, we performed an analysis of the broadband spectral energy distribution (SED) of the star together with the Gaia EDR3 parallax (with no systematic offset applied; see, e.g., Stassun & Torres 2021) in order to determine an empirical measurement of the stellar radius, following the procedures described in Stassun & Torres (2016) and Stassun et al. (2017, 2018). We took the BT and VT magnitudes from Tycho-2, the JHKS magnitudes from 2MASS, the W1–W4 magnitudes from WISE, the GBP and GRP magnitudes from Gaia, and the NUV magnitude from GALEX. Together, the available photometry spans the full stellar SED over the wavelength range 0.2–22 μm (see Fig. 6).

We performed a fit using NExtGen stellar atmosphere models, with the free parameters being the effective temperature (Teff) and metallicity ([Fe/H]). The remaining free parameter is the extinction AV, which we fixed at zero due to the proximity of the star. The resulting fit (Fig. 6) has a reduced χ2${\chi ^2}$ of 1.1, excluding the GALEX NUV flux which indicates a moderate level of activity (see below). Integrating the (unreddened) model SED gives the bolometric flux at Earth, Fbol = 9.45 ± 0.11 × 10−10 erg s−1 cm−2. We estimate the stellar radius by taking the Fbol and Teff together with the Gaia parallax, and the stellar mass using the empirical relations of Torres (2010) and Mann et al. (2020). All these values agree with the spectroscopically derived parameters within the uncertainties and are reported in Table 2 (third column).

thumbnail Fig. 5

Periodogram of the radial velocity and activity indicator time series of the HARPS dataset. Panel a corresponds to the original RV time series, while panel b corresponds to the time series after a simple linear detrending. In panel c we include an additional linear detrend-ing corresponding to the activity indicator most correlated with the RVs (in this case the BIS-span). Panels d–e include the periodogram of the FWHM of the cross-correlation function and the bisector span, respectively. Finally, the bottom panel shows the window function corresponding to the HARPS sampling. The periodicity corresponding to TOI-969 b and TOI-969 c are shown as vertical dotted lines in all panels. The stellar rotation period and its first two harmonics are shown as dashed vertical lines.

4 Analysis and results

Due to the complexity of the system and the large amount of data in hand, we first perform an analysis of the individual datasets, namely the TESS light curve (Sect. 4.1) and the RV (Sect. 4.2). From those independent analyses, we obtain relevant information about the models to be tested in the final joint (light curves and RV) analysis in Sect. 4.3, and parameter ranges to appropriately set up the priors in this more complex and computationally expensive joint model.

Table 1

General and dynamical properties of TOI-969 obtained in this work.

Table 2

Stellar physical properties of TOI-969.

4.1 Light-curve analysis

We first inspect the photometric data through an independent analysis of the TESS and LCOGT light curves. Our model is composed of a transit model assuming a single planet in the system (the one with the detected transits) defined by the period (Pb), a time of inferior conjunction (T0,b), the orbital inclination the planet radius (Rb), and the stellar mass (M⋆) and radius (R⋆). The latter are used to check for possible phase-curve-variation effects, such as the ellipsoidal, Doppler beaming, and reflection effects (which in turn are all negligible in this case). We use the limb-darkening parametrization described in Kipping (2013) and apply the quadratic law with a pair of parameters per instrument and bandpass(u1.j,u2.j). In addition, we assume a photometric jitter and offset for each instrument (σj and F0j, respectively). Finally, in the case of the ground-based observations, we use the pre-detrended time series and use a linear detrending with the airmass within our model by including a slope per bandpass for the LCOGT observations (mX,j).

In order to account for the stellar variability seen in the TESS sectors, we also include a Gaussian process (GP) with a kernel composed of a mixture of harmonic oscillators designed to model stellar rotation and implemented in the celerite2 (Foreman-Mackey et al. 2017, 2018) package as the term named RotationTerm. One of the harmonic oscillators is devoted to the rotation period and the other one to half of the rotation period. This kernel includes five hyperparameters, namely an amplitude (ησ,LC), a periodicity (ηP), a quality factor for the second oscillator (ηQ0), a difference between the two quality factors of both oscillators (ηδQ), and a factor representing the relevance of one oscillator as compared to that of the other (ηf).

Based on this model, which includes 35 parameters, we explore the parameter space informed by the observed time series using the Monte-Carlo Markov chain (MCMC) affine-invariant ensemble sampler emcee (Foreman-Mackey et al. 2013). We launch our MCMC with 140 walkers (four times the number of parameters) and a total of 100 000 steps per walker. The priors for all parameters are detailed in Table A.5. We use Gaussian priors for the orbital period and time of mid-transit as published by the TESS alert and using a broad standard deviation corresponding to five times the published uncertainty on these values. Gaussian priors are also used for the limb-darkening coefficients, with the mean value corresponding to the derived value using the limb-darkening code by Espinoza & Jordan (2015) using the ATLAS models and the stellar properties derived in Sect. 3; and a standard deviation of 0.1. The stellar mass and radius priors are also included as Gaussian distributions with the parameters described in Sect. 3. The priors for the remaining parameters are set to uniform distributions within the ranges stated in Table A.5.

We perform a first burn-in phase with the previously mentioned setup. Subsequently, we use the maximum a posteriori set of parameters from this first phase to start a second chain initialized at those values and running for half the number of steps as in the first phase (i.e., 50 000 steps per walker). The posterior distributions are computed using the full chains from the second phase. The median and 68.7% credible interval for each parameter are shown in Table A.5. The phase-folded light curves showing the median transit models from this analysis are shown in Fig. 7.

thumbnail Fig. 6

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

4.2 Radial velocity

We now use the radial velocity dataset described in Sect. 2.3 as well as the input information from the previous TESS analysis to scrutinize the Doppler information. To this end, we first build our model, which is composed of N Keplerians, each of them potentially including an orbital period (Pi), a time of inferior conjunction (T0,i), a radial velocity semi-amplitude (Ki), and, in the case of noncircular orbit models, we also add eccentricity(ei) and the corresponding argument of the periastron (ωi) Additionally, our models can include a slope (m) to account for possible linear trends (when this slope is assumed zero, we include an “F” (for “flat") in the model labels) and a quadratic trend (when this is added, a “Q” is added to the model label). We also include a radial velocity offset per instrument (γj) and a jitter per instrument to account for all systematic errors not taken into account in the model (σj). In our analysis, we test models that aim to account for the radial-velocity variations caused by the stellar activity (potentially relevant as seen in the previous section in the TESS light curve) through GP informed by an activity indicator. Based on the radial velocity periodogram presented in Sect. 2.3, the bisector span seems to match the radial velocity variations in a similar manner, showing strong periodicities at similar frequencies. This does not happen with the FWHM, which does not show significant variations other than a long-term periodicity. Hence, we use the BIS-span indicator to inform the GP about the activity-induced component of the radial velocity. Correspondingly, for each instrument, we add an offset (γBis,j) and a jitter (σBIS,j).

We include a GP for the RVs and a GP for the BIS-span, both sharing all parameters except the amplitude of the kernel. We use the celerite2 (Foreman-Mackey et al. 2017, 2018) implementation to build the GPs and we choose the RotationTerm kernel (see Sect. 4.1) to account for the stellar activity. All hyperparam-eters but the amplitude one (ησ) are shared between the RV and BIS-span GPs.

The prior distributions used in this analysis are summarised in Table A.6. Uninformative priors were generally used in most of the parameters involved in our model, except in some cases. We used Gaussian priors for the orbital period and mid-transit time for TOI-969 b centered at the values coming from the transit analysis from the previous section. This is because the precision in those parameters achievable with the transit method is higher than with RVs. We used log-uniform priors for the orbital period of the long-period candidate and the RV and BIS jitter, as well as for the amplitude hyperparameters of the GP (ησ,FWHM and ησ). Finally, given the discussion presented in Sect. 3, we opted for a Gaussian prior on the hyperparameter related to the rotation period of the star.

The parameter space is explored by using emcee. We use a number of walker equal to four times the amount of free parameters and 100 000 steps per walkers in a first burnout phase. The maximum a posteriori model from this first phase is selected as the initial guess for the final production, which now contains 100 000 steps and the same number of walkers as the initial phase. The convergence of the chains is checked by estimating the autocorrelation time and the corresponding chain length, with the latest being at least 20 times longer than the autocorrelation time to consider convergence. We then use 15% of the final flattened chain (typically composed of 105–106 elements) to estimate the Bayesian evidence of each model (In 𝒵i) and its corresponding uncertainty through the perrakis implementation7 (Díaz et al. 2016).

The preliminary analysis of the data already revealed the presence of different signals (see Sect. 2.3). Consequently, we test models with different numbers of planets, each potentially having either circular (fixed values for ei = 0 and ωi = 90°) or eccentric orbits. We also tested different types of trends, from the simplest linear trend (including a slope parameter) to a quadratic one (including a quadratic term) for the cases of zero-planet (labelled “Op*”) and one-planet (“lp*”) models. In total, 21 models with different levels of complexity are tested (from 0 to 3 planets). The labelling for the models follows the nomenclature Xp[Pic] [Z], where X is the total number of planets assumed in the system, Pic are the identifiers of the planets with assumed circular orbits (e.g., “2plc2c” if both planets “1” and “2” are assumed in circular orbits; or “2plc” if only planet “1” is assumed in circular orbit), and the latter letter (“[Z]”) indicates whether we are using a flat (“F”) slope (i.e., assuming slope = 0) or a quadratic (“Q”) trend.

The Bayesian evidence is then used to select the simplest model that best represents our data. This is done by sorting the models by complexity and selecting the one with the highest evidence, assuming that a difference of 6 in logarithmic space (i.e., ℬ = Δ ln 𝒵 > 6, Trotta 2008) corresponds to strong evidence in favor of the more complex model. Figure 8 shows this log-evidence for each of the 21 models. From this, we see that the model with a flat slope and two components in eccentric orbits (“2pF”) has the largest evidence. However, the odds against the simpler model with the inner component in a circular orbit “2plcF” are still not significant and we choose the latter as the preferred model informed by the current dataset. The evidence of this model against the no-planet model and the one-planet model is significantly larger (𝒵2p1cF – 𝒵0p✱ > +46 and 𝒵2p1cF – 𝒵1p✱ > +33), thus supporting the existence of another Keplerian signal in the system other than the transiting planet. On the other side, the “2plcF” model also has significantly larger evidence than the more complex three-planet models (𝒵2p1cF – > +31). The phase-folded radial velocity curves of both components for the preferred model “2plcF” are shown in Fig. 9. The posterior distributions of all parameters involved in this analysis are presented in Table A.6.

thumbnail Fig. 7

Phase-folded light curves from TOI-969 after removing the GP for the TESS data and the linear dependency for the LCOGT dataset. Left panels show the individual phase-folded light curves per instrument and the corresponding joint model. Black symbols correspond to 9-min bins. Right panels show the residuals from the joint model.
thumbnail Fig. 8

Comparison of the logarithm of the Bayesian evidence (ℬ = In 𝒵) for the different tested models (labeled in the X axis) in the case of the RV-only analysis (Sect. 4.2). The labels for the models correspond to Xp[Pic][Z], where X is the number of planets assumed in the system, Pic are the identifiers of the planets with assumed circular orbits, and Z indicates whether we are using a flat (F) slope (i.e., assuming slope = 0) or a quadratic (Q) trend. The gray-shaded region displays the region of Δ𝒵 = −6 from the larger evidence model showing its strength against simpler models.

4.3 Joint analysis

We use the previous individual analysis of the light curve and radial velocity datasets as exploratory studies of the parameter space and for model comparison purposes. Based on the results obtained in these preliminary analyses, we now run a joint model including both the radial velocity and photometric light curves, from which we obtain the final parameters of the system. In this case, the planetary and orbital parameters are parameterized to include the shared information between both techniques. That is, the Keplerians are modeled by a planet period, time of mid-transit, planet mass (mi) and radius (Ri), orbital inclination (ii), eccentricity, and argument of periastron. In addition, as detailed in Sects. 4.1 and 4.2, we use the RotationTerm kernel with shared hyperparameters (except for the amplitude) for the TES S light curves, radial velocity, and BIS-span time series.

We proceed in a similar fashion as explained in Sect. 4.2 to explore the parameter space of this joint analysis. We reach convergence of the chains under this strategy and confirm the planetary nature of the inner component and place the second component in the planet-to-brown-dwarf transition. Unfortunately, the TESS observations do not cover the conjunction time of TOI-969 c and consequently we cannot come to a conclusion as to its transiting nature. The final set of parameters (including the absolute mass for TOI-969 b) are presented in Table A.7 and a summary of the main properties of the two components is presented in Table 3.

thumbnail Fig. 9

Phase-folded radial velocity data from TOI-969 b (left) and TOI-969 c (right) after removing the signal from the other component in the system and the GP noise for the selected model (labeled “2plcF") from the radial velocity analysis (see Sect. 4.2). For TOI-969 b, we only include HARPS and PFS data, while removing CORALIE for improved clarity given the larger uncertainties for that dataset. The gray-shaded regions correspond to the region covered by models within 68.7% (dark gray) and 95% (light gray) of the confidence interval.

4.4 Transit-timing variations for TOI-969 b

Given the relatively large number of transits of TOI-969 b in our photometric dataset, we can study the presence of additional planets in the system based on the gravitational pull on this planet. We estimate the central times of individual transit events of TOI-969 b based on deviations from a single-planet model. To this end, we use the results of the joint modeling of Sect. 4.3 (Table A.7) as priors in this analysis. As in the global analysis, we also use GPs to detrend the TESS light curves. On the contrary, we use the already detrended LCOGT/MusCat light curves to estimate a single transit mid-time out of the simultaneous modeling of the four-filter light curves. We use the Juliet package8 (Espinoza et al. 2019) to perform this analysis. By doing this, we estimate the central times of the 25 TESS transits and 1 of the LCOGT light curves (Table 4). We compare these observed times with the ephemeris reported in Table A.7. The resulting observed minus calculated (O–C) diagram of the transit times of TOI-969 b are shown in Fig. 10. The RMS of the O-C points is 7.1 min with a maximum deviation of 14.13 min. We note that the amplitude of the timing deviations are larger in the transits observed in the last TESS sector (bottom panel of Fig. 10). However, it is unlikely that these deviations are related to the stellar activity; as can be seen in Fig. 2, the stellar activity is lower in this TESS sector.

Based on the constraints of these transit-timing results, we explore the limits that can be placed on the mass of a potential unseen companion to TOI-969 b. For this purpose, we used the TTV2Fast2Furious code9 (Hadden et al. 2019). The resulting upper limits on the companion are not very constraining, ruling out planets with masses above >100 M⊕ for periods larger than twice the orbital period of TOI-969 b. Better constraints can be placed with an extended timing monitoring of additional TOI-969 b transits.

5 Discussion

5.1 TOI-969 b: a new hot-Neptune in the desert

The derived properties of planet TOI-969 b shown in the previous sections place this planet in the so-called mini-Neptune regime in terms of mass (mb=9.1−1.0+1.1M⊕)$\left( {{m_b} = 9.1_{ - 1.0}^{ + 1.1}{M_ \oplus }} \right)$ and radius (Rb2.765−0.097+0.088R⊕)$\left( {{R_{\rm{b}}}2.765_{ - 0.097}^{ + 0.088}{R_ \oplus }} \right)$ see Fig.11 Ks relatively short period also places it in the boundary of the hot-Neptune desert (Mazeh et al. 2016) and suggests that it possesses similar properties to the iconic GJ 1214 b (Charbonneau et al. 2009).

As for most Neptune-like planets, for a planet of the mass and radius of TOI-969 b and a H/He-dominated atmosphere, the atmospheric loss by blowoff is energy limited (Owen & Alvarez 2016). In this case, the atmospheric mass-loss rate is (Erkaev et al. 2007; Owen & Wu 2013): M˙=∈πFXUVRp3GMp,$\dot M = \in {{\pi {F_{{\rm{XUV}}}}R_{\rm{p}}^3} \over {G{M_{\rm{p}}}}},$(1)

where FXUV is the XUV flux received by the planet, G the gravitational constant, and ε is an efficiency parameter. We approximate the XUV luminosity by the analytical fit obtained by Sanz-Forcada et al. (2011), and we estimate ε ≃ 0.07 from Owen & Jackson (2012). This yields a present-day mass-loss rate of 0.03 M⊕ Gyr−1, assuming the age of the star is 2.03 Gyr. Following the approach of Aguichine et al. (2021), the total mass lost by the planet can be estimated by integrating m˙$\dot m$ over time assuming that only the XUV luminosity changes, and that the planet properties have remained roughly constant during its evolution. In this case, we find that TOI-969 b would have lost ~0.42 M⊕ of H/He during its evolution. The possibility that TOI-969 b formed with a H/He envelope of mass greater than ~0.42 M⊕ (~ 4% by mass) cannot be excluded. It is therefore possible that TOI-969 b is still subject to atmospheric escape, and could be entirely stripped of its H/He envelope. However, the presence of an atmosphere made of heavier volatiles, such as He, O2, H2O, and others (see Hu et al. 2015; Bolmont et al. 2017; Aguichine et al. 2021; Ito & Ikoma 2021, respectively), provides a more natural explanation as to how TOI-969 b retained its volatile envelope, as these molecules have smaller escape rates (Owen & Jackson 2012; Ito & Ikoma 2021).

TOI-969 b is therefore another excellent target for testing theories of atmospheric evaporation (Lecavelier Des Etangs 2007; Owen & Adams 2019). Interestingly, its transmission spectroscopy metric (TSM, Kempton et al. 2018) corresponds to TSM = 93 ± 18, making it one of the best targets in this regime for atmospheric purposes, especially with the JWST.

Table 3

Summary of the posterior distributions of the main physical and orbital properties of TOI-969 b and TOI-969 c from the joint analysis (Sect. 4.3).

5.2 Internal structure of TOI-969 b

With the precise mass and radius determinations for TOI-969 b, we can study its internal structure. We performed an MCMC Bayesian analysis (Dorn et al. 2015) of its internal composition using the interior structure model introduced in Mousis et al. (2020) and Brugger et al. (2017), which comprises three layers: a Fe-rich core, a silicate-rich mantle, and a water layer. With an irradiance temperature10 of the order of 103 K, TOI-969 b can present vapor and supercritical phases if water is found at its surface. Water is the second-most abundant volatile after H/He. The supercritical phase can be reached for planets that are highly irradiated and present a volatile envelope that reaches high temperatures (around 103 K) at the bottom of this layer, which corresponds to high pressures of approximately 300 bar (see its phase diagram in, e.g., Fig. 1 in Mousis et al. 2020).

Therefore, we couple an atmosphere-interior model that calculates the surface conditions and the contribution of the atmosphere to the total radius, which in the case of low-density, warm planets is significant, as shown by Acuña et al. (2021). We consider two scenarios to obtain the interior structure of TOI-969 b: scenario 1, where only the mass and radius of the planet (from Table A.7) are considered as inputs to the MCMC analysis; and scenario 2, where the planetary mass and radius, and the stellar Fe/Si mole ratios (see Table 2) are the input data11. To compute the Fe/Si and Mg/Si mole ratios with the stellar abundances, we follow the approach described in Brugger et al. (2017) and Sotin et al. (2007), and obtain Fe/Si = 0.77±0.23 and Mg/Si = 1.01 ±0.44. The outputs of the MCMC analysis are the posterior distribution functions (PDF) of the core mass fraction (CMF), the water mass fraction (WMF), and the atmospheric parameters, which are the temperature at 300 bar, the planetary albedo, and the atmospheric thickness from transit pressure (20 mbar, see Mousis et al. 2020; Grimm et al. 2018) to 300 bar12. We assume a water-rich atmosphere. Table 5 (upper part) shows the 1D, 1σ confidence intervals of the MCMC output parameters and Fig. 12 presents the ternary diagram and the 1σ confidence region for the internal structure of TOI-969 b. Under the assumption of the absence of an H/He layer, TOI-969 b is a volatile-rich planet that could have up to 60% of its mass in a water-rich volatile layer in both scenarios. The CMF is compatible with an Earth-like value in scenario 1, whereas in scenario 2 the CMF is lower due to a lower Fe/Si mole ratio of the host star compared to the solar value (Fe/Si⊙ = 0.96).

H/He and water are the most abundant volatiles. In interior modeling, when including a volatile layer, it is widely acknowledged to assume that the main component of the envelope is either hydrogen, which is representative of a primordial atmosphere, or water, whose density better represents that of a secondary atmosphere. The interior structure model described above presents the implementation of a water layer. However, as we describe in Sect. 5.1 the presence of H/He in the current atmosphere of TOI-969 b cannot be discarded based on atmospheric loss. Moreover, the water layer model described above yields a maximum WMF of almost 60% within 1σ (WMF = 50 ± 9%). This is less than the maximum WMF found in Solar System bodies (i.e., comets), which is 70–80% (McKay et al. 2019). However, it is still a high WMF compared to other low-mass exoplanets. Hence, we also explore the possibility of a H/He atmosphere instead of a water-dominated volatile layer.

To do this, we combine our interior structure model (only the Fe-rich core and the mantle) with the mass-radius relations of Zeng et al. (2019). We obtain the atmospheric thickness by subtracting the mass-radius relationship of an Earth-like bulk with a H/He envelope minus that of a bare Earth-like core, both as provided by Zeng et al. (2019)13. This allows us to express the atmospheric thickness as a function of the surface gravity, g0 = M/R2, and the H/He mass fraction, zatm = zatm(g0,xH/He). The boundary surface conditions for our interior model are Tsurf = 2000 K and Psurf = 1 bar, because Zeng et al. (2019) consider a maximum irradiation temperature of 2000 K and an isothermal profile in the atmosphere. Table 5 (bottom part) shows the observables and the compositional parameters retrieved from this second analysis assuming a H/He atmosphere. We observe that the H/He volatile mass fraction is four orders of magnitude lower that the water mass fraction obtained in the previous analysis, which is expected because water at high pressures is significantly more dense than H and He. The 1σ confidence intervals of the H/He volatile mass fractions are situated between 0.1 and 0.3%, as expected from the position of TOI-969b in the mass–radius diagram in Fig. 13. To reproduce the density of TOI-969 b, the H/He case does not need a volatile layer as massive as that of the water case because a H/He atmosphere is more expanded than a secondary atmosphere for a similar atmospheric mass. This results in a larger portion of the total mass being constituted by the Fe-rich core when we consider a H/He atmosphere, yielding higher CMFs in both scenarios. In scenario 2, we can observe that the CMF is compatible within uncertainties with the Earth value (CMF⊕ = 0.32).

Finally, we applied the stoichiometric model of Santos et al. (2017) to determine the iron-mass fraction (which can be translated into core-mass fraction) of the planet building blocks in the protoplanetary disks of TOI-969. This model is based on the chemical abundances of this star listed in Table 2, and assuming that C and O abundances for the two stars scale with metallicity. We find an iron mass fraction of the planet building blocks of 33.1 ± 5.1%. This value is compatible with that predicted by this model for the Solar System planet-building blocks (33.2%, Santos et al. 2015).

Table 4

Transit times of TOI-969 b estimated from the TESS and LCOGT light curves.

thumbnail Fig. 10

Observed–calculated diagram of the transit times of TOI-969 b. LCOGT transit corresponds to transit number #400. The gray regions represent the 1σ uncertainties of the ephemeris reported in Table A.7.

5.3 Prospects for atmospheric studies with JWST

We simulated synthetic spectra for a range of atmospheric scenarios and instrumental setups. We adopted Tau-REx III (Al-Refaie et al. 2021) to compute the model-atmosphere spectra, using the stellar parameters from Table 2 (arithmetic average) and the planetary parameters from Table A.7. Our test-bench models assume atmospheric chemical equilibrium (ACE, Agúndez et al. 2012) with an isothermal profile at the equilibrium temperature, scaled (1× and 100×) solar abundances, collisionally induced absorption (CIA) by H2-H1 and H2-He (Abel et al. 2011, 2012; Fletcher et al. 2018), and Rayleigh scattering. We show two cloud-free models with different metallicities, and a third model with solar metallicity and an optically thick cloud deck with top pressure of 100 Pa. The cloud-free spectra present molecular absorption features of 100–300 ppm, which are slightly smaller by a factor of about 2 in the case with enhanced metallicity. The cloudy model also presents smaller absorption features due to the suppression of contributions from deeper atmospheric layers. A flat spectrum due to higher altitude (lower top pressure) clouds cannot be excluded for TOI-969 b.

We used ExoTETHyS14 (Morello et al. 2021) to compute bin-averaged spectra, taking into account the spectral response of the JWST instruments, noise scatter, and error bars. We simulated JWST spectra for the NIRISS-SOSS (0.6–2.8 µm), NIRSpec-G395M (2.88–5.20µm), and MIRI-LRS (5–12µm) instrumental modes. The wavelength bins were specifically determined so as to have similar counts, leading to nearly uniform error bars per spectral point. In particular, we set a median resolving power of R ~ 50 for the NIRISS-SOSS and NIRSpec-G395M modes, and bin sizes of ~0.l–02µm for the MIRI-LRS. The error bars have been calculated for a single visit of twice the transit duration in each instrumental mode, scaled by the inverse of the observing efficiency estimated with the Exoplanet Characterization Toolkit (ExoCTK15, Bourque et al. 2021), and a factor 1.2 to account for correlated noise. We obtained error bars of 35–0 ppm per spectral point for the NIRISS-SOSS and NIRSpec-G395M modes, and 110–14 for the MIRI-LRS bins. These numbers suggest that a single transit observation is sufficient to sample the molecular absorption features in case of a clear atmosphere, even with 100× solar metallicity. A single NIRSpec-G395M observation could also be sufficient to detect the absorption features in a cloudy scenario, if the top pressure is higher than 100 Pa (see Fig. 14).

thumbnail Fig. 11

Comparative properties of TOI-969b. Left: mass–radius diagram of known exoplanets with measured mass precisions better than 30% (gray symbols) and the location of TOI-969b. The bulk density lines corresponding to different compositions from Zeng et al. (2019) are shown as solid traces, and the dashed lines correspond to iso-densities of 0.3, 1, and 3 gm−3. Right: period-radius diagram of known exoplanets (gray crosses) highlighting those with measured masses (open black circles). The light-gray-shaded region represents the mean boundaries of the hot-Neptune desert derived by Mazeh et al. (2016) and the dark-grey-shaded region corresponds to the interior envelope of the \cr boundaries. The location of TOI-969 b is marked in red, as are the locations two other planets in the hot-Neptune desert, namely TOI-849 b (Armstrong et al. 2020) and LTT 9779 b (Jenkins et al. 2020).
Table 5

Median and 1σ confidence intervals of the interior and atmosphere modeling.

thumbnail Fig. 12

Ternary diagram with the 1σ confidence region of TOI-969b in scenario 1 (black), where only the mass and the radius are considered as input data, and in scenario 2 (red), where the stellar abundances are also included as input data in the MCMC interior structure analysis. The mantle mass fraction (MMF) is defined as MMF = 1 - CMF - WMF. The green dot and red square indicate the position of Earth and Mercury in the ternary diagram, respectively.
thumbnail Fig. 13

Mass-radius relationships for supercritical water (SW) planets (Acuña et al. 2021; Mousis et al. 2020), planets with Earth-like cores and H/He atmospheres (Zeng et al. 2019), and dry planets with different Fe contents (Brugger et al. 2017). We assume high chemical equilibrium temperatures for supercritical and H/He volatile layers of 1200 K and 2000 K, respectively. The position of TOI-969 b is indicated as a magenta star.

5.4 TOI-969 c: a massive eccentric cold component in the planetary-to-brown-dwarf transition

Our RV analysis strongly favors the presence of an outer companion at 2.4 AU in a high eccentric orbit that we call TOI-969 c. With a minimum mass of 11.3−0.9+1.1$11.3_{ - 0.9}^{ + 1.1}$ MJup and a derived eccentricity of ec=0.628−0.036+0.043$11.3_{ - 0.9}^{ + 1.1}$, our analysis presents this planetary system as an ideal case for dynamical studies of planet migration and formation. Figure 15 displays the orbital shape of the system. These high eccentricities are expected for massive planets as demonstrated in Adibekyan et al. (2013) and previously predicted by numerical simulations (e.g., Papaloizou et al. 2001; Bitsch et al. 2013). In particular, one of the main scenarios to explain the current population of eccentric cold giant planets is that of planet-disk interactions during the planet formation phase (Papaloizou et al. 2001; Kley & Dirksen 2006; Bitsch et al. 2013). As discussed in Bitsch et al. (2020), in the case of very massive gas giants with more than 5 MJup, the gaps opened by these massive bodies are so deep that they prevent the Lindblad resonances that allow the damping of high eccentricities. This is one possibility to explain the large eccentricity of TOI-969 c in the absence of additional massive bodies in the system. The alternative scenario consists of planet-planet scattering. This has been tested by different authors using different approaches and configurations of the forming planetary system (e.g., Juric & Tremaine 2008; Raymond et al. 2009; Bitsch et al. 2020), leading to the conclusion that planetary systems including a massive eccentric giant should not harbor any inner super-Earths. By contrast, Zhu & Wu (2018) and Bryan et al. (2019) found that there is a positive correlation between the presence of hot super-Earths and cold Jupiters in the same system. This is therefore still an open question.

The presence of the hot mini-Neptune TOI-969 b at such a close-in orbit is therefore intriguing in this regard, as it seems to be an exception to the planet-planet scattering scenario, which should have destabilized the inner planet and removed it from the system. The alternative would then be that TOI-969 b formed inner to the cold giant and subsequently migrated toward its current location through type-I migration. However, the expected location of the snow line for the protoplanetary disk of this late-K dwarf star is at around 1 AU, which is too close to the current location of TOI-969 c to have had enough dynamical room to form a Neptune-like planet such as TOI-969 b. Consequently, the formation of this system and its evolution until reaching its present configuration remain unclear, and a complete re-arrangement of the initial planetary configuration cannot be discarded.

6 Conclusions

In this paper, we confirm the planetary nature of TOI-969 b, a planet candidate detected by the TESS mission through the transit technique. Additional transits of TOI-969 b were also observed from the ground using LCOGT and MEarth telescopes. We used high-spatial resolution images to discard possible contaminant sources within the TESS aperture as well as an intense RV monitoring with different instruments (namely HARPS, PFS and CORALIE). These RV time series reveal the presence of a long-term signal compatible with a new planet in the system.

We analyzed this large dataset and conclude that TOI-969 b is a hot (Pb ~ 1 82 days) mini-Neptune with a radius of 2.765−0.097+0.088$2.765_{ - 0.097}^{ + 0.088}$ R⊕ and a mass of 9.1−1.0+1.1$9.1_{ - 1.0}^{ + 1.1}$ M⊕. These properties place this planet within the lower boundary of the hot-Neptune desert. The internal structure analysis also suggests a core mass fraction of 12 ± 4% and a large fraction of its mass contained in a potential atmospheric layer. The simulations we performed suggest this atmosphere could be characterized with JWST in a single visit, thus putting this target in the priority list for sub-Neptune planets. The transit times of this planet also display some deviations from a strictly periodic signal at the 10-min level. However, the precision from the TESS data is not sufficient to properly assess a TTV analysis. These data therefore cannot offer further insight into the possible presence of additional planets beyond that provided by the RV data.

The long time-span radial velocity time series indicate the presence of a second component in the system (TOI-969 c) in a long-period orbit with Pc=1700−280+290${P_c} = 1700_{ - 280}^{ + 290}$ days. Our two-planet model is here favored against a simpler one-planet model plus a quadratic trend. This second body has a minimum mass of 11.3−0.9+1.1$11.3_{ - 0.9}^{ + 1.1}$MJup (in the transition between the planetary and brown-dwarf domain) and a large eccentricity of ec=0.628−0.036+0.043${e_c} = 0.628_{ - 0.036}^{ + 0.043}$, constituting a relevant piece of the system in terms of dynamical stability, and placing it as a key actor of TOI-969 in molding the current architecture of the planetary system. TOI-969 thus becomes a test bench for planet migration theories and in particular for planet-planet scattering studies.

thumbnail Fig. 14

Synthetic transmission spectra of TOI-969 b. Fiducial models assuming a cloud-free atmosphere with solar abundances (bottom), cloud-free atmosphere with enhanced metallicity by a factor of 100 (middle), and atmosphere with solar abundances and an optically thick cloud deck with top pressure of 100 Pa (top). Estimated uncertainties are shown for the observation of one transit with JWST NIRISS-SOSS, NIRSpec-G395M, and MIRI-LRS configurations. Simulated data are color coded according to wavelength.
thumbnail Fig. 15

Schematic view of the planetary system TOI-969, with the Earth direction toward the bottom of the panel. A total of 500 randomly selected orbits from the posterior distribution of the joint model (see Sect. 4.3) are drawn for TOI-969 b (orange) and TOI-969 c (blue). The orbit corresponding to the median solution for each parameter is displayed as a thicker solid line of the corresponding color for each planet.

Acknowledgements

We thank the anonymous referee for their though revision of this manuscript that has improved its final quality. J.L-B. acknowledges financial support received from “la Caixa” Foundation (ID 100010434) and from the European Unions Horizon 2020 research and innovation programme under the Marie Slodowska-Curie grant agreement No 847648, with fellowship code LCF/BQ/PI20/11760023. This research has also been partly funded by the Spanish State Research Agency (AEI) Projects No.PID2019-107061GB-C6l and No. MDM-2017-0737 Unidad de Excelencia “Maria de Maeztu” – Centro de Astrobiología (INTA-CSIC). R.L. acknowledges financial support from the Spanish Ministerio de Ciencia e Innovación, through project PID2019-109522GB-C52, and the Centre of Excellence “Severo Ochoa” award to the Instituto de Astrofísica de Andalucía (SEV-2017-0709). DJ.A. acknowledges support from the STFC via an Ernest Rutherford Fellowship (ST/R00384X/1). S.G.S acknowledges the support from FCT through Estimulo FCT contract nr.CEECIND/00826/2018 and POPH/FSE (EC). G.M. has received funding from the European Union’s Horizon 2020 research and innovation programme under the Marie Skłodowska-Curie grant agreement No. 895525. S.H. acknowledges CNES funding through the grant 837319. The French group acknowledges financial support from the French Programme National de Planétologie (PNP, INSU). This work is partly financed by the Spanish Mnistry of Economics and Competitiveness through grants PGC2018-098153-B-C31. We acknowledge the support by FCT – Fundação para a Ciência e a Tecnologia through national funds and by FEDER through COMPETE2020 – Programa Operacional Competitividade e Internacionalização by these grants: UID/FIS/04434/2019; UIDB/04434/2020; UIDP/04434/2020; PTDC/FIS-AST/32113/2017 & POCI-01-0145-FEDER-032113; PTDC/FISAST/28953/2017 & POCI-01-0145-FEDER-028953. P.J.W is supported by an STFC consolidated grant (ST/T000406/1). F.H. is funded by an STFC studentship. T.A.S acknowledges support from the Fundação para a Ciência e a Tecnologia (FCT) through the Fellowship PD/BD/150416/2019 and POCH/FSE (EC). C.M.P. acknowledges support from the SNSA (dnr 65/19P). This work has been carried out within the framework of the National Centre of Competence in Research (NCCR) PlanetS supported by the Swiss National Science Foundation. This project has received funding from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation programme (grant agreement SCORE No 851555). O.D.S.D. is supported in the form of work contract (DL 57/2016/CP1364/CT0004) funded by national funds through Fundação para a Ciência e a Tecnologia (FCT). M.E. acknowledges the support of the DFG priority programSPP 1992 “Exploring the Diversity of Extrasolar Planets” (HA 3279/12-1). A.O. is funded by an STFC studentship. J.K. gratefully acknowledge the support of the Swedish National Space Agency (SNSA; DNR 2020-00104). This work makes use of observations from the LCOGT network. This paper is based on observations made with the MuSCAT3 instrument, developed by the Astrobiology Center and under financial supports by ISPS KAKENHI (IP18H05439) and 1ST PRESTO (IPMIPR1775), at Faulkes Telescope North on Maui, HI, operated by the Las Cumbres Observatory. Some of the observations in the paper made use of the High-Resolution Imaging instrument Zorro obtained under Gemini LLP Proposal Number: GN/S-2021A-LP-105. Zorro was funded by the NASA Exoplanet Exploration Program and built at the NASA Ames Research Center by Steve B. Howell, Nie Scott, Elliott P. Horch, and Emmett Quigley. Zorro was mounted on the Gemini North (and/or South) telescope of the international Gemini Observatory, a program of NSF’s OIR Lab, which is managed by the Association of Universities for Research in Astronomy (AURA) under a cooperative agreement with the National Science Foundation. on behalf of the Gemini partnership: the National Science Foundation (United States), National Research Council (Canada), Agencia Nacional de Investigación y Desarrollo (Chile), Ministerio de Ciencia, Tecnología e Innovación (Argentina), Ministério da Ciência, Tecnologia, Inovações e Comunicações (Brazil), and Korea Astronomy and Space Science Institute (Republic of Korea). We acknowledge the use of public TESS data from pipelines at the TESS Science Office and at the TESS Science Processing Operations Center. Resources supporting this work were provided by the NASA High-End Computing (HEC) Program through the NASA Advanced Supercomputing (NAS) Division at Ames Research Center for the production of the SPOC data products. The MEarth Team gratefully acknowledges funding from the David and Lucile Packard Fellowship for Science and Engineering (awarded to D.C.). This material is based upon work supported by the National Science Foundation under grants AST-0807690, AST-1109468, AST-1004488 (Alan T. Waterman Award), and AST-1616624, and upon work supported by the National Aeronautics and Space Administration under Grant No. 80NSSC18K0476 issued through the XRP Program. This work is made possible by a grant from the John Templeton Foundation. The opinions expressed in this publication are those of the authors and do not necessarily reflect the views of the John Templeton Foundation. This research made use of Astropy, (a community-developed core Python package for Astronomy, Astropy Collaboration 2013, 2018), SciPy (Virtanen et al. 2020), matplotlib (a Python library for publication quality graphics Hunter 2007), and numpy (Harris et al. 2020). This research has made use of NASA’s Astrophysics Data System Bibliographic Services. This research has made use of the SIMBAD database, operated at CDS, Strasbourg, France. The authors wish to recognize and acknowledge the very significant cultural role and reverence that the summit of Mauna Kea has always had within the indigenous Hawaiian community. We are most fortunate to have the opportunity to conduct observations from this mountain.

Appendix A Additional tables

Table A.1

TESS light curves of TOI-969 (see Sect. 2.1.1).

Table A.2

Ground-based photometric time series of TOI-969 from the LCOGT and MEarth facilities (see Sects. 2.1.2 and 2.1.3).

Table A.3

Photometric time series from WASP (see Sect. 2.1.4).

Table A.4

Radial velocity and activity indicators for all instruments.

Table A.5

Inferred and derived parameters from the light-curve analysis (Sect. 4.1) for the one-planet model.

Table A.6

Inferred and derived parameters from the RV-only analysis (Sect. 4.2) for the preferred model (labeled “2p1cF”) with the inner planet in a circular orbit and the external component in an eccentric orbit.

Table A.7

Inferred and derived parameters from the joint radial velocity and light-curves analysis (Sect. 4.3) for the preferred model.

References

  1. Abel, M., Frommhold, L., Li, X., & Hunt, K. L. C. 2011, J. Phys. Chem. A, 115, 6805 [NASA ADS] [CrossRef] [Google Scholar]
  2. Abel, M., Frommhold, L., Li, X., & Hunt, K. L. C. 2012, J. Chem. Phys., 136, 044319 [NASA ADS] [CrossRef] [Google Scholar]
  3. Acuña, L., Deleuil, M., Mousis, O., et al. 2021, A&A, 647, A53 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  4. Adibekyan, V. Z., Sousa, S. G., Santos, N. C., et al. 2012, A&A, 545, A32 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  5. Adibekan, V. Z., Figueira, P., Santos, N. C., et al. 2013, A&A, 560, A51 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  6. Adibekyan, V., Figueira, P., Santos, N. C., et al. 2015, A&A, 583, A94 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  7. Adibekyan, V., Dorn, C., Sousa, S. G., et al. 2021, Science, 374, 330 [NASA ADS] [CrossRef] [Google Scholar]
  8. Aguichine, A., Mousis, O., Deleuil, M., & Marcq, E. 2021, ApJ, 914, 84 [NASA ADS] [CrossRef] [Google Scholar]
  9. Agúndez, M., Venot, O., Iro, N., et al. 2012, A&A, 548, A73 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  10. Akeson, R. L., Chen, X., Ciardi, D., et al. 2013, PASP, 125, 989 [Google Scholar]
  11. Al-Refaie, A. F., Changeat, Q., Waldmann, I. P., & Tinetti, G. 2021, ApJ, 917, 37 [NASA ADS] [CrossRef] [Google Scholar]
  12. Aller, A., Lillo-Box, J., Jones, D., Miranda, L. F., & Barceló Forteza, S. 2020, A&A, 635, A128 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  13. Ambikasaran, S., Foreman-Mackey, D., Greengard, L., Hogg, D. W., & O’Neil, M. 2015, IEEE Transactions on Pattern Analysis and Machine Intelligence, 38, 2 [Google Scholar]
  14. Anderson, D. R., Collier Cameron, A., Delrez, L., et al. 2017, A&A, 604, A110 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  15. Angus, R., Morton, T. D., Foreman-Mackey, D., et al. 2019, AJ, 158, 173 [Google Scholar]
  16. Armstrong, D. J., Lopez, T. A., Adibekyan, V., et al. 2020, Nature, 583, 39 [Google Scholar]
  17. Astropy Collaboration, (Robitaille, T. P., et al.) 2013, A&A, 558, A33 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  18. Astropy Collaboration, (Price-Whelan, A. M., et al.) 2018, AJ, 156, 123 [NASA ADS] [CrossRef] [Google Scholar]
  19. Baranne, A., Queloz, D., Mayor, M., et al. 1996, A&AS, 119, 373 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  20. Bensby, T., Feltzing, S., & Lundström, I. 2003, A&A, 410, 527 [CrossRef] [EDP Sciences] [Google Scholar]
  21. Berta, Z. K., Irwin, J., Charbonneau, D., Burke, C. J., & Falco, E. E. 2012, AJ, 144, 145 [Google Scholar]
  22. Bitsch, B., Crida, A., Libert, A. S., & Lega, E. 2013, A&A, 555, A124 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  23. Bitsch, B., Trifonov, T., & Izidoro, A. 2020, A&A, 643, A66 [EDP Sciences] [Google Scholar]
  24. Bolmont, E., Selsis, F., Owen, J. E., et al. 2017, MNRAS, 464, 3728 [NASA ADS] [CrossRef] [Google Scholar]
  25. Borucki, W. J., Koch, D., Basri, G., et al. 2010, Science, 327, 977 [Google Scholar]
  26. Bourque, M., Espinoza, N., Filippazzo, J., et al. 2021, The Exoplanet Characterization Toolkit (USA: ExoCTK) [Google Scholar]
  27. Brown, T. M., Baliber, N., Bianco, F. B., et al. 2013, PASP, 125, 1031 [Google Scholar]
  28. Brugger, B., Mousis, O., Deleuil, M., & Deschamps, F. 2017, ApJ, 850, 93 [NASA ADS] [CrossRef] [Google Scholar]
  29. Bryan, M. L., Knutson, H. A., Lee, E. J., et al. 2019, AJ, 157, 52 [Google Scholar]
  30. Butler, R. P., & Marcy, G. W. 1996, ApJ, 464, L153 [NASA ADS] [CrossRef] [Google Scholar]
  31. Caldwell, D. A. et al. 2020, RNAAS, 4, 201 [NASA ADS] [Google Scholar]
  32. Charbonneau, D., Berta, Z. K., Irwin, J., et al. 2009, Nature, 462, 891 [NASA ADS] [CrossRef] [Google Scholar]
  33. Ciardi, D. R., Beichman, C. A., Horch, E. P., & Howell, S. B. 2015, ApJ, 805, 16 [NASA ADS] [CrossRef] [Google Scholar]
  34. Cohen, M., Wheaton, W. A., & Megeath, S. T. 2003, AJ, 126, 1090 [Google Scholar]
  35. Crane, J. D., Shectman, S. A., & Butler, R. P. 2006, SPIE Conf. Ser., 6269, 626931 [NASA ADS] [Google Scholar]
  36. Crane, J. D., Shectman, S. A., Butler, R. P., Thompson, I. B., & Burley, G. S. 2008, SPIE Conf. Ser., 7014, 701479 [NASA ADS] [Google Scholar]
  37. Crane, J. D., Shectman, S. A., Butler, R. P., et al. 2010, SPIE Conf. Ser., 7735, 773553 [NASA ADS] [Google Scholar]
  38. da Silva, L., Girardi, L., Pasquini, L., et al. 2006, A&A, 458, 609 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  39. Delrez, L., Ehrenreich, D., Alibert, Y., et al. 2021, Nat. Astron., 5, 775 [NASA ADS] [CrossRef] [Google Scholar]
  40. Díaz, R. F., Ségransan, D., Udry, S., et al. 2016, A&A, 585, A134 [Google Scholar]
  41. Dorn, C., Khan, A., Heng, K., et al. 2015, A&A, 577, A83 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  42. Erkaev, N. V., Kulikov, Y. N., Lammer, H., et al. 2007, A&A, 472, 329 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  43. Espinoza, N., & Jordán, A. 2015, MNRAS, 450, 1879 [Google Scholar]
  44. Espinoza, N., Kossakowski, D., & Brahm, R. 2019, MNRAS, 490, 2262 [Google Scholar]
  45. Fletcher, L. N., Gustafsson, M., & Orton, G. S. 2018, ApJS, 235, 24 [NASA ADS] [CrossRef] [Google Scholar]
  46. Flock, M., Turner, N. J., Mulders, G. D., et al. 2019, A&A, 630, A147 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  47. Foreman-Mackey, D. 2018, Res. Notes Am. Astron. Soc., 2, 31 [NASA ADS] [Google Scholar]
  48. Foreman-Mackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, PASP, 125, 306 [Google Scholar]
  49. Foreman-Mackey, D., Agol, E., Ambikasaran, S., & Angus, R. 2017, AJ, 154, 220 [Google Scholar]
  50. Furlan, E., Ciardi, D. R., Everett, M. E., et al. 2017, AJ, 153, 71 [NASA ADS] [CrossRef] [Google Scholar]
  51. Gaia Collaboration (Prusti, T., et al.) 2016, A&A, 595, A1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  52. Gaia Collaboration (Brown, A. G. A., et al.) 2021, A&A, 649, A1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  53. Gaia Collaboration, (Vallenari, A., et al.) 2022, A&A, in press https://doi.org/10.1051/0004-6361/202243940 [Google Scholar]
  54. Gardner, J. P., Mather, J. C., Clampin, M., et al. 2006, Space Sci. Rev., 123, 485 [Google Scholar]
  55. Gavel, D., Kupke, R., Dillon, D., et al. 2014, SPIE Conf. Ser., 9148, 914805 [NASA ADS] [Google Scholar]
  56. Gillon, M., Jehin, E., Lederer, S. M., et al. 2016, Nature, 533, 221 [Google Scholar]
  57. Grimm, S. L., Demory, B.-O., Gillon, M., et al. 2018, A&A, 613, A68 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  58. Hadden, S., Barclay, T., Payne, M. J., & Holman, M. J. 2019, AJ, 158, 146 [NASA ADS] [CrossRef] [Google Scholar]
  59. Harris, C. R., Millman, K. J., van der Walt, S. J., et al. 2020, Nature, 585, 357 [Google Scholar]
  60. Høg, E., Fabricius, C., Makarov, V. V., et al. 2000, A&A, 355, L27 [Google Scholar]
  61. Howell, S. B., Everett, M. E., Sherry, W., Horch, E., & Ciardi, D. R. 2011, AJ, 142, 19 [Google Scholar]
  62. Howell, S. B., Everett, M. E., Horch, E. P., et al. 2016, ApJ, 829, L2 [NASA ADS] [CrossRef] [Google Scholar]
  63. Hu, R., Seager, S., & Yung, Y. L. 2015, ApJ, 807, 8 [Google Scholar]
  64. Hunter, J. D. 2007, Comput. Sci. Eng., 9, 90 [Google Scholar]
  65. Ida, S., & Lin, D. N. C. 2008, ApJ, 685, 584 [Google Scholar]
  66. Irwin, J., Irwin, M., Aigrain, S., et al. 2007, MNRAS, 375, 1449 [Google Scholar]
  67. Irwin, J. M., Berta-Thompson, Z. K., Charbonneau, D., et al. 2015, in Cambridge Workshop on Cool Stars, Stellar Systems, and the Sun (Berlin: Springer), 18, 767 [NASA ADS] [Google Scholar]
  68. Ito, Y., & Ikoma, M. 2021, MNRAS, 502, 750 [NASA ADS] [CrossRef] [Google Scholar]
  69. Jenkins, J. M., Twicken, J. D., McCauliff, S., et al. 2016, SPIE Conf. Ser., 9913, 99133E [Google Scholar]
  70. Jenkins, J. S., Díaz, M. R., Kurtovic, N. T., et al. 2020, Nat. Astron., 4, 1148 [Google Scholar]
  71. Jensen, E. 2013, Astrophysics Source Code Library [record ascl:1306.007] [Google Scholar]
  72. Juric, M., & Tremaine, S. 2008, ApJ, 686, 603 [NASA ADS] [CrossRef] [Google Scholar]
  73. Kempton, E. M. R., Bean, J. L., Louie, D. R., et al. 2018, PASP, 130, 114401 [Google Scholar]
  74. Kipping, D. M. 2013, MNRAS, 434, L51 [Google Scholar]
  75. Kley, W., & Dirksen, G. 2006, A&A, 447, 369 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  76. Kreidberg, L. 2015, PASP, 127, 1161 [Google Scholar]
  77. Kupke, R., Gavel, D., Roskosi, C., et al. 2012, SPIE Conf. Ser., 8447, 84473G [NASA ADS] [Google Scholar]
  78. Kurucz, R. L. 1993, SYNTHE Spectrum Synthesis Programs and Line Data (Cambridge, Mass.: Smithsonian Astrophysical Observatory) [Google Scholar]
  79. Lecavelier Des Etangs, A. 2007, A&A, 461, 1185 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  80. Lillo-Box, J., Barrado, D., & Bouy, H. 2012, A&A, 546, A10 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  81. Lillo-Box, J., Barrado, D., & Bouy, H. 2014, A&A, 566, A103 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  82. Lillo-Box, J., Santos, N. C., Santerne, A., et al. 2022, A&A, 667, A102 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  83. Lindegren, L., Klioner, S. A., Hernández, J., et al. 2021, A&A, 649, A2 [EDP Sciences] [Google Scholar]
  84. Lopez, E. D., & Fortney, J. J. 2014, ApJ, 792, 1 [Google Scholar]
  85. Luger, R., Agol, E., Foreman-Mackey, D., et al. 2019, AJ, 157, 64 [Google Scholar]
  86. Mann, A. W., & von Braun, K. 2015, PASP, 127, 102 [NASA ADS] [CrossRef] [Google Scholar]
  87. Mann, A. W., Johnson, M. C., Vanderburg, A., et al. 2020, AJ, 160, 179 [Google Scholar]
  88. Marcy, G. W., & Butler, R. P. 1992, PASP, 104, 270 [Google Scholar]
  89. Maxted, P. F. L., Anderson, D. R., Collier Cameron, A., et al. 2011, PASP, 123, 547 [NASA ADS] [CrossRef] [Google Scholar]
  90. Mayor, M., Pepe, F., Queloz, D., et al. 2003, The Messenger, 114, 20 [NASA ADS] [Google Scholar]
  91. Mazeh, T., Holczer, T., & Faigler, S. 2016, A&A, 589, A75 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  92. McGurk, R., Rockosi, C., Gavel, D., et al. 2014, SPIE Conf. Ser., 9148, 91483A [NASA ADS] [Google Scholar]
  93. McKay, A. J., DiSanti, M. A., Kelley, M. S. P., et al. 2019, AJ, 158, 128 [NASA ADS] [CrossRef] [Google Scholar]
  94. Middelkoop, F. 1982, A&A, 107, 31 [NASA ADS] [Google Scholar]
  95. Morello, G., Zingales, T., Martin-Lagarde, M., Gastaud, R., & Lagage, P.-O. 2021, AJ, 161, 174 [NASA ADS] [CrossRef] [Google Scholar]
  96. Mousis, O., Deleuil, M., Aguichine, A., et al. 2020, ApJ, 896, L22 [NASA ADS] [CrossRef] [Google Scholar]
  97. Narita, N., Fukui, A., Yamamuro, T., et al. 2020, SPIE Conf. Ser., 11447, 114475K [NASA ADS] [Google Scholar]
  98. Noyes, R. W., Weiss, N. O., & Vaughan, A. H. 1984, ApJ, 287, 769 [Google Scholar]
  99. Owen, J. E., & Adams, F. C. 2019, MNRAS, 490, 15 [Google Scholar]
  100. Owen, J. E., & Alvarez, M. A. 2016, ApJ, 816, 34 [Google Scholar]
  101. Owen, J. E., & Jackson, A. P. 2012, MNRAS, 425, 2931 [Google Scholar]
  102. Owen, J. E., & Lai, D. 2018, MNRAS, 479, 5012 [Google Scholar]
  103. Owen, J. E., & Wu, Y. 2013, ApJ, 775, 105 [Google Scholar]
  104. Papaloizou, J. C. B., Nelson, R. P., & Masset, F. 2001, A&A, 366, 263 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  105. Pepe, F., Mayor, M., Galland, F., et al. 2002, A&A, 388, 632 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  106. Perrakis, K., Ntzoufras, I., & Tsionas, E. G. 2014, Comput. Stat. Data Anal., 77, 54 [Google Scholar]
  107. Piaulet, C., Benneke, B., Rubenzahl, R. A., et al. 2021, AJ, 161, 70 [NASA ADS] [CrossRef] [Google Scholar]
  108. Pollacco, D. L., Skillen, I., Collier Cameron, A., et al. 2006, PASP, 118, 1407 [NASA ADS] [CrossRef] [Google Scholar]
  109. Queloz, D., Mayor, M., Udry, S., et al. 2001, The Messenger, 105, 1 [NASA ADS] [Google Scholar]
  110. Raymond, S. N., Barnes, R., Veras, D., et al. 2009, ApJ, 696, L98 [NASA ADS] [CrossRef] [Google Scholar]
  111. Ricker, G. R., Winn, J. N., Vanderspek, R., et al. 2014, SPIE Conf. Ser., 9143, 20 [Google Scholar]
  112. Santos, N. C., Sousa, S. G., Mortier, A., et al. 2013, A&A, 556, A150 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  113. Santos, N. C., Adibekyan, V., Mordasini, C., et al. 2015, A&A, 580, L13 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  114. Santos, N. C., Adibekyan, V., Dorn, C., et al. 2017, A&A, 608, A94 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  115. Sanz-Forcada, J., Micela, G., Ribas, I., et al. 2011, A&A, 532, A6 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  116. Savel, A. B., Dressing, C. D., Hirsch, L. A., et al. 2020, AJ, 160, 287 [NASA ADS] [CrossRef] [Google Scholar]
  117. Scott, N. J., Howell, S. B., Gnilka, C. L., et al. 2021, Front. Astron. Space Sci., 8, 138 [NASA ADS] [CrossRef] [Google Scholar]
  118. Smith, J. C., Stumpe, M. C., Van Cleve, J. E., et al. 2012, PASP, 124, 1000 [Google Scholar]
  119. Sneden, C. A. 1973, PhD thesis, The University of Texas at Austin, USA [Google Scholar]
  120. Sotin, C., Grasset, O., & Mocquet, A. 2007, Icarus, 191, 337 [Google Scholar]
  121. Sousa, S. G. 2014, ARES + MOOG: A Practical Overview of an Equivalent Width (EW) Method to Derive Stellar Parameters (Berlin: Springer), 297 [Google Scholar]
  122. Sousa, S. G., Santos, N. C., Israelian, G., Mayor, M., & Monteiro, M. J. P. F. G. 2007, A&A, 469, 783 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  123. Sousa, S. G., Santos, N. C., Adibekyan, V., Delgado-Mena, E., & Israelian, G. 2015, A&A, 577, A67 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  124. Sousa, S. G., Adibekyan, V., Delgado-Mena, E., et al. 2021, A&A 656, A53 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  125. Speagle, J. S. 2020, MNRAS, 493, 3132 [Google Scholar]
  126. Stassun, K. G., & Torres, G. 2016, AJ, 152, 180 [Google Scholar]
  127. Stassun, K. G., & Torres, G. 2021, ApJ, 907, L33 [NASA ADS] [CrossRef] [Google Scholar]
  128. Stassun, K. G., Collins, K. A., & Gaudi, B. S. 2017, AJ, 153, 136 [Google Scholar]
  129. Stassun, K. G., Corsaro, E., Pepper, J. A., & Gaudi, B. S. 2018, AJ, 155, 22 [Google Scholar]
  130. Stumpe, M. C., Smith, J. C., Van Cleve, J. E., et al. 2012, PASP, 124, 985 [Google Scholar]
  131. Stumpe, M. C., Smith, J. C., Catanzarite, J. H., et al. 2014, PASP, 126, 100 [Google Scholar]
  132. Suárez Mascareño, A., Rebolo, R., & González Hernández, J. I. 2016, A&A, 595, A12 [Google Scholar]
  133. Szabó, G. M., & Kiss, L. L. 2011, ApJ, 727, L44 [Google Scholar]
  134. Torres, G. 2010, AJ, 140, 1158 [Google Scholar]
  135. Trifonov, T., Kürster, M., Zechmeister, M., et al. 2018, A&A, 609, A117 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  136. Trotta, R. 2008, Contemp. Phys., 49, 71 [Google Scholar]
  137. Tsantaki, M., Sousa, S. G., Adibekyan, V. Z., et al. 2013, A&A, 555, A150 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  138. Vaughan, A. H., Preston, G. W., & Wilson, O. C. 1978, PASP, 90, 267 [Google Scholar]
  139. Virtanen, P., Gommers, R., Oliphant, T. E., et al. 2020, Nat. Methods, 17, 261 [Google Scholar]
  140. Zeng, L., Jacobsen, S. B., Sasselov, D. D., et al. 2019, Proc. Natl. Acad. Sci., 116, 9723 [Google Scholar]
  141. Zhu, W., & Wu, Y. 2018, AJ, 156, 92 [Google Scholar]
2

This code is publicly available in Github through the following URL: https://github.com/jlillo/tpfplotter

4

The last version of ARES code (ARES v2) can be downloaded at https://github.com/sousasag/ARES

7

https://github.com/exord/bayev. A python implementation by R. Diaz of the formalism explained in Perrakis et al. (2014).

8

Juliet core tools used in this analysis are: batman (Kreidberg 2015), starry (Luger et al. 2019), dynesty (Speagle 2020), and Gaussian Processes with george (Ambikasaran et al. 2015) and celerite 1 (Foreman-Mackey et al. 2017).

10

Irradiance temperature is the term used in atmospheric physics to name the equilibrium temperature at zero albedo.

11

Adibekyan et al. (2021) recently found a relation between the composition of low-mass planets and their host stars, but the relation was not one-to-one.

12

The critical point of water is located at 220 bar, this is the pressure at which water transitions from vapor to supercritical. 300 bar is close enough to the pressure of the critical point to prevent the atmospheric model from taking over pressures where the opacity is very high.

All Tables

Table 1

General and dynamical properties of TOI-969 obtained in this work.

In the text
Table 2

Stellar physical properties of TOI-969.

In the text
Table 3

Summary of the posterior distributions of the main physical and orbital properties of TOI-969 b and TOI-969 c from the joint analysis (Sect. 4.3).

In the text
Table 4

Transit times of TOI-969 b estimated from the TESS and LCOGT light curves.

In the text
Table 5

Median and 1σ confidence intervals of the interior and atmosphere modeling.

In the text
Table A.1

TESS light curves of TOI-969 (see Sect. 2.1.1).

In the text
Table A.2

Ground-based photometric time series of TOI-969 from the LCOGT and MEarth facilities (see Sects. 2.1.2 and 2.1.3).

In the text
Table A.3

Photometric time series from WASP (see Sect. 2.1.4).

In the text
Table A.4

Radial velocity and activity indicators for all instruments.

In the text
Table A.5

Inferred and derived parameters from the light-curve analysis (Sect. 4.1) for the one-planet model.

In the text
Table A.6

Inferred and derived parameters from the RV-only analysis (Sect. 4.2) for the preferred model (labeled “2p1cF”) with the inner planet in a circular orbit and the external component in an eccentric orbit.

In the text
Table A.7

Inferred and derived parameters from the joint radial velocity and light-curves analysis (Sect. 4.3) for the preferred model.

In the text

All Figures

thumbnail Fig. 1

Target pixel file corresponding to TOI-969 observations by the TESS mission in sectors 7 (upper panel) and 34 (lower panel). In both panels, the sources from Gaia DR3 are shown as red circles, with their size corresponding to the magnitude contrast against TOI-969 (marked with the label “1” and a white cross). The aperture used by the TESS-SPOC pipeline is shown as a shaded red region in each panel. The figures are computed using the tpfplotter algorithm by Aller et al. (2020).
In the text
thumbnail Fig. 2

TESS photometric time series from sectors #7 (left panels) and #34 (right panels). Upper panels display the PDCSAP flux and a simple median filter computed trend (blue solid line). The detrended light curves are shown in the lower panels and the transits from TOI-969 b are marked by vertical red lines at the bottom of each panel.
In the text
thumbnail Fig. 3

High-spatial-resolution contrast curves for TOI-969 from the Zorro speckle instrument on Gemini South (blue curve for the 562 nm band and green curve for the 832 nm band), and the 3 m Shane telescope in the J-band (orange) and Ks-band (red). The two dotted lines correspond to 10% (upper) and 5% (lower) contamination levels in the light curve. The dashed line corresponds to the maximum contrast that a blended binary could have to mimic the transit depth of TOI-969 b.
In the text
thumbnail Fig. 4

Radial velocity time series and model inference from the joint radial velocity and light curve analysis in Sect. 4.3, including the GP noise (blue dashed line) and full RV model (black solid line for the median posterior model and gray-shaded regions for 68.7 (dark gray) and 99.7% (light gray) confidence intervals).
In the text
thumbnail Fig. 5

Periodogram of the radial velocity and activity indicator time series of the HARPS dataset. Panel a corresponds to the original RV time series, while panel b corresponds to the time series after a simple linear detrending. In panel c we include an additional linear detrend-ing corresponding to the activity indicator most correlated with the RVs (in this case the BIS-span). Panels d–e include the periodogram of the FWHM of the cross-correlation function and the bisector span, respectively. Finally, the bottom panel shows the window function corresponding to the HARPS sampling. The periodicity corresponding to TOI-969 b and TOI-969 c are shown as vertical dotted lines in all panels. The stellar rotation period and its first two harmonics are shown as dashed vertical lines.
In the text
thumbnail Fig. 6

Spectral energy distribution of TOI-969. Red symbols represent the observed photometric measurements, where the horizontal bars represent the effective width of the passband. Blue symbols are the model fluxes from the best-fit NextGen atmosphere model (black).
In the text
thumbnail Fig. 7

Phase-folded light curves from TOI-969 after removing the GP for the TESS data and the linear dependency for the LCOGT dataset. Left panels show the individual phase-folded light curves per instrument and the corresponding joint model. Black symbols correspond to 9-min bins. Right panels show the residuals from the joint model.
In the text
thumbnail Fig. 8

Comparison of the logarithm of the Bayesian evidence (ℬ = In 𝒵) for the different tested models (labeled in the X axis) in the case of the RV-only analysis (Sect. 4.2). The labels for the models correspond to Xp[Pic][Z], where X is the number of planets assumed in the system, Pic are the identifiers of the planets with assumed circular orbits, and Z indicates whether we are using a flat (F) slope (i.e., assuming slope = 0) or a quadratic (Q) trend. The gray-shaded region displays the region of Δ𝒵 = −6 from the larger evidence model showing its strength against simpler models.
In the text
thumbnail Fig. 9

Phase-folded radial velocity data from TOI-969 b (left) and TOI-969 c (right) after removing the signal from the other component in the system and the GP noise for the selected model (labeled “2plcF") from the radial velocity analysis (see Sect. 4.2). For TOI-969 b, we only include HARPS and PFS data, while removing CORALIE for improved clarity given the larger uncertainties for that dataset. The gray-shaded regions correspond to the region covered by models within 68.7% (dark gray) and 95% (light gray) of the confidence interval.
In the text
thumbnail Fig. 10

Observed–calculated diagram of the transit times of TOI-969 b. LCOGT transit corresponds to transit number #400. The gray regions represent the 1σ uncertainties of the ephemeris reported in Table A.7.
In the text
thumbnail Fig. 11

Comparative properties of TOI-969b. Left: mass–radius diagram of known exoplanets with measured mass precisions better than 30% (gray symbols) and the location of TOI-969b. The bulk density lines corresponding to different compositions from Zeng et al. (2019) are shown as solid traces, and the dashed lines correspond to iso-densities of 0.3, 1, and 3 gm−3. Right: period-radius diagram of known exoplanets (gray crosses) highlighting those with measured masses (open black circles). The light-gray-shaded region represents the mean boundaries of the hot-Neptune desert derived by Mazeh et al. (2016) and the dark-grey-shaded region corresponds to the interior envelope of the \cr boundaries. The location of TOI-969 b is marked in red, as are the locations two other planets in the hot-Neptune desert, namely TOI-849 b (Armstrong et al. 2020) and LTT 9779 b (Jenkins et al. 2020).
In the text
thumbnail Fig. 12

Ternary diagram with the 1σ confidence region of TOI-969b in scenario 1 (black), where only the mass and the radius are considered as input data, and in scenario 2 (red), where the stellar abundances are also included as input data in the MCMC interior structure analysis. The mantle mass fraction (MMF) is defined as MMF = 1 - CMF - WMF. The green dot and red square indicate the position of Earth and Mercury in the ternary diagram, respectively.
In the text
thumbnail Fig. 13

Mass-radius relationships for supercritical water (SW) planets (Acuña et al. 2021; Mousis et al. 2020), planets with Earth-like cores and H/He atmospheres (Zeng et al. 2019), and dry planets with different Fe contents (Brugger et al. 2017). We assume high chemical equilibrium temperatures for supercritical and H/He volatile layers of 1200 K and 2000 K, respectively. The position of TOI-969 b is indicated as a magenta star.
In the text
thumbnail Fig. 14

Synthetic transmission spectra of TOI-969 b. Fiducial models assuming a cloud-free atmosphere with solar abundances (bottom), cloud-free atmosphere with enhanced metallicity by a factor of 100 (middle), and atmosphere with solar abundances and an optically thick cloud deck with top pressure of 100 Pa (top). Estimated uncertainties are shown for the observation of one transit with JWST NIRISS-SOSS, NIRSpec-G395M, and MIRI-LRS configurations. Simulated data are color coded according to wavelength.
In the text
thumbnail Fig. 15

Schematic view of the planetary system TOI-969, with the Earth direction toward the bottom of the panel. A total of 500 randomly selected orbits from the posterior distribution of the joint model (see Sect. 4.3) are drawn for TOI-969 b (orange) and TOI-969 c (blue). The orbit corresponding to the median solution for each parameter is displayed as a thicker solid line of the corresponding color for each planet.
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.