EDP Sciences logo
Subscriber Authentication Point
Search
Menu
Home All issues Volume 687 (July 2024) A&A, 687 (2024) A48 Full HTML
Open Access
Issue
A&A
Volume 687, July 2024
Article Number A48
Number of page(s) 27
Section Planets and planetary systems
DOI https://doi.org/10.1051/0004-6361/202347981
Published online 28 June 2024

© The Authors 2024

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

The Transiting Exoplanet Survey Satellite (TESS, Ricker et al. 2015) has already added about 400 confirmed planets to the known sample1. More than half of these new objects are smaller than Neptune, and transit stars that are bright enough for a detailed characterization with high-precision spectroscopy. Additionally, about 70 of these planets are hosted by M dwarfs. While the Kepler mission revealed a large abundance of such planets with no equivalent in the Solar System (Borucki et al. 2011; Dressing & Charbonneau 2013), most of the Kepler planets are currently out of reach for a detailed characterization due to the faintness of their host star.

The study of the bulk composition of these small planets shows two distinct populations with radii smaller than that of Neptune and bigger than that of the Earth. First, “super-Earths” are thought to be rocky planets with either thin or no atmosphere, while the other “mini-Neptunes” show smaller densities consistent with an extended atmosphere or a significant water fraction (e.g., Rogers 2015; Adams et al. 2008). The formation pathways of these planets are not fully understood (for a comprehensive review, see Bean et al. 2021), especially the paucity of planets found between ~1.5 and ~2.5 R⊕ (Fulton et al. 2017) for FGK stars. Several theories such as gas-poor formation (Lee et al. 2014), and atmospheric loss either by photoevaporation (Owen & Wu 2013) or core-powered mass loss (Gupta & Schlichting 2019) could explain this so-called radius valley, which may drift toward smaller radii for low-mass stars (Berger et al. 2020; Cloutier & Menou 2020), depending on the model. Another approach considers a density valley rather than a radius valley for M-type stars (Luque & Pallé 2022). The study distinguishes between three populations: rocky planets, water worlds, and puffy mini-Neptunes, with a common formation history for these last two. To illuminate these different formation theories, an increase in the sample of thoroughly characterized sub-Neptune-sized planets is critical. Small and cool M dwarfs are particularly interesting targets for transmission spectroscopy as they present high signal-to-noise ratios (S/Ns) even for smaller transiting planets. Extended atmospheres, such as those expected for mini-Neptunes, should further increase the S/N for transmission spectroscopy.

TOI-4336 A is part of a hierarchical triple M-dwarf system (M3.5-, M3.5-, and M4-type stars) located at 22 pc from the Sun. The host star is the brightest component of the inner binary pair that has a minimum orbital separation of over a hundred au, and is at a projected angular separation of 6.25″. We adopted them as TOI-4336 A and B (TIC 166184428 and TIC 166184426, respectively) according to their brightness. The third star, which we refer to as TOI-4336 C (TIC 166184390), is the latest of the system and is on a wider orbit at a distance of over 2900 au, and at an angular separation of 98.44″. Here we report the detection of TOI-4336 A b, a planet that lies at the inner edge of the empirical habitable zone (HZ) (Kopparapu et al. 2013, 2017; Kaltenegger 2017; Zsom et al. 2013) of the system. It receives less irradiation than a young Venus; however, since M-dwarf irradiation warms planets more effectively than Sun-like stars (see, e.g., Kasting et al. 1993), that flux moves TOI-4336 A b just closer to the star than the inner edge of the empirical HZ, making the planet an intriguing example of a planet near a boundary of the HZ.

We describe in Sect. 2 the observations and methods used to characterize the triple star system. The TESS and ground-based observations used to validate the planetary nature of TOI-4336 A b are detailed in Sect. 3, and the statistical validation is reported in Sect. 4. The joint analysis of all the photometric data is presented in Sect. 5 and the results are discussed in Sect. 6. Finally, we present our conclusions in Sect. 7.

2 Stellar characterization

2.1 Spectroscopic reconnaissance

We gathered near-infrared spectra of TOI-4336 A and the two resolved, co-moving stars (TOI-4336 B and TOI-4336 C) with the SpeX spectrograph (Rayner et al. 2003) on the 3.2-m NASA Infrared Telescope Facility (IRTF) on 2021 Jun 27 (UT) and again on 2021 Jun 28 (UT). We used the short-wavelength cross-dispersed (SXD) mode of SpeX and the 0.″5 × 15″ slit aligned to the parallactic angle, which yielded spectra covering 0.75–2.42 μm at R~1200 (Fig. 1). We collected 6 exposures of each target, nodding in an ABBA pattern. Our total exposure times were 360 s each for TOI-4336 A and TOI-4336 B, and 540 s for TIC 166184390 TOI-4336 C. We collected the standard set of SXD flat-field and arc-lamp exposures immediately after the science frames, followed by the A0 V standard HD 130163. We reduced the data with SPEXTOOL v4.1 (Cushing et al. 2004). We used SPLAT to compare the spectra to standards in the IRTF Spectral Library (Cushing et al. 2005; Rayner et al. 2009), focusing on the 0.9–1.4 μm region for the spectral classification (Kirkpatrick et al. 2010), and to estimate metallicity [Fe/H] via the Mann et al. (2013) relation (see Delrez et al. 2022, for details). Between nights, the spectral-type determinations of each target are identical and their metallicity estimates are consistent at <1σ. We estimate spectral types of M3.5 ± 0.5 for TOI-4336 A and TOI-4336 B and M4.0 ± 0.5 for TOI-4336 C. Combining measurements from both nights, we obtain final [Fe/H] estimates of −0.20 ± 0.12, −0.21 ± 0.12, and −0.17 ± 0.12 for TOI-4336 A, B, and C, respectively (see Table 1).

We acquired an optical spectrum of TOI-4336 A on 2022 Jan 07 (UT) using the Low Dispersion Survey Spectrograph (LDSS3-C, Stevenson et al. 2016) on the 6.5-m Magellan II (Clay) Telescope under clear and stable conditions. We used LDSS-3C in long-slit mode with the standard setup (fast readout speed, low gain, 1 × 1 binning) and the VPH-Red grism, OG-590 blocking filter, and the 0.″75 × 4′ center slit, which provides spectra covering 6000–10 000 Å at R ~ 1810. We collected eight, 60-s exposures of the target, followed by a 1-s arc-lamp exposure and three 10-s flat field exposures with the quartz high lamp and reduced the data with a custom Python-based pipeline (Dransfield et al. 2024). We used the ratio of the spectrum of the G2 V star HR 5325, observed at a similar airmass, to a G2V template from Pickles (1998) to compute a relative flux calibration of the TOI-4336 spectrum. No correction was made to address telluric absorption. The reduced spectrum is shown in Fig. 2. It was then analyzed using tools in the kastredux package2. Comparison with the Sloan Digital Sky Survey templates from Kesseli et al. (2017) shows an excellent match to an M4 dwarf template. This is confirmed by analysis of spectral classification indices from Reid et al. (1995); Lépine et al. (2003); and Riddick et al. (2007). We also computed the ζ metallicity index (Lépine et al. 2007, 2013), determining a value of 1.026 ± 0.002, corresponding to a metallicity of [Fe/H] = +0.04+0.20 based on the empirical calibration of Mann et al. (2013).

We obtained high-resolution spectroscopic observations with the CHIRON spectrograph, located on the SMARTS 1.5-meter telescope at Cerro Tololo Inter-American Observatory, (Tokovinin et al. 2013) for both TOI-4336 A and the co-moving companion TOI-4336 B. We used CHIRON’s ‘slicer’ mode, which employs an image slicer to achieve a resolving power of R ~ 80 000 from 4100 to 8700 Å. Our observations yield per-pixel S/Ns that range from 10–12 in the TiO bands around 7100 Å. Spectra were extracted using the official CHIRON pipeline (Paredes et al. 2021), and we derived radial velocities (RVs) and spectral line profiles as described by Pass et al. (2023)3. This reduction is specifically optimized for mid-to-late M dwarfs and produces carefully calibrated relative RVs. For our relative RVs, the error budget is dominated by spectrograph stability. For our absolute RVs, the dominant uncertainty is a 0.5 km s−1 error in the RV scale; this error stems from the absolute RV uncertainty in a comparison spectrum of Barnard’s Star, derived using 17 yr of measurements from the CfA Digital Speedometer. The RV measurements are given in Table 1. TOI-4336 A was observed on 2021 Feb. 06 and 2021 Jul. 12, and TOI-4336 B was observed on 2021 Feb. 06.

We find no evidence for significant velocity variation in TOI-4336 A, and the difference in absolute radial velocity between the two stars is consistent with orbital motion in a wide binary given their separation. We detect no significant rotational broadening for either star, setting upper limits of v sin i* < 1.9 km s−1. This limit corresponds to half a resolution element of the CHIRON spectrograph. Hα is seen in absorption for both stars, and following Newton et al. (2017), we measure EWHα = 0.1644 ± 0.0025Å for TOI-4336 A, which places it among the sequence of quiescent M dwarfs. Taken together, we interpret the lack of detectable activity and rotational broadening as an indication that the star is not young, which is in agreement with the results of the LDSS3 spectral analysis, this is discussed in Sect. 2.3.

2.2 Spectral energy distribution

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 DR3 parallax (with no systematic offset applied; see, e.g., Stassun & Torres 2021), following the procedures described in Stassun & Torres (2016); Stassun et al. (2017, 2018). We pulled the JHKS magnitudes from 2MASS (Cutri et al. 2003), the W1–W4 magnitudes from WISE (Cutri et al. 2021), the GBP and GRP magnitudes from Gaia (Gaia Collaboration 2022), and the NUV flux from GALEX (Bianchi et al. 2017). We also used the Gaia spectrophotometry spanning 0.4–1.0 μm, providing an especially strong constraint on the overall absolute flux calibration. Together, the available photometry spans the stellar SED over the wavelength range 0.2–22 μm (see Fig. 3).

We performed a fit using PHOENIX stellar atmosphere models (Husser et al. 2013), with the effective temperature (Teff) and metallicity ([Fe/H]) as free parameters (the surface gravity, log g, has very little influence on the broadband SED). We set the extinction, AV, to zero given the close proximity of the system. The resulting fit (Fig. 3) has a reduced χ2 of 1.3 (excluding the NUV measurement, which suggests some chromospheric activity), with best-fit Teff = 3300 ± 75 K and [Fe/H] = 0.0 ± 0.2. The derived parameters are shown in Table 1.

Integrating the model SED gives the bolometric flux at Earth, Fbol = 7.40 ± 0.17 × 10−10 erg s−1 cm−2. Taking the Fbol together with the Gaia parallax directly gives the luminosity, Lbol = 0.01163 ± 0.00027 L⊙. Similarly, the Fbol together with the Teff and the parallax gives the stellar radius, R* = 0.330 ± 0.015 R⊙. The stellar mass can also be estimated via the empirical MK based relations of Mann et al. (2019), giving M* = 0.331 ± 0.010 M⊙. All uncertainties are propagated in the usual manner, except for the MK-based mass estimate for which we adopt the systematic uncertainties quoted for the Mann et al. (2019) relations as the dominant source of error. These parameters are summarized in Table 2.

For completeness, we applied the same SED-fitting procedures to the other two stars in the system, with the results shown in Fig. 3 and summarized in Table 1. We placed the three stars of the TOI-4336 system in a color-magnitude diagram to compare their properties to nearby M dwarfs (see Fig. A.1). TOI-4336 C appears less luminous than the other two, which is consistent with their spectral types.

thumbnail Fig. 1

SpeX spectra of TOI-4336 (TOI-4336 A, top row), TIC 166184426 (TOI-4336 B, middle row), and TIC 166184390 (TOI-4336 C, bottom row) from 2021 Jun. 27 (left column) and 2021 Jun. 28 (right column). The target spectrum (blue) is shown alongside the best-fit spectral template from the IRTF Spectral Library (gray). All spectra are normalized to their flux in the 0.9-1.4 μπι region. Wavelengths with strong telluric absorption are shaded (and largely masked from the spectra), and prominent M-dwarf absorption features are highlighted.

2.3 Age estimation

To estimate the age of the system, we compared its kinematics and metallicity to a large sample of stars from GALAH survey DR3 (Buder et al. 2018). Ages for stars in the GALAH survey are estimated based on the Bayesian method by Sharma et al. (2018), which is designed to fit the measured distribution of stellar parameters (sky positions, Teff, log g, [Fe/H], magnitudes) with a combination of model isochrones and a Galactic population synthesis model. We used the same reference frame as the GALAH DR3 catalog and the Gaia radial velocities to compute the UVW velocities of the system, given in Table 1. We selected objects with distances < 300 pc, and total 3D-velocity |V – Vsystem| < 10 km s−1 and metallicities within 1-standard deviation of the system, our sample is shown in Fig. B.1. These constraints yield an age of 6.7−3+2.7$6.7_{ - 3}^{ + 2.7}$ Gyr. Additionally, the lack of Hα (6563 Å) and Li I (6708 Å) emission in any of the optical spectra further confirms that the system is likely old (Newton et al. 2017; Kiman et al. 2021). In particular, examining the LDSS3 spectra, we find no evidence of significant Hα emission with a 3σ equivalent width limit of 0.12 Å, suggesting an age ≳4.5 Gyr according to the age–activity relation of West et al. (2008). However, TOI-4336 A has a mass placing it near the fully convective boundary, and more recent studies have found the average age of transition to the inactive mode is 2.4 ± 0.3 Gyr for fully convective M dwarfs (Medina et al. 2022) with some of them transitioning as early as 600 Myr (Pass et al. 2022).

Finally, using the latest version of the BANYAN tool (Gagné et al. 2018)4, we find that the system has a 99.9 percent probability of being a field object, with no association with nearby young moving groups. This further confirms that the triple system is likely old.

Table 1

Properties of the TOI-4336 system.

thumbnail Fig. 2

LDSS3 red optical spectrum of TOI-4336 A (black line), compared to its best-fit M4 template (Kesseli et al. 2017, magenta line). Spectra are normalized in the 7400–7500 Å region, and major absorption features are labeled, including regions of strong telluric absorption (⊕). The inset box shows a close-up of the region encompassing Hα (6563 Å) and Li I (6708 Å) features, neither of which is detected.

2.4 Triple M-dwarf system

TOI-4336 A is the primary star of a resolved triple system, with angular separations of 6.25″ and 98.44″ from the Β and C components respectively, identified as WDS J13444-4020 in the Washington Double Stars Catalog (WDS, Mason et al. 2001). We used the astrometric measurements from WDS, the stellar masses from SED analysis and the Gaia astrometry, distances and proper motions to get a first-order estimation of the orbital elements of the system with LOFTI-Gaia package5 (Pearce et al. 2020). We obtained a posterior sample of 10 000 accepted orbits for each system (see Appendix C). The posterior parameters are given in Table C.2. The results show that the close binary has a median semi-major axis of 133 au with a 68% confidence interval from 72 to 175 au. This results correspond to a median periastron separation of 34.7 au. with a 68% confidence interval from 5 to 118 au. Although the short phase coverage does not allow a detailed orbital characterization of the system, these results are used as a point of reference to evaluate the effect of the triple system configuration on the formation of TOI-4336 A b (see Sect. 6.2).

thumbnail Fig. 3

Spectral energy distribution of TOI-4336 A (top panel) and its companion stars TOI-4336 Β (middle panel) and TOI-4336 C (bottom panel). Red symbols represent the observed photometric measurements, the horizontal bars represent the effective width of the passband. Blue symbols are the model fluxes from the best-fit PHOENIX atmosphere model (black). The Gaia spectrophotometry is represented as a grey swathe; a closeup view is shown in the inset.

3 Photometric observations

3.1 TESS

The two close components of the system, TOI-4336 A and B, have an angular separation of 6″ and thus are located in the same TESS pixel, given the pixel size of 21″. The tertiary component is, however, resolved in TESS. They were observed in Sector 11 (2019 Apr. 26 to 2019 May 20), Sector 38 (2021 Apr. 29 to 2021 May 26) and Sector 64 (2023 Apr. 06 to 2023 May 04), with both short- and long- cadence modes. The 2-min cadence image data were reduced and analyzed using the Science processing Operations Center (SPOC, Jenkins et al. 2016) pipeline at NASA Ames Research Center. 30-min cadence data were obtained by the TESS-SPOC pipeline (Caldwell et al. 2020) and 10-min cadence data were obtained by the Quick Look Pipeline (QLP, Huang et al. 2020).

We extracted the photometry for the global analysis from the TESS target pixel files (TPFs) for all three sectors. We retrieved the 2-min TESS TPFs from the Mikulski Archive for Space Telescopes using the lightkurve6 package (Lightkurve Collaboration et al. 2018) using the hard quality bitmask to remove the images affected by scattered light or sub-optimal attitude control. At the beginning of each orbit of Sector 11, Camera 1 was subject to scattered light, and the attitude control was disabled for a short period. Using the hard-quality bitmask allows removing this affected bit of data (quality flag 7407). We extracted the photometry from the images using a set of custom apertures inscribed within a square of 4 × 4 pixels centered on the target, these are shown in Fig. 4. The TESS TPFs are obtained with tpfplotter7 (Aller et al. 2020). After the removal of the 5σ outliers, we used the cotrending basis vectors (CBVs) obtained by the presearch and data conditioning (PDC) pipeline module of SPOC, a method initially developed to remove low and high-frequency systematic trends in Kepler data (Kinemuchi et al. 2012). We used the multi-scale correction method, which is preferred to preserve transit signals, combined with a correction of short spikes. The correction was calculated using a linear regression approach with an L2 regularization implemented in lightkurve. To prevent under- or over-fitting of the data, the value of the L2 penalty is optimized according to the PDC goodness metrics. We compared the corrected light curves using the combined differential photometric precision (CDPP) metric as implemented in lightkurve and selected the one for which the metric was minimal. The extracted light curves are shown in Fig. 5.

Data Validation Reports were produced by the TESS pipelines each time a new sector became available. The validation tests systematically performed for a planet candidate were passed, although in the first instance the period found was half the true period and the centroid shift analysis did match to TOI-4336 A, TOI-1955.01 was in fact issued for TOI-4336 B. Our ground-based photometric follow-up confirmed the source to be TOI-4336 A. On the TESS side, the SPOC conducted a transit search of Sectors 38 and 64 on 2021 Jul. 2 and 2023 Jun. 20 respectively. The transit search was conducted with an adaptive, noise-compensating matched filter (Jenkins et al. 2002, 2010, 2020), producing a threshold crossing event for which an initial limb-darkened transit model was fitted (Li et al. 2019) and a suite of diagnostic tests were conducted to help make or break the planetary nature of the signal (Twicken et al. 2018). The TESS Science Office reviewed the vetting information and issued an alert for TOI-4336.01 on 2021 Jul. 28 (Guerrero et al. 2021). According to the difference image centroiding tests, the host star is located within 4.0 ± 2.7″ of the transit signal source of TOI-4336 A.

Table 2

Properties of TOI-4336 A and TOI-4336 A b based on our global transit model (see Sect. 5).

thumbnail Fig. 4

TESS target pixel files showing the custom apertures used to extract the photometry of TOI-4336 A b in this work for Sectors 11, 38, and 64.

3.2 Ground-based photometry

We collected photometric data using four sets of facilities as part of the TESS follow-up observing program (TFOP) sub-group 1 (SG1) for seeing-limited photometry: SPECULOOS-Southern Observatory (SSO; Delrez et al. 2018; Sebastian et al. 2021), TRAPPIST-South (Gillon et al. 2011; Jehin et al. 2011), Las Cumbres Observatory Global Telescope (LCOGT; Brown et al. 2013), and ExTrA (Bonfils et al. 2015). The observations are summarized in Table D.1, and the quality of the light curves is quantified by their RMS value, it is given in Table E.1.

All but the ExTrA data reduction and analysis was done using a custom pipeline for image processing and photometric extraction built with the prose8 package (Garcia et al. 2022, 2021). First, the image calibration was performed, and the images aligned using twirl9. We detected the stars using an implementation of DAOPhot from Photutils (Bradley et al. 2023), and used ballet10 as centroiding algorithm. The photometric extraction was done on a set of 30 circular apertures ranging between 0.5 and twice the value of the full width at half maximum (FWHM) of the target’s point-spread function (PSF), and the position of the background annulus was selected on the basis of the lack of contaminant in the vicinity of the target. We performed differential photometry on all data sets to retrieve our light curves. In doing so we selected comparison stars in the field as well as apertures minimizing the white and red noise in the target light curve as calculated by Pont et al. (2006). We treated each observation separately, and the number of comparison stars varied between 3 and 9. Given the proximity of stars A and B, the flux in the aperture of the TRAPPIST-South observations is contaminated as it contains fully the PSFs of the two stars. This is taken into account in the global analysis where we fit the dilution factor (see Sect. 5). For the observations with the 1-m class telescopes (SSO and LCO), the stars are completely resolved with a negligible overlap of the PSFs. In that case, we systematically selected TOI-4336 B as a comparison star.

3.2.1 TRAPPIST-South

TRAPPIST-South is a 0.6m Ritchey-Chrétien telescope located at La Silla Observatory in Chile (Gillon et al. 2011; Jehin et al. 2011). It is equipped with an FLI ProLine CCD camera with a pixel scale of 0.64″, providing a field of view of 22′ × 22′. We observed four transits of TOI-4336 A b with this facility in the Sloan-z′ filter with 15 s exposures on 2021 Apr. 30, 2021 Jun. 18, 2021 Aug. 6, and 2023 May 27.

3.2.2 SPECULOOS-South

The SPECULOOS Southern Observatory, located at ESO Paranal Observatory in Chile, is composed of four 1.0m F/8 Ritchey-Chrétien telescopes (Delrez et al. 2018; Sebastian et al. 2021), named after the Galilean moons Io, Europa, Ganymede, and Callisto. All telescopes are equipped with a deep-depletion Andor iKon-L 2k × 2k CCD camera with a total field of view of 12′ for a pixel scale of 0.35″ (Burdanov et al. 2018). We collected data on the nights of 2021 Jun. 18, 2022 Apr. 8, 2022 May 27, and 2023 Mar. 17 with SSO in the Sloan-r′ filter with 10s exposures and the Sloan-g′ filter with 24s exposures.

thumbnail Fig. 5

TESS 2-min cadence photometry obtained using custom apertures for Sectors 11, 38, and 64. The transits are highlighted in purple.

3.2.3 LCOGT

A total of six transits were obtained with LCOGT in the Sloan-g′ with 150s exposures and Pan-STARRS zs with 45s exposures on the nights of 2021 Jun. 18, 2022Feb. 18, 2022 Apr. 8, and 2023 Mar. 17. The telescopes are equipped with 4096 × 4096 SINISTRO Cameras, having an image scale of0.389″ per pixel, resulting a field of view of 26′ × 26′. The raw data were calibrated by the standard LCOGT BANZAI pipeline (McCully et al. 2018) and photometric measurements were extracted using the prose package, similarly to the TRAPPIST and SPECULOOS data.

3.2.4 ExTrA

The ExTrA facility, located at ESO’s La Silla Observatory, is composed of three 0.6m telescopes feeding a near-infrared multi-object spectrograph. We simultaneously obtained one transit with two ExTrA telescopes (Bonfils et al. 2015) on 2022 May 27 in a bandpass centered on 1.21 μm with 60 s exposures and the spectrograph’s low resolution mode (R ~ 20) and 8″ diameter aperture fibers. The resulting ExTrA data were analyzed with custom data reduction software, detailed in Cointepas et al. (2021).

4 Validation of the planet

4.1 Archival imaging

Investigating archival imaging of the field of view of the considered target is common practice to look for a possible blend with a background object in the current images. The large proper motion of the TOI-4336 system (166 mas yr−1; Gaia Collaboration 2022) makes this possible. We inspected DSS/POSS II images in blue and IR (Lasker et al. 1996; Reid et al. 1991) spanning 47 yr compared to the SSO observation, as shown in Fig. 6. We conclude there is no background star with a limiting magnitude of B~22.5 at the current position of TOI-4336 A which could potentially be the source of the transit events or impact our conclusions.

4.2 High-angular resolution imaging

A search for proper motion blend candidates does not eliminate the possibility that the star itself is an unresolved binary. Close-bound stellar companions can confound exoplanet discovery and parameter derivation. If a close companion does exist, the detected transit signal will yield incorrect stellar and exo-planet parameters (Ciardi et al. 2015; Furlan & Howell 2017). Additionally, the presence of a close companion star could lead to the non-detection of small planets residing within the same exoplanetary system (Lester et al. 2021).

The renormalized unit weight error (RUWE) quantifies the possible presence of an unresolved companion or acceleration during the astrometric solution estimation from Gaia measurements (Gaia Collaboration 2022). As all components of the systems have values above the expected value of ~1.0 for a single star (see Table 1), we use high-resolution imaging to constrain the presence of unresolved companions.

TOI-4336 A was observed three times with the Zorro instrument on the Gemini-South 8-m telescope (Scott et al. 2021; Howell & Furlan 2022): 2022 Mar. 19, 2022 May 17, and 2023 May 27. Zorro provides simultaneous speckle imaging in two bands (562 nm and 832 nm) with output data products including a reconstructed image with robust contrast limits on companion detections (e.g., Howell et al. 2016). All three observations were processed with our standard reduction pipeline (Howell et al. 2011), Fig. 7 shows the final contrast curves and the reconstructed speckle images for the 2022 May 17 observation. TOI-4336 A is an isolated star with no companion brighter than 5–7 magnitudes below that of the target star from the 8-m telescope diffraction limit of 0.2″ out to 1.2″. This excludes the presence of companion stars with spectral types between M4 and early-L at these angular limits. At the distance of the TOI-4336 system (d = 22.5 pc), they correspond to spatial limits of 0.45 to 27 au.

We acquired a second set of speckle imaging for the A and Β components of the TOI-4336 system with the HRCam instrument of the SOAR telescope (Tokovinin 2018), the data were analyzed following the method outlined in Ziegler et al. (2020). The observations were obtained on 2021 Jul. 14, and 2022 Mar. 20 for Β and A respectively, both in the Cousins-I filter, and the contrast curves yielded no companions within 1″ with a contrast of 6.7 and 7 magnitudes. The 5σ sensitivity and speckle autocorrelation functions are shown in Fig. 8.

thumbnail Fig. 6

Archival images of the field around the TOI-4336 system. From left to right: 1974 image taken with the blue plate of DSS/POSS-II, 1993 image taken with the InfraRed plate of DSS/POSS-II, and 2021 image taken with SSO in Sloan-r′. The white circles indicate the position of the stars on the 2021 images.
thumbnail Fig. 7

Contrast curve obtained with Zorro for TOI-4336 A in two bands (562 nm and 832 nm) and the reconstructed speckle image of the observation from 2022 May 17.

4.3 Statistical validation

We made use of the statistical validation package TRICERAT0PS11 (Giacalone et al. 2021; Giacalone & Dressing 2020) using the same procedures detailed in Dransfield et al. (2024); Pozuelos et al. (2023); Barkaoui et al. (2023), to validate the planetary nature of TOI-4336.01. The threshold for validation is set at 0.015, and initially we find the false positive probability (FPP) to be 0.444. However, upon inspection of the probability breakdown, we found that the nearby transiting planet (ΝΤΡ) probability on TOI-4336 Β is 0.444, while probabilities for all other false positive scenarios are of order 10−8 or smaller. The reason for the high FPP is therefore that the two stars are blended within the TESS aperture. While we do fold in ground-based observations to the light curve used for the statistical validation to improve the precision, TRICERATOPS only makes use of TESS apertures for the calculation of scenario probabilities. However, in our ground-based observations, the close binary components of the triple system are resolved and we were able to confirm that the transits are on TOI-4336 A. Therefore, by eliminating this scenario, the FPP for TOI-4336.01 falls to ~10−8, placing it well below the threshold for validation.

5 Global analysis

5.1 Transit analysis

We performed a global analysis of the photometric data using a custom Fortran code described in Gillon et al. (2012); Delrez et al. (2022) in which we made use of emcee (Foreman-Mackey et al. 2013), a Markov Chain Monte Carlo affme invariant ensemble sampler proposed by Goodman & Weare (2010). We used a combination of quadratic limb-darkening transit model (Mandel & Agol 2002) and baseline model to fit the data. The baseline model represents the combination of systematic effects producing correlated (red) noise, such as atmospheric conditions (airmass, FWHM of the PSF, sky background) or instrumental effects (the variation of the position of the star on the detector along the x and y directions). To remove these effects, we select a linear combination of low-order polynomials with respect to these five parameters, in addition to the time as a parameter, which minimizes the Bayesian information criterion (BIC, Schwarz 1978). Fitting simultaneously the transit signal and the correlated noise allows a good propagation of the uncertainties to the derived parameters. The TS light curves are also affected by an additional offset coming from a constraint of the German equatorial mount it is equipped with. As it crosses the meridian, the telescope mount rotates by 180° and the stars fall onto different pixels with varying sensitivity, affecting the flux measurements. The code also produces ßw and ßr, two scaling factors to apply to the photometric errors to account for an under-or over-estimation of white and red noise in each light curve. The baseline models and the error scaling factors are shown in Table E.1.

Given that the TRAPPIST-South and TESS photometric data sets are contaminated by the second fainter star of the system, we computed the dilution factor as the flux ratio in the TESS band following the definition used in the code: Dil=FblendFsource${\rm{Dil}} = {{{F_{{\rm{blend}}}}} \over {{F_{{\rm{source}}}}}}$ with Fblend the flux of the contaminant and Fsource the flux of the target star. We used a normal prior with a conservative 3% error to account for possible faint stars contaminating the apertures. The quadratic limb darkening coefficients u1, u2 are taken from Claret (2018) for TESS and the ground-based observations from Claret et al. (2012), except for the ExTrA 1.2 μm and Sloan-zs filters for which the priors were obtained with PyLDTK (Parviainen & Aigrain 2015) and the PHOENIX model atmospheres (Husser et al. 2013). We used a conservative value of 0.05 for the uncertainty in the normal prior distributions of these parameters. We also assumed a normal prior probability distribution function (PDF) for the effective temperature, and the stellar mass and radius based on the results described in Sect. 2.2 (see Table 1). The jump parameters sampled in our analysis are: the effective temperature, the metallicity, the transit epoch, the log of the orbital period, the transit depth as defined by dF=Rp2/R⋆2${\rm{d}}F = R_p^2/R_ \star ^2$, the cosine of the orbital inclination, the log of the stellar density, and the log of the stellar mass. The limb darkening coefficients are also taken as jump parameters following the parametrization of Kipping (2013) for triangular sampling with the quadratic limb darkening law: q1 = (u1 + u2)2 and q2 = 0.5u1(u1 + u2)−1. Finally, we also fitted for the dilution parameters of the TRAPPIST-South and TESS data. We ran the emcee fit using two repeats of 100 walkers with 1000 steps per walker to explore efficiently the full parameter space. We monitored the convergence of the fit using the Gelman-Rubin statistic (Gelman & Rubin 1992) which allowed to check that the two independent emcee analyses had produced consistent posterior PDFs for the jump parameters. Table 2 shows the results of the fit for a fixed depth across all filters. This includes the jump parameters as well as the derived parameters, for each we give the median of the posterior distribution and the 1 – σ interval. The posterior distributions of the jump parameters are given in the form of a corner plot (Foreman-Mackey 2016) in Fig. F.1.

We then performed a chromaticity check, using the same priors but allowing the depth to vary as a function of wavelength. All the depths found agree within 1-σ, they are shown in Fig. 9. We also performed an eccentric fit, adding epcosωp$\sqrt {{e_p}} \cos {\omega _p}$ and epsinωp$\sqrt {{e_p}} \sin {\omega _p}$ as jump parameters with ep the eccentricity and ωp the argument of periapsis, both of these quantities left as free parameters with no priors. We find a value of 0.12−0.09+0.18$0.12_{ - 0.09}^{ + 0.18}$ for the eccentricity. Following the Bayesian Information Criterion approximation to the Bayes Factor outlined in Wagenmakers (2007), we find a Bayes Factor of 10 654 which heavily favors the circular fit. Finally, we performed an analysis allowing the transit timings to vary to check for the existence of Transit Timing Variations (TTVs) which could be indicative of a third body in the system. We found the fitted timings to be consistent with no TTVs down to a few minutes, as shown on Fig. G.1. The lack of detected TTVs does not allow us to conclude on the presence of an additional planet orbiting TOI-4336 A.

thumbnail Fig. 8

Speckle autocorrelation functions (ACF) and 5σ sensitivity functions of the SOAR observations obtained for A (on the left) and B (on the right).

5.2 Search for additional candidates

Given the high occurrence rates of small planets around M stars (Dressing & Charbonneau 2013; Hsu et al. 2020), a system hosting one confirmed planet is likely to host additional ones. As the threshold for a possible TOI detection is set at S/N = 7 by the TESS pipelines, we then used the SHERLOCK12 package to determine if there could be other candidates in the system that would have been missed (Pozuelos et al. 2020; Demory et al. 2020). SHERLOCK is a community pipeline built on robust and deeply tested astrophysical tools that performs an iterative search for signals on a given star for missions such as TESS. To carry on the task, SHERLOCK follows the next steps: (1) Downloading and preparing the light curves from online databases using the lightkurve package. (2) Computing a Lomb-Scargle periodogram (Lomb 1976; Scargle 1982; VanderPlas 2018) to identify stellar variability, and the field of view plots using tpfplotter. Other preprocessing steps include correction of stellar variability, and adding a high RMS mask to remove outliers. (3) Following a multi-detrend approach, a bi-weight filter provided by the WōtanS13 package (Hippke et al. 2019) is applied on the light curves using a range of window sizes. (4) The transit search is performed iteratively on the nominal light curve as well as all the detrended light curves using the Transit Least Squares14 package (Hippke & Heller 2019). Once a signal is found above a certain threshold of S/N, it is masked for the next iterations, and this goes on until no more signals with sufficient S/N are found or the search reaches a certain number of iterations. (5) Vetting reports can be created for interesting signals in PDF format, where some metrics are computed and flagged in red when they are considered problematic. (6) A statistical validation can be performed using TRICERATOPS. (7) A Bayesian fit using the nested sampling algorithm of Allesfitter15 (Günther & Daylan 2019, 2021), Dynesty16 (Speagle 2020; Koposov et al. 2023), can be run for a set of selected signals to refine the system parameters. (8) An observation plan can be created based on the results of the fit for a chosen ground-based observatory.

For TOI-4336 A, we performed the transit search on the short-cadence data we extracted with custom apertures for Sectors 11, 38, and 64 of TESS, as described in Sect. 3.1. We first corrected the light curves for dilution using the value of 86.3% we obtained in the global analysis. We selected a range of ten window sizes between 0.2 and 1.3 days to generate the detrended light curves. We explored periods between 0.5 and 20 days, and selected a threshold signal of S/N ≥5, with a maximum number of five runs. We ran the transit search algorithm considering all three sectors combined in addition to individual sectors. We found the signal of TOI-4336 A b to be recovered easily in the first run of all our analyses with an S/N of 18.70. We also found a second potential candidate in the combined search with a period of 7.59 days, a duration of 1.74 h, and a depth of 1.18 ppt. We denote it as TIC 166184428.0217, assuming it corresponds to a transiting planet orbiting TOI-4336 A as well. The candidate is recovered in four out of the ten detrended light curves in the second transit search run with an S/N of 5.35. We performed the vetting stage and TIC 166184428.02 passed all the tests except for two: (1) The transit source offset computed from difference imaging, a method introduced in Bryson et al. (2017). (2) The per-pixel Box Least Squares (BLS, Kovács et al. 2002) search showed some deviation from the target. This method uses a fixed epoch and period for the BLS run on each pixels and computes the S/N of the detections. The S/N are then normalized and the centroid is found, then compared to the one found by difference imaging. These metrics are affected by the significant contamination from nearby sources which shifts the centroid on the TPFs. Given the presence of TOI-4336 Β in the aperture, we do not consider these failed tests to be critical. Finally, we performed the fitting stage of SHERLOCK to recover the system parameters. The results of the nested sampling fit are found in Table H.1, the folded light curve and the posterior distributions of the fitted parameters can be found in Figs. H.1 and H.2. Allesfitter uses the parameters found in the transit search as uniform priors for the Nested Sampling fit and a GP with a Matern 3/2 kernel for detrending.

TIC 166184428.02 is at the limit of the detection threshold we set in our transit search with SHERLOCK. The shallow depth of the signal makes it challenging to check whether the single events look consistent in shape with the transit of a planet. We also performed an independent transit search on the TESS data to assess the detection limits of the data, and we did not recover the candidate (see Sect. 5.3).

thumbnail Fig. 9

Chromaticity check of TOI-4336 A b. On the left: phase folded transits from TESS and ground-based observations. The filters above the dashed line correspond to the TESS and TRAPPIST-South observations which are diluted by TOI-4336 B. On the right: comparison of the depths obtained for the chromatic fit. All the bands agree within 1-σ, and the shaded region corresponds to the depth obtained from the achromatic fit for comparison.
thumbnail Fig. 10

Results of the injection-recovery tests on TESS light curves, described in Sect. 5.3. Almost all planets larger than 1.5 R⊕ were recovered successfully up to 30 d.

5.3 Detection limits and injection-recovery tests

We used the SPECULOOS-Southern Observatory to collect ground-based data to check if monitoring of TOI-4336 A from 1-m class telescopes would allow the detection of smaller planets thanks to an improved precision compared to TESS, as it has been successfully done in Delrez et al. (2022). We gathered twelve nights of photometric observations without any transits, in the Sloan-r′ filter, the details are shown in Table D.1. The photometric light curves for each night were extracted using the prose package, following the method described in Sect. 3.2 as for other ground-based data, and combined. We ran transit-search pipeline occultence18 to search for additional transiting planets around TOI-4336 A with both SPECULOOS and TESS (Sect. 3.1). occultence consists of several steps: cleaning (cosmic ray and spurious data point removal), an initial search to mask transit structures (for detrending) using Box Least Squares (BLS, Kovács et al. 2002), a detrending method, followed by a final transit search on the detrended light curve using BLS. For SPECULOOS detrending, we chose ridge regression (Hoerl & Kennard 1970) to fit polynomials of airmass, FWHM, sky background, δx and δy (change in the position of the target on the CCD). For each night, we fit for all combinations of these parameters with different polynomial orders (up to cubic) and selected the combination with the lowest AIC (Akaike information criterion, Akaike et al. 1973), which assesses a model’s likelihood of describing the data while penalizing larger numbers of parameters. For TESS’s detrending step instead we performed GP-detrending to capture the remaining correlated noise. We did not detect any transit structures in either SPECULOOS or TESS detrended light curve with S/N > 3.

To assess the detection efficiency of our transit-search pipeline we ran injection-recovery tests on both the SPECU-LOOS and TESS light curves (with real transits masked).

For each instrument, we generated transits for 3000 artificial planets using PYTRANSIT (Parviainen 2015) and injected each in turn into TOI-4336 A’s light curves. We used the Sloan–r’ and TESS limb-darkening coefficients from Table 2 for SPECULOOS and TESS, respectively. Each planet’s remaining parameters (radius Rp, period P, and inclination i) were drawn from the following uniform distributions: Rp ~ U(0.5,5.0) R⊕, cos i ~ U(cos imin, cos imax), log P ~ U(log 0.5, log 30.5) d, where U(a, b) represents a value drawn from a uniform distribution between a and b. imin and imax are the minimum and maximum inclinations for a transiting planet. Due to the shorter baseline of SPECULOOS data, we limited the explored period range up to 10.5 d. The host mass and radius were taken from Table 1. However, the inclination limits depend on the orbital period; therefore, when drawing the planetary parameters, we drew each inclination from a range set by the period. Only circu1ar orbits we considered. The time at which the first transit was injected was also drawn from a uniform distribution, φ ~ U(0,1), where φ is the phase of the period, such that the first transit was injected at φΡ from the start of observations.

We used the transit-search method described at the start of this section to search for the injected planets. A planet was recovered if at least one epoch from the highest likelihood period from BLS was within 1 h of an injected transit and the S/N of that transit is >3. This recovery criterion allows us to detect even single transits. The results from our injection-recovery for TOI-4336 A with TESS are shown in Fig. 10. Due to the day-night cycle and the complications in dealing with ground-based systematics, the injection-recovery results for TESS exceed SPECULOOS for all periods > 1 d, though the results are comparable for small (<1.5 R⊕) planets and very short periods. From the TESS results it is unlikely that there exist additional planets orbiting TOI-4336 A with R>1.5 R⊕ and P<30 d, though we cannot rule out smaller, or longer period, planets.

6 Discussion

6.1 Habitability

The stellar irradiation from its host, 1.50−0.17+0.18 S⊕$1.50_{ - 0.17}^{ + 0.18}\,{S_ \oplus }$, puts TOI-4336 A b very close at a distance consistent with the inner edge of the empirical HZ of the system (about 1.488 Sθ) (Kopparapu et al. 2013), as shown in Fig. 11. The limits of the HZ and the environments of planets orbiting close to these limits are not fully understood. Considering the error bars on the stellar irradiation, the planet could receive less or more irradiation than the irradiation at the empirical HZ limit. Probing the atmosphere of planets at the limit of the HZ will provide fur-therinsights into the environments at these stellarirradiations. A planetary radius of 2.12 R⊕ places the planet beyond the radius valley (e.g., Gupta et al. 2022), in the realm of mini-Neptunes where planets are likely to have an extended gaseous atmosphere. This expected atmosphere should allow for ease of atmospheric characterization. This makes TOI-4336 A b an interesting target to explore the environment at and around the empirical inner edge limit of the HZ and assess whether it is the likely case of a Mini-Neptune or the less likely case of a rocky planet at the inner edge of the HZ.

TOI-4336 A b orbits a host star in a triple system. While multiple host stars can influence the HZ boundaries (Kaltenegger & Haghighipour 2013; Haghighipour & Kaltenegger 2013; Kane & Hinkel 2013), the two other mid-M stars in the system orbit far enough apart to have no significant influence on the HZ limits. Using the bolometric luminosity from the SED fit (see Sect. 2.2) and for a median semi-major axis of 133 au (see Sect. 2.4), the added flux from the second star would only add 5.76710−7 S⊕ to the overall irradiation. Still, TOI-4336 A b provides a very interesting, as well as accessible, target to explore the region around the empirical HZ limits.

6.2 Formation in the triple system

Hundreds of planets have been discovered in binary systems (e.g., Raghavan et al. 2010; Matson et al. 2018). Yet binaries – especially on close or eccentric orbits – shrink the orbital real estate available for planets (Holman & Wiegert 1999). In addition, an inclined binary can trigger Kozai oscillations (Takeda et al. 2008), although this is impeded during the disk phase (Batygin et al. 2011) and in multiple-planet systems (Innanen et al. 1997; Kaib et al. 2011). In general, wider binaries have a progressively weaker influence. There is an observed deficit of exoplanets in systems with binaries closer than ~50 au (Kraus et al. 2016). Yet, in some cases, binaries wider than ~1000 au can reach very eccentric orbits due to external Galactic perturbations and destabilize the orbits of gas giants orbiting at Jupiter- to Saturn-like distances (Kaib et al. 2013). However, there is no evidence that the overall occurrence rate of exoplanets is strongly affected by wide binary companions (e.g., Kraus et al. 2016; Ziegler et al. 2020).

For the case of TOI-4336 A b, we do not expect the triple stellar configuration to have played a strong role in the planet’s formation or evolution. The maximum orbital radius that remains stable in the face of perturbations from a binary is a function of the companion’s mass and the binary orbital parameters (Holman & Wiegert 1999; Pilat-Lohinger et al. 2003). The closer stellar companion of TOI-4336 A has a semi-major axis of 133 au (see Sect. 2.4). Assuming a binary orbital eccentricity of 0.7 (the median for a thermal distribution) and equal stellar masses, the outermost stable orbit around TOI-4336 A would be about 3% of the binary semi-major axis, or ~4 au (assuming coplanarity between the stellar and planetary orbits; Holman & Wiegert 1999), more than an order of magnitude larger than the orbital radius of 0.09 au of TOI-4336 A b. If the binary’s periastron distance were 5 au, which is the minimum value in the 1-σ contours from the analysis in Sect. 2.4, then the binary’s eccentricity would have been so high that it could potentially have disrupted planet formation entirely. Given the existence of the planet, we expect that the binary’s eccentricity is likely no higher than 0.7–0.8, to avoid any drastic consequences during planetary growth. The second companion, with a semi-major axis of 2915 au, would have a much smaller effect. It therefore seems unlikely that the companion stars played a direct role in the formation or evolution of the sub-Neptune TOI-4336 A b.

6.3 Prospects for detailed characterization

We explore the possibility to constrain the planet mass using high-resolution spectroscopy. Following the mass–radius relationship by Chen & Kipping (2017), we estimate the planet mass to be 5.4−2.2+4.1 M⊕$5.4_{ - 2.2}^{ + 4.1}\,{M_ \oplus }$. Measuring the implied semi-amplitude of the spectrocopic orbit, 2.9 ± 1.6 m s−1, would require high-precision spectroscopy using stabilized spectrographs such as ESPRESSO (Pepe et al. 2010). Thanks to its apparent brightness (Vmag~12.9), the star is within reach of ESPRESSO with a 3.2 ± 1.8σ detection of the planet with 15 spectra. Because of its infrared brightness (Hmag=8.9), we expect detection to be also possible using NIRPS (Bouchy et al. 2017), which has been proven to show a long-term stability of 2–3 m s−1 for bright targets. There is already an ongoing radial velocity campaign with both these facilities to measure the mass of TOI-4336 A b.

Given the infrared brightness of the host (Kmag ~ 8.6) and the favorable planet-to-star radius ratio due to the low-mass host star, we assess the suitability of TOI-4336 A b for the characterization of an atmosphere. We first calculated the transmission spectroscopy metric (TSM, Kempton et al. 2018) as it is a convenient metric to compare the amenability for characterization of different planets. We find that TOI-4336 A b has an exceptionally high TSM of 83−4+5$83_{ - 4}^{ + 5}$ in the temperate sub-Neptune category, making it even more favourable than some of the best-studied sub-Neptunes for a detailed atmospheric characterization by transit transmission spectroscopy with HST and JWST (e.g., LHS 1140 b, Edwards et al. 2021).

We then used ExoTransmit (Kempton et al. 2017) to calculate two simulated transmission spectra for TOI-4336 A b, assuming a typical H/He-rich atmosphere with an isothermal P-T profile at the planet’s equilibrium temperature (~300K). We also assume equilibrium chemistry and include all the opacities available in ExoTransmit. The spectra are shown in the bottom right panel of Fig. 11. The first spectrum assumes a mass equal to the estimate derived from empirical mass-radius relations presented in Chen & Kipping (2017), while the second allows for a case where the planet has a similar density to one of the densest known mini-Neptune (Kepler-231 b, Hadden & Lithwick 2014). A lower mass would simply increase the scale height making atmospheric investigations easier. We find that for the mass estimate of ~5.4 M⊕, we could detect atmospheric features with > 16σ significance in the case of a transparent atmosphere with just three transits of HST. Even in the higher mass limit of ~10.1 M⊕, our simulations show that we could still detect a putative atmosphere to > 5σ with the same observing setup. This makes TOI-4336 A b one of the most promising sub-Neptunes for atmospheric investigations. However, the presence of clouds and/or hazes could hinder the detection of atmospheric signals. Atmospheric exploration with HST is then crucial to determine whether this is the case or not. In fact, a proposal for these observations has already been accepted for HST and the campaign to measure the atmosphere of this exquisite sub-Neptune is underway. We expect this planet to become one of the best studied in its category in the coming years.

7 Conclusions

Our analyses show that TOI-4336 A b is a temperate subNeptune orbiting a nearby star part of a triple M-dwarf system. We tested the chromaticity of the transit with multicolor photometry and found the transit depths to be consistent within 1-σ in all the bands. We found no significant variations in the transit depth between odd and even transits, which could be an indicator of binarity or a blended source. In addition, we statistically validate the planetary nature using archival imaging spanning 47 yr of observations, and high angular resolution imaging, which excludes unresolved companions with spectral type between M4 and early-L between 0.2 and 1.2″ in separation. Our global model yields a planetary radius of 2.12−0.09+0.08 R⊕$2.12_{ - 0.09}^{ + 0.08}\,{R_ \oplus }$ and an orbital period of 16.33 days, resulting in an incident irradiation of 1.50−0.17+0.18 S⊕$1.50_{ - 0.17}^{ + 0.18}\,{S_ \oplus }$. This incident irradiation would lead to an equilibrium temperature estimate of 308 ± 9 K, assuming a Bond albedo of zero and perfect heat redistribution (f = 1/4). We find the host to be an M3.5 star (Teff=3298−73+75 K${T_{{\rm{eff}}}} = 3298_{ - 73}^{ + 75}\,{\rm{K}}$, M* = 0.33 ± 0.02 M⊙), with estimated semi-major axes of 133 au and 2915 au with respect to the other two mid-M stars of the triple system.

The radius of this new planet puts it most likely in the miniNeptune category, and is thus likely to have retained an extended atmosphere. The incident radiation places TOI-4336 A b at the inner edge of the empirical HZ, which makes it a good candidate to explore this region and its consequences on habitability. We investigated the implications of the triple star system on the formation of the planet and found that the eccentricity of the closer pair should be no higher than 0.8. In addition, the orbital configuration of the system implies TOI-4336 B did not have any effect on the formation of TOI-4336 A b. Indeed, the planet would need to be at more than 40 times the orbital distance from its host star for the companion to disrupt the formation of the planet.

TOI-4336 A b shows similar properties to the widely studied K2-18 b (e.g., Benneke et al. 2017, 2019; Cloutier et al. 2017; Tsiaras et al. 2019). We examined the suitability of this temperate sub-Neptune for detailed characterization and found it is a prime target for atmospheric studies with HST and JWST. The determination of its mass is a key ingredient for the interpretation of transmission spectra. We find both ESPRESSO and NIRPS are appropriate to measure the radial velocities of the planet and constrain its mass.

thumbnail Fig. 11

Top left panel: complete sample of all known exoplanets with measured masses (Data from NASA Exoplanet Archive, Apr 2024, https://exoplanetarchive.ipac.caltech.edu/). Following Kempton et al. (2018), terrestrial planets correspond to planetary radii below 1.25 R⊕, small sub-Neptunes radii between 1.25 and 2.75 R⊕, large sub-Neptunes radii between 2.75 and 4 R⊕, and sub-Jovians are planets with radii between 4 and 10 R⊕. TOI-4336 A b is shown by the star symbol. Top right panel: zoomed in view that puts TOI-4336 A b in the context of small sub-Neptunes only. Bottom left panel: stellar effective temperature as a function of insolation of transiting exoplanets orbiting hosts cooler than 4000 K. The solid black line denotes the empirical HZ boundary and the dashed line the conservative one (Kopparapu et al. 2013). The size of points scales with the planetary radius. The points are colored according to their equilibrium temperature. TOI-4336 A b is highlighted with errorbars, and it is well placed at the inner edge of the HZ of its host star. Bottom right panel: synthetic transmission spectra for TOI-4336 A b, assuming a typical H/He-rich atmosphere with isothermal temperature profile at the planet’s equilibrium temperature (~300 K). The simulated data points are for the Wide Field Camera 3 instrument of HST.

Acknowledgements

Funding for the TESS mission is provided by NASA’s Science Mission Directorate. We acknowledge the use of public TESS data from pipelines at the TESS Science Office and at the TESS Science Processing Operations Center. This research has made use of the Exoplanet Follow-up Observation Program website, which is operated by the California Institute of Technology, under contract with the National Aeronautics and Space Administration under the Exoplanet Exploration Program. This paper includes data collected by the TESS mission that are publicly available from the Mikulski Archive for Space Telescopes (MAST). Based on data collected by the SPECULOOS-South Observatory at the ESO Paranal Observatory in Chile.The ULiege s contribution to SPECULOOS has received funding from the European Research Council under the European Union’s Seventh Framework Programme (FP/2007-2013) (grant Agreement n° 336480/SPECULOOS), from the Balzan Prize and Francqui Foundations, from the Belgian Scientific Research Foundation (F.R.S.-FNRS; grant n° T.0109.20), from the University of Liege, and from the ARC grant for Concerted Research Actions financed by the Wallonia-Brussels Federation. This work is supported by a grant from the Simons Foundation (PI Queloz, grant number 327127). This research is in part funded by the European Union s Horizon 2020 research and innovation programme (grants agreements n° 803193/BEBOP), and from the Science and Technology Facilities Council (STFC; grant n° ST/S00193X/1, and ST/W000385/1). The material is based upon work supported by NASA under award number 80GSFC21M0002. Based on data collected by the TRAPPIST-South telescope at the ESO La Silla Observatory. TRAPPIST is funded by the Belgian Fund for Scientific Research (Fond National de la Recherche Scientifique, FNRS) under the grant PDR T.0120.21, with the participation of the Swiss National Science Fundation (SNF). Based on data collected under the ExTrA project at the ESO La Silla Paranal Observatory. ExTrA is a project of Institut de Planétologie et d’Astrophysique de Grenoble (IPAG/CNRS/UGA), funded by the European Research Council under the ERC Grant Agreement n. 337591-ExTrA. This work has been supported by a grant from Labex OSUG@2020 (Investissements d’avenir — ANR10 LABX56). This work has been carried out within the framework of the NCCR PlanetS supported by the Swiss National Science Foundation. The Digitized Sky Surveys were produced at the Space Telescope Science Institute under U.S. Government grant NAG W-2166. The images of these surveys are based on photographic data obtained using the Oschin Schmidt Telescope on Palomar Mountain and the UK Schmidt Telescope. The plates were processed into the present compressed digital form with the permission of these institutions. This work makes use of observations from the LCOGT network. Part of the LCOGT telescope time was granted by NOIRLab through the Mid-Scale Innovations Program (MSIP). MSIP is funded by NSF. Based in part on observations obtained at the Southern Astro-physical Research (SOAR) telescope, which is a joint project of the Ministério da Ciência, Tecnologia e Inovaçôes (MCTI/LNA) do Brasil, the US National Science Foundation NOIRLab, the University of North Carolina at Chapel Hill (UNC), and Michigan State University (MSU). IRAF was distributed by the National Optical Astronomy Observatory, which was managed by the Association of Universities for Research in Astronomy (AURA) under a cooperative agreement with the National Science Foundation. 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, Nic Scott, Elliott P. Horch, and Emmett Quigley. Zorro was mounted on the Gemini 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 Investigacion y Desarrollo (Chile), Ministerio de Ciencia, Tecnologia e Innovacion (Argentina), Ministério da Ciência, Tecnologia, Inovaçôes e Comunicaçôes (Brazil), and Korea Astronomy and Space Science Institute (Republic of Korea). The postdoctoral fellowship of K.B. is funded by F.R.S.-FNRS grant T.0109.20 and by the Francqui Foundation. This publication benefits from the support of the French Community of Belgium in the context of the FRIA Doctoral Grant awarded to M.T. E.D. acknowledges support from the innovation and research Horizon 2020 program in the context of the Marie Sklodowska-Curie subvention 945298. B.V.R. thanks the HeisingSimons Foundation for support. M.G. is F.R.S.-FBRS Research Director. F.J.P. acknowledges financial support from the grant CEX2021-001131-S funded by MCIN/AEI/10.13039/501100011033 and through projects PID2019-109522GB-C52 and PID2022-137241NB-C43. KAC and SNQ acknowledge support from the TESS mission via subaward s3449 from MIT. B.-O.D. acknowledges support from the Swiss State Secretariat for Education, Research and Innovation (SERI) under contract number MB22.00046. E.J. is a Belgian FNRS Senior Research Associate. YGMC acknowledges support from UNAM-PAPIIT-IG101321. L.D. is an F.R.S.-FNRS Postdoctoral Researcher. SNR acknowledges support from the French Programme National de Planétologie (PNP). This research made use of Lightkurve, a Python package for Kepler and TESS data analysis (Lightkurve Collaboration, 2018). This research has made use of the NASA Exoplanet Archive, which is operated by the California Institute of Technology, under contract with the National Aeronautics and Space Administration under the Exo-planet Exploration Program. 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. This work made use of tpfplotter by J. Lillo-Box (publicly available in www.github.com/jlillo/tpfplotter), which also made use of the python packages astropy, lightkurve, matplotlib and numpy. This work has made use of data from the European Space Agency (ESA) mission Gaia (https://www.cosmos.esa.int/gaia), processed by the Gaia Data Processing and Analysis Consortium (DPAC, https://www.cosmos.esa.int/web/gaia/dpac/consortium). Funding for the DPAC has been provided by national institutions, in particular the institutions participating in the Gaia Multilateral Agreement. This research has made use of the Washington Double Star Catalog maintained at the U.S. Naval Observatory. This publication makes use of data products from the Two Micron All Sky Survey, which is a joint project of the University of Massachusetts and the Infrared Processing and Analysis Center/California Institute of Technology, funded by the National Aeronautics and Space Administration and the National Science Foundation.

Appendix A Color-magnitude diagram

thumbnail Fig. A.1

Color-magnitude diagram for nearby M dwarfs.The TOI-4336 system is marked with a star symbol and color-coded to distinguish the three stars. TOI-4336 C appear redder than A and B, which is consistent with a later spectral type (M4±0.5 compared to M3.5±0.5. We used the extended target list of the SPECULOOS survey which gathers over 14 000 M-dwarf stars within 40 parsecs to generate the diagram (Sebastian et al. 2021).

Appendix B Age estimation

thumbnail Fig. B.1

Age distribution of the GALAH survey (light grey) and the distribution of stars within our selection criteria (dark grey). The yellow line indicates the median age and the yellow bands are the 16 percentile and 84 percentile regions. We estimate an age of 6.7−3.1+2.7$6.7_{ - 3.1}^{ + 2.7}$ Gyr for the system.

Appendix C Triple system orbital analysis

thumbnail Fig. C.1

Top row: AB (left) and AC (right) orbits with ten random orbital solutions from posterior of the orbit fit. Middle row: separation with WDS astrometry for AB (left) and AC (right) orbits with a hundred random orbital solutions from posterior of the fit.Bottom row: the same for the position angle.
Table C.1

Astrometry measurements from WDS (Mason et al. 2001) for AC (left) and AB (right) systems used in the orbital fit.

Table C.2

Summary of orbital parameters for the TOI-4336 system.

thumbnail Fig. C.2

Posterior samples of the AB system. Histogram sub-panels show the posterior distribution, with the median and 68% confidence intervals marked by dashed lines, with titles quantifying those ranges.

Appendix D Ground-based photometric observations

Table D.1

Summary of the ground-based follow-up observations obtained for the validation of TOI-4336 A b.

Appendix E Global model errors and scaling factors

Table E.1

Global model errors and scaling factors for the MCMC analysis of the photometric data.

For each transit observation, the chosen baseline model is denoted as p(αn) with α a parameter such as t = time, a = airmass, f = FWHM of the PSF, x or y = position of the target on the detector along the x− and y− axes, and n the order of the considered polynomial. The residual RMS is also given with the corresponding exposure time and the error scaling factors ßw and ßr. The column T0 represents the expected mid-transit times.

Appendix F Posterior distributions of the transit model

thumbnail Fig. F.1

Posterior distributions for the main jump parameters used in the emcee fit (see Sect. 5): the effective temperature (Teff), the metallicity ([Fe/H]), the log of the stellar density (log ρ*), the log of the stellar mass (log M*), the dilution in the TESS and z′ bands, the epoch (T0), the transit depth (dF), the log of the orbital period (log P) and the cosine of the orbital inclination of the planet (cos i).

Appendix G Transit timing variations

thumbnail Fig. G.1

O-C diagram obtained from a global fit allowing transit timing variations. We used as reference the timing obtained in the global transit model of Table 2.

Appendix H Search for additional candidates in the system

thumbnail Fig. H.1

Transit model of TIC 166184428.02. Top panel: phase folded photometry of TIC 166184428.02. The best-fit model is shown in solid black. The grey points are the raw flux and the red points are the 15-minutes binned flux. Bottom panel: residuals of the best-fit model.
Table H.1

Properties of TIC 166184428.02 obtained from the fitting stage of SHERLOCK, see Sect. 5.2.

thumbnail Fig. H.2

Posterior distributions for the main parameters fitted in the Nested Sampling fit for the second candidate found in the TOI-4336 system, TIC 166184428.02 (see Sect. 5.2), with Rb the candidate radius, R* the stellar radius, ib the inclination, T0;b the epoch, Pb the period, q1:1c and q2:1c the Kipping parametrisation (Kipping 2013) for the limb darkening coefficients, logσlc the error scaling parameter, and lnσ(lc) and lnp(lc) the two GP parameters related to the amplitude and length scale.

References

  1. Adams, E. R., Seager, S., & Elkins-Tanton, L. 2008, ApJ, 673, 1160 [Google Scholar]
  2. Akaike, H., Petrov, B. N., & Csaki, F. 1973, Second International Symposium on Information Theory (Budapest: Akadémiai Kiadó) [Google Scholar]
  3. 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]
  4. Bailer-Jones, C. A. L., Rybizki, J., Fouesneau, M., Demleitner, M., & Andrae, R. 2021, AJ, 161, 147 [Google Scholar]
  5. Barkaoui, K., Timmermans, M., Soubkiou, A., et al. 2023, A&A, 677, A38 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  6. Batygin, K., Morbidelli, A., & Tsiganis, K. 2011, A&A, 533, A7 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  7. Bean, J. L., Raymond, S. N., & Owen, J. E. 2021, J. Geophys. Res. (Planets), 126, e06639 [NASA ADS] [Google Scholar]
  8. Benneke, B., Werner, M., Petigura, E., et al. 2017, ApJ, 834, 187 [Google Scholar]
  9. Benneke, B., Wong, I., Piaulet, C., et al. 2019, ApJ, 887, L14 [Google Scholar]
  10. Berger, T. A., Huber, D., Gaidos, E., van Saders, J. L., & Weiss, L. M. 2020, AJ, 160, 108 [NASA ADS] [CrossRef] [Google Scholar]
  11. Bianchi, L., Shiao, B., & Thilker, D. 2017, ApJS, 230, 24 [Google Scholar]
  12. Bonfils, X., Almenara, J. M., Jocou, L., et al. 2015, SPIE Conf. Ser., 9605, 96051L [NASA ADS] [Google Scholar]
  13. Borucki, W. J., Koch, D. G., Basri, G., et al. 2011, ApJ, 736, 19 [Google Scholar]
  14. Bouchy, F., Doyon, R., Artigau, É., et al. 2017, The Messenger, 169, 21 [NASA ADS] [Google Scholar]
  15. Bradley, L., Sipőcz, B., Robitaille, T., et al. 2023, https://doi.org/10.5281/zenodo.7946442 [Google Scholar]
  16. Brown, T. M., Baliber, N., Bianco, F. B., et al. 2013, PASP, 125, 1031 [Google Scholar]
  17. Bryson, S. T., Jenkins, J. M., Klaus, T. C., et al. 2017, Kepler Data Processing Handbook: Target and Aperture Definitions: Selecting Pixels for Kepler Downlink, Kepler Science Document KSCI-19081-002, 3, ed. J. M. Jenkins [Google Scholar]
  18. Buder, S., Asplund, M., Duong, L., et al. 2018, MNRAS, 478, 4513 [Google Scholar]
  19. Burdanov, A., Delrez, L., Gillon, M., & Jehin, E. 2018, in Handbook of Exoplanets, eds. H. J. Deeg, & J. A. Belmonte (Springer International Publishing AG), 130 [Google Scholar]
  20. Caldwell, D. A., Tenenbaum, P., Twicken, J. D., et al. 2020, RNAAS, 4, 201 [NASA ADS] [Google Scholar]
  21. Chen, J., & Kipping, D. 2017, ApJ, 834, 17 [Google Scholar]
  22. Ciardi, D. R., Beichman, C. A., Horch, E. P., & Howell, S. B. 2015, ApJ, 805, 16 [NASA ADS] [CrossRef] [Google Scholar]
  23. Claret, A. 2018, A&A, 618, A20 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  24. Claret, A., Hauschildt, P. H., & Witte, S. 2012, A&A, 546, A14 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  25. Cloutier, R., & Menou, K. 2020, AJ, 159, 211 [NASA ADS] [CrossRef] [Google Scholar]
  26. Cloutier, R., Astudillo-Defru, N., Doyon, R., et al. 2017, A&A, 608, A35 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  27. Cointepas, M., Almenara, J. M., Bonfils, X., et al. 2021, A&A, 650, A145 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  28. Cushing, M. C., Vacca, W. D., & Rayner, J. T. 2004, PASP, 116, 362 [NASA ADS] [CrossRef] [Google Scholar]
  29. Cushing, M. C., Rayner, J. T., & Vacca, W. D. 2005, ApJ, 623, 1115 [Google Scholar]
  30. Cutri, R. M., Skrutskie, M. F., van Dyk, S., et al. 2003, VizieR Online Data Catalog II/246 [Google Scholar]
  31. Cutri, R. M., Wright, E. L., Conrow, T., et al. 2021, VizieR Online Data Catalog II/328 [Google Scholar]
  32. Delrez, L., Gillon, M., Queloz, D., et al. 2018, SPIE Conf. Ser., 10700, 107001I [Google Scholar]
  33. Delrez, L., Murray, C. A., Pozuelos, F. J., et al. 2022, A&A, 667, A59 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  34. Demory, B. O., Pozuelos, F. J., Gomez Maqueo Chew, Y., et al. 2020, A&A, 642, A49 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  35. Dransfield, G., Timmermans, M., Triaud, A. H. M. J., et al. 2024, MNRAS, 527, 35 [Google Scholar]
  36. Dressing, C. D., & Charbonneau, D. 2013, ApJ, 767, 95 [Google Scholar]
  37. Edwards, B., Changeat, Q., Mori, M., et al. 2021, AJ, 161, 44 [Google Scholar]
  38. Foreman-Mackey, D. 2016, J. Open Source Softw., 1, 24 [Google Scholar]
  39. Foreman-Mackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, PASP, 125, 306 [Google Scholar]
  40. Fulton, B. J., Petigura, E. A., Howard, A. W., et al. 2017, AJ, 154, 109 [Google Scholar]
  41. Furlan, E., & Howell, S. B. 2017, AJ, 154, 66 [NASA ADS] [CrossRef] [Google Scholar]
  42. Gagné, J., Mamajek, E. E., Malo, L., et al. 2018, ApJ, 856, 23 [Google Scholar]
  43. Gaia Collaboration 2022, VizieR Online Data Catalog: I/355 [Google Scholar]
  44. Garcia, L. J., Timmermans, M., Pozuelos, F. J., et al. 2021, Astrophysics Source Code Library, [record ascl:2111.006] [Google Scholar]
  45. Garcia, L. J., Timmermans, M., Pozuelos, F. J., et al. 2022, MNRAS, 509, 4817 [Google Scholar]
  46. Gelman, A., & Rubin, D. B. 1992, Statist. Sci., 7, 457 [NASA ADS] [Google Scholar]
  47. Giacalone, S., & Dressing, C. D. 2020, Astrophysics Source Code Library, [record ascl:2002.004] [Google Scholar]
  48. Giacalone, S., Dressing, C. D., Jensen, E. L. N., et al. 2021, AJ, 161, 24 [Google Scholar]
  49. Gillon, M., Jehin, E., Magain, P., et al. 2011, EPJ Web Conf., 11, 06002 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  50. Gillon, M., Triaud, A. H. M. J., Fortney, J. J., et al. 2012, A&A, 542, A4 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  51. Goodman, J., & Weare, J. 2010, Commun. Appl. Math. Computat. Sci., 5, 65 [Google Scholar]
  52. Guerrero, N. M., Seager, S., Huang, C. X., et al. 2021, ApJS, 254, 39 [NASA ADS] [CrossRef] [Google Scholar]
  53. Günther, M. N., & Daylan, T. 2019, Astrophysics Source Code Library [record ascl:1903.003] [Google Scholar]
  54. Günther, M. N., & Daylan, T. 2021, ApJS, 254, 13 [CrossRef] [Google Scholar]
  55. Gupta, A., & Schlichting, H. E. 2019, MNRAS, 487, 24 [Google Scholar]
  56. Gupta, A., Nicholson, L., & Schlichting, H. E. 2022, MNRAS, 516, 4585 [NASA ADS] [CrossRef] [Google Scholar]
  57. Hadden, S., & Lithwick, Y. 2014, ApJ, 787, 80 [Google Scholar]
  58. Haghighipour, N., & Kaltenegger, L. 2013, ApJ, 777, 166 [Google Scholar]
  59. Henden, A. A., Levine, S., Terrell, D., & Welch, D. L. 2015, in Am. Astron. Soc. Meeting Abstracts, 225, 336.16 [NASA ADS] [Google Scholar]
  60. Hippke, M., & Heller, R. 2019, Astrophysics Source Code Library [record ascl:1910.007] [Google Scholar]
  61. Hippke, M., David, T. J., Mulders, G. D., & Heller, R. 2019, AJ, 158, 143 [Google Scholar]
  62. Hoerl, A. E., & Kennard, R. W. 1970, Technometrics, 12, 55 [Google Scholar]
  63. Holman, M. J., & Wiegert, P. A. 1999, AJ, 117, 621 [Google Scholar]
  64. Howell, S. B., & Furlan, E. 2022, Front. Astron. Space Sci., 9, 871163 [NASA ADS] [CrossRef] [Google Scholar]
  65. Howell, S. B., Everett, M. E., Sherry, W., Horch, E., & Ciardi, D. R. 2011, AJ, 142, 19 [Google Scholar]
  66. Howell, S. B., Everett, M. E., Horch, E. P., et al. 2016, ApJ, 829, L2 [NASA ADS] [CrossRef] [Google Scholar]
  67. Hsu, D. C., Ford, E. B., & Terrien, R. 2020, MNRAS, 498, 2249 [Google Scholar]
  68. Huang, C. X., Vanderburg, A., Pál, A., et al. 2020, RNAAS, 4, 204 [Google Scholar]
  69. Husser, T. O., Wende-von Berg, S., Dreizler, S., et al. 2013, A&A, 553, A6 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  70. Innanen, K. A., Zheng, J. Q., Mikkola, S., & Valtonen, M. J. 1997, AJ, 113, 1915 [NASA ADS] [CrossRef] [Google Scholar]
  71. Jehin, E., Gillon, M., Queloz, D., et al. 2011, The Messenger, 145, 2 [NASA ADS] [Google Scholar]
  72. Jenkins, J. M., Caldwell, D. A., & Borucki, W. J. 2002, ApJ, 564, 495 [NASA ADS] [CrossRef] [Google Scholar]
  73. Jenkins, J. M., Caldwell, D. A., Chandrasekaran, H., et al. 2010, ApJ, 713, L87 [Google Scholar]
  74. Jenkins, J. M., Twicken, J. D., McCauliff, S., et al. 2016, SPIE Conf. Ser., 9913, 99133E [Google Scholar]
  75. Jenkins, J. M., Tenenbaum, P., Seader, S., et al. 2020, Kepler Data Processing Handbook: Transiting Planet Search, Kepler Science Document KSCI-19081-003, 9, ed. J. M. Jenkins [Google Scholar]
  76. Kaib, N. A., Raymond, S. N., & Duncan, M. J. 2011, ApJ, 742, L24 [NASA ADS] [CrossRef] [Google Scholar]
  77. Kaib, N. A., Raymond, S. N., & Duncan, M. 2013, Nature, 493, 381 [Google Scholar]
  78. Kaltenegger, L. 2017, ARA&A, 55, 433 [NASA ADS] [CrossRef] [Google Scholar]
  79. Kaltenegger, L., & Haghighipour, N. 2013, ApJ, 777, 165 [NASA ADS] [CrossRef] [Google Scholar]
  80. Kane, S. R., & Hinkel, N. R. 2013, ApJ, 762, 7 [Google Scholar]
  81. Kasting, J. F., Whitmire, D. P., & Reynolds, R. T. 1993, Icarus, 101, 108 [Google Scholar]
  82. Kempton, E. M. R., Lupu, R., Owusu-Asare, A., Slough, P., & Cale, B. 2017, PASP, 129, 044402 [NASA ADS] [CrossRef] [Google Scholar]
  83. Kempton, E. M. R., Bean, J. L., Louie, D. R., et al. 2018, PASP, 130, 114401 [Google Scholar]
  84. Kesseli, A. Y., West, A. A., Veyette, M., et al. 2017, ApJS, 230, 16 [NASA ADS] [CrossRef] [Google Scholar]
  85. Kiman, R., Faherty, J. K., Cruz, K. L., et al. 2021, AJ, 161, 277 [NASA ADS] [CrossRef] [Google Scholar]
  86. Kinemuchi, K., Barclay, T., Fanelli, M., et al. 2012, PASP, 124, 963 [Google Scholar]
  87. Kipping, D. M. 2013, MNRAS, 435, 2152 [Google Scholar]
  88. Kirkpatrick, J. D., Looper, D. L., Burgasser, A. J., et al. 2010, ApJS, 190, 100 [Google Scholar]
  89. Koposov, S., Speagle, J., Barbary, K., et al. 2023, https://doi.org/10.5281/zenodo.8408702 [Google Scholar]
  90. Kopparapu, R. K., Ramirez, R., Kasting, J. F., et al. 2013, ApJ, 765, 131 [NASA ADS] [CrossRef] [Google Scholar]
  91. Kopparapu, R. K., Wolf, E. T., Arney, G., et al. 2017, ApJ, 845, 5 [NASA ADS] [CrossRef] [Google Scholar]
  92. Kovács, G., Zucker, S., & Mazeh, T. 2002, A&A, 391, 369 [Google Scholar]
  93. Kraus, A. L., Ireland, M. J., Huber, D., Mann, A. W., & Dupuy, T. J. 2016, AJ, 152, 8 [NASA ADS] [CrossRef] [Google Scholar]
  94. Lasker, B. M., Doggett, J., McLean, B., et al. 1996, in ASP Conf. Ser., 101, Astronomical Data Analysis Software and Systems V, eds. G. H. Jacoby, & J. Barnes, 88 [NASA ADS] [Google Scholar]
  95. Lee, E. J., Chiang, E., & Ormel, C. W. 2014, ApJ, 797, 95 [Google Scholar]
  96. Lépine, S., Rich, R. M., & Shara, M. M. 2003, AJ, 125, 1598 [CrossRef] [Google Scholar]
  97. Lépine, S., Rich, R. M., & Shara, M. M. 2007, ApJ, 669, 1235 [CrossRef] [Google Scholar]
  98. Lépine, S., Hilton, E. J., Mann, A. W., et al. 2013, AJ, 145, 102 [Google Scholar]
  99. Lester, K. V., Matson, R. A., Howell, S. B., et al. 2021, AJ, 162, 75 [NASA ADS] [CrossRef] [Google Scholar]
  100. Li, J., Tenenbaum, P., Twicken, J. D., et al. 2019, PASP, 131, 024506 [NASA ADS] [CrossRef] [Google Scholar]
  101. Lightkurve Collaboration (Cardoso, J. V. d. M., et al.) 2018, Astrophysics Source Code Library [record ascl:1812.013] [Google Scholar]
  102. Lomb, N. R. 1976, Ap&SS, 39, 447 [Google Scholar]
  103. Luque, R., & Pallé, E. 2022, Science, 377, 6611 [Google Scholar]
  104. Mandel, K., & Agol, E. 2002, ApJ, 580, L171 [Google Scholar]
  105. Mann, A. W., Brewer, J. M., Gaidos, E., Lépine, S., & Hilton, E. J. 2013, AJ, 145, 52 [Google Scholar]
  106. Mann, A. W., Dupuy, T., Kraus, A. L., et al. 2019, ApJ, 871, 63 [Google Scholar]
  107. Mason, B. D., Wycoff, G. L., Hartkopf, W. I., Douglass, G. G., & Worley, C. E. 2001, AJ, 122, 3466 [Google Scholar]
  108. Matson, R. A., Howell, S. B., Horch, E. P., & Everett, M. E. 2018, AJ, 156, 31 [NASA ADS] [CrossRef] [Google Scholar]
  109. McCully, C., Volgenau, N. H., Harbeck, D.-R., et al. 2018, SPIE Conf. Ser., 10707, 107070K [Google Scholar]
  110. Medina, A. A., Winters, J. G., Irwin, J. M., & Charbonneau, D. 2022, ApJ, 935, 104 [NASA ADS] [CrossRef] [Google Scholar]
  111. Newton, E. R., Irwin, J., Charbonneau, D., et al. 2017, ApJ, 834, 85 [Google Scholar]
  112. Owen, J. E., & Wu, Y. 2013, ApJ, 775, 105 [Google Scholar]
  113. Paredes, L. A., Henry, T. J., Quinn, S. N., et al. 2021, AJ, 162, 176 [NASA ADS] [CrossRef] [Google Scholar]
  114. Parviainen, H. 2015, MNRAS, 450, 3233 [Google Scholar]
  115. Parviainen, H., & Aigrain, S. 2015, MNRAS, 453, 3821 [Google Scholar]
  116. Pass, E. K., Charbonneau, D., Irwin, J. M., & Winters, J. G. 2022, ApJ, 936, 109 [NASA ADS] [CrossRef] [Google Scholar]
  117. Pass, E. K., Winters, J. G., Charbonneau, D., et al. 2023, AJ, 166, 11 [NASA ADS] [CrossRef] [Google Scholar]
  118. Pearce, L. A., Kraus, A. L., Dupuy, T. J., et al. 2020, ApJ, 894, 115 [CrossRef] [Google Scholar]
  119. Pepe, F. A., Cristiani, S., Rebolo Lopez, R., et al. 2010, SPIE Conf. Ser., 7735, 77350F [Google Scholar]
  120. Pickles, A. J. 1998, PASP, 110, 863 [NASA ADS] [CrossRef] [Google Scholar]
  121. Pilat-Lohinger, E., Funk, B., & Dvorak, R. 2003, A&A, 400, 1085 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  122. Pont, F., Zucker, S., & Queloz, D. 2006, MNRAS, 373, 231 [NASA ADS] [CrossRef] [Google Scholar]
  123. Pozuelos, F. J., Suárez, J. C., de Elía, G. C., et al. 2020, A&A, 641, A23 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  124. Pozuelos, F. J., Timmermans, M., Rackham, B. V., et al. 2023, A&A, 672, A70 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  125. Raghavan, D., McAlister, H. A., Henry, T. J., et al. 2010, ApJS, 190, 1 [Google Scholar]
  126. Rayner, J. T., Toomey, D. W., Onaka, P. M., et al. 2003, PASP, 115, 362 [NASA ADS] [CrossRef] [Google Scholar]
  127. Rayner, J. T., Cushing, M. C., & Vacca, W. D. 2009, ApJS, 185, 289 [Google Scholar]
  128. Reid, I. N., Brewer, C., Brucato, R. J., et al. 1991, PASP, 103, 661 [NASA ADS] [CrossRef] [Google Scholar]
  129. Reid, I. N., Hawley, S. L., & Gizis, J. E. 1995, AJ, 110, 1838 [Google Scholar]
  130. Ricker, G. R., Winn, J. N., Vanderspek, R., et al. 2015, J. Astron. Telescopes Instrum. Syst., 1, 014003 [Google Scholar]
  131. Riddick, F. C., Roche, P. F., & Lucas, P. W. 2007, MNRAS, 381, 1067 [Google Scholar]
  132. Rogers, L. A. 2015, ApJ, 801, 41 [Google Scholar]
  133. Scargle, J. D. 1982, ApJ, 263, 835 [Google Scholar]
  134. Schwarz, G. 1978, Ann. Statist., 6, 461 [Google Scholar]
  135. Scott, N. J., Howell, S. B., Gnilka, C. L., et al. 2021, Front. Astron. Space Sci., 8, 138 [NASA ADS] [CrossRef] [Google Scholar]
  136. Sebastian, D., Gillon, M., Ducrot, E., et al. 2021, A&A, 645, A100 [EDP Sciences] [Google Scholar]
  137. Sharma, S., Stello, D., Buder, S., et al. 2018, MNRAS, 473, 2004 [NASA ADS] [CrossRef] [Google Scholar]
  138. Speagle, J. S. 2020, MNRAS, 493, 3132 [Google Scholar]
  139. Stassun, K. G., & Torres, G. 2016, AJ, 152, 180 [Google Scholar]
  140. Stassun, K. G., & Torres, G. 2021, ApJ, 907, L33 [NASA ADS] [CrossRef] [Google Scholar]
  141. Stassun, K. G., Collins, K. A., & Gaudi, B. S. 2017, AJ, 153, 136 [Google Scholar]
  142. Stassun, K. G., Corsaro, E., Pepper, J. A., & Gaudi, B. S. 2018, AJ, 155, 22 [Google Scholar]
  143. Stassun, K. G., Oelkers, R. J., Paegert, M., et al. 2019, AJ, 158, 138 [Google Scholar]
  144. Stevenson, K. B., Bean, J. L., Seifahrt, A., et al. 2016, ApJ, 817, 141 [NASA ADS] [CrossRef] [Google Scholar]
  145. Takeda, G., Kita, R., & Rasio, F. A. 2008, ApJ, 683, 1063 [NASA ADS] [CrossRef] [Google Scholar]
  146. Tokovinin, A. 2018, PASP, 130, 035002 [Google Scholar]
  147. Tokovinin, A., Fischer, D. A., Bonati, M., et al. 2013, PASP, 125, 1336 [NASA ADS] [CrossRef] [Google Scholar]
  148. Tsiaras, A., Waldmann, I. P., Tinetti, G., Tennyson, J., & Yurchenko, S. N. 2019, Nat. Astron., 3, 1086 [Google Scholar]
  149. Twicken, J. D., Catanzarite, J. H., Clarke, B. D., et al. 2018, PASP, 130, 064502 [Google Scholar]
  150. VanderPlas, J. T. 2018, ApJS, 236, 16 [Google Scholar]
  151. Wagenmakers, E.-J. 2007, Psychon. Bull. Rev., 14, 779 [Google Scholar]
  152. West, A. A., Hawley, S. L., Bochanski, J. J., et al. 2008, AJ, 135, 785 [Google Scholar]
  153. Zacharias, N., Finch, C. T., Girard, T. M., et al. 2013, AJ, 145, 44 [Google Scholar]
  154. Ziegler, C., Tokovinin, A., Briceflo, C., et al. 2020, AJ, 159, 19 [Google Scholar]
  155. Zsom, A., Seager, S., de Wit, J., & Stamenkovic, V. 2013, ApJ, 778, 109 [NASA ADS] [CrossRef] [Google Scholar]
1

NASA Exoplanet Archive, Jun 2023, https://exoplanetarchive.ipac.caltech.edu/

3

This analysis builds upon the tres-tools package: https://github.com/mdwarfgeek/tres-tools.

17

The candidate has been submitted as a Community TESS Object of Interest https://exofop.ipac.caltech.edu/tess/target.php?id=166184428

All Tables

Table 1

Properties of the TOI-4336 system.

In the text
Table 2

Properties of TOI-4336 A and TOI-4336 A b based on our global transit model (see Sect. 5).

In the text
Table C.1

Astrometry measurements from WDS (Mason et al. 2001) for AC (left) and AB (right) systems used in the orbital fit.

In the text
Table C.2

Summary of orbital parameters for the TOI-4336 system.

In the text
Table D.1

Summary of the ground-based follow-up observations obtained for the validation of TOI-4336 A b.

In the text
Table E.1

Global model errors and scaling factors for the MCMC analysis of the photometric data.

In the text
Table H.1

Properties of TIC 166184428.02 obtained from the fitting stage of SHERLOCK, see Sect. 5.2.

In the text

All Figures

thumbnail Fig. 1

SpeX spectra of TOI-4336 (TOI-4336 A, top row), TIC 166184426 (TOI-4336 B, middle row), and TIC 166184390 (TOI-4336 C, bottom row) from 2021 Jun. 27 (left column) and 2021 Jun. 28 (right column). The target spectrum (blue) is shown alongside the best-fit spectral template from the IRTF Spectral Library (gray). All spectra are normalized to their flux in the 0.9-1.4 μπι region. Wavelengths with strong telluric absorption are shaded (and largely masked from the spectra), and prominent M-dwarf absorption features are highlighted.
In the text
thumbnail Fig. 2

LDSS3 red optical spectrum of TOI-4336 A (black line), compared to its best-fit M4 template (Kesseli et al. 2017, magenta line). Spectra are normalized in the 7400–7500 Å region, and major absorption features are labeled, including regions of strong telluric absorption (⊕). The inset box shows a close-up of the region encompassing Hα (6563 Å) and Li I (6708 Å) features, neither of which is detected.
In the text
thumbnail Fig. 3

Spectral energy distribution of TOI-4336 A (top panel) and its companion stars TOI-4336 Β (middle panel) and TOI-4336 C (bottom panel). Red symbols represent the observed photometric measurements, the horizontal bars represent the effective width of the passband. Blue symbols are the model fluxes from the best-fit PHOENIX atmosphere model (black). The Gaia spectrophotometry is represented as a grey swathe; a closeup view is shown in the inset.
In the text
thumbnail Fig. 4

TESS target pixel files showing the custom apertures used to extract the photometry of TOI-4336 A b in this work for Sectors 11, 38, and 64.
In the text
thumbnail Fig. 5

TESS 2-min cadence photometry obtained using custom apertures for Sectors 11, 38, and 64. The transits are highlighted in purple.
In the text
thumbnail Fig. 6

Archival images of the field around the TOI-4336 system. From left to right: 1974 image taken with the blue plate of DSS/POSS-II, 1993 image taken with the InfraRed plate of DSS/POSS-II, and 2021 image taken with SSO in Sloan-r′. The white circles indicate the position of the stars on the 2021 images.
In the text
thumbnail Fig. 7

Contrast curve obtained with Zorro for TOI-4336 A in two bands (562 nm and 832 nm) and the reconstructed speckle image of the observation from 2022 May 17.
In the text
thumbnail Fig. 8

Speckle autocorrelation functions (ACF) and 5σ sensitivity functions of the SOAR observations obtained for A (on the left) and B (on the right).
In the text
thumbnail Fig. 9

Chromaticity check of TOI-4336 A b. On the left: phase folded transits from TESS and ground-based observations. The filters above the dashed line correspond to the TESS and TRAPPIST-South observations which are diluted by TOI-4336 B. On the right: comparison of the depths obtained for the chromatic fit. All the bands agree within 1-σ, and the shaded region corresponds to the depth obtained from the achromatic fit for comparison.
In the text
thumbnail Fig. 10

Results of the injection-recovery tests on TESS light curves, described in Sect. 5.3. Almost all planets larger than 1.5 R⊕ were recovered successfully up to 30 d.
In the text
thumbnail Fig. 11

Top left panel: complete sample of all known exoplanets with measured masses (Data from NASA Exoplanet Archive, Apr 2024, https://exoplanetarchive.ipac.caltech.edu/). Following Kempton et al. (2018), terrestrial planets correspond to planetary radii below 1.25 R⊕, small sub-Neptunes radii between 1.25 and 2.75 R⊕, large sub-Neptunes radii between 2.75 and 4 R⊕, and sub-Jovians are planets with radii between 4 and 10 R⊕. TOI-4336 A b is shown by the star symbol. Top right panel: zoomed in view that puts TOI-4336 A b in the context of small sub-Neptunes only. Bottom left panel: stellar effective temperature as a function of insolation of transiting exoplanets orbiting hosts cooler than 4000 K. The solid black line denotes the empirical HZ boundary and the dashed line the conservative one (Kopparapu et al. 2013). The size of points scales with the planetary radius. The points are colored according to their equilibrium temperature. TOI-4336 A b is highlighted with errorbars, and it is well placed at the inner edge of the HZ of its host star. Bottom right panel: synthetic transmission spectra for TOI-4336 A b, assuming a typical H/He-rich atmosphere with isothermal temperature profile at the planet’s equilibrium temperature (~300 K). The simulated data points are for the Wide Field Camera 3 instrument of HST.
In the text
thumbnail Fig. A.1

Color-magnitude diagram for nearby M dwarfs.The TOI-4336 system is marked with a star symbol and color-coded to distinguish the three stars. TOI-4336 C appear redder than A and B, which is consistent with a later spectral type (M4±0.5 compared to M3.5±0.5. We used the extended target list of the SPECULOOS survey which gathers over 14 000 M-dwarf stars within 40 parsecs to generate the diagram (Sebastian et al. 2021).
In the text
thumbnail Fig. B.1

Age distribution of the GALAH survey (light grey) and the distribution of stars within our selection criteria (dark grey). The yellow line indicates the median age and the yellow bands are the 16 percentile and 84 percentile regions. We estimate an age of 6.7−3.1+2.7$6.7_{ - 3.1}^{ + 2.7}$ Gyr for the system.
In the text
thumbnail Fig. C.1

Top row: AB (left) and AC (right) orbits with ten random orbital solutions from posterior of the orbit fit. Middle row: separation with WDS astrometry for AB (left) and AC (right) orbits with a hundred random orbital solutions from posterior of the fit.Bottom row: the same for the position angle.
In the text
thumbnail Fig. C.2

Posterior samples of the AB system. Histogram sub-panels show the posterior distribution, with the median and 68% confidence intervals marked by dashed lines, with titles quantifying those ranges.
In the text
thumbnail Fig. F.1

Posterior distributions for the main jump parameters used in the emcee fit (see Sect. 5): the effective temperature (Teff), the metallicity ([Fe/H]), the log of the stellar density (log ρ*), the log of the stellar mass (log M*), the dilution in the TESS and z′ bands, the epoch (T0), the transit depth (dF), the log of the orbital period (log P) and the cosine of the orbital inclination of the planet (cos i).
In the text
thumbnail Fig. G.1

O-C diagram obtained from a global fit allowing transit timing variations. We used as reference the timing obtained in the global transit model of Table 2.
In the text
thumbnail Fig. H.1

Transit model of TIC 166184428.02. Top panel: phase folded photometry of TIC 166184428.02. The best-fit model is shown in solid black. The grey points are the raw flux and the red points are the 15-minutes binned flux. Bottom panel: residuals of the best-fit model.
In the text
thumbnail Fig. H.2

Posterior distributions for the main parameters fitted in the Nested Sampling fit for the second candidate found in the TOI-4336 system, TIC 166184428.02 (see Sect. 5.2), with Rb the candidate radius, R* the stellar radius, ib the inclination, T0;b the epoch, Pb the period, q1:1c and q2:1c the Kipping parametrisation (Kipping 2013) for the limb darkening coefficients, logσlc the error scaling parameter, and lnσ(lc) and lnp(lc) the two GP parameters related to the amplitude and length scale.
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.