A PSF-based Approach to TESS High quality data Of Stellar clusters (PATHOS) – IV. Candidate exoplanets around stars in open clusters: frequency and age–planetary radius distribution Free

ABSTRACT

The knowledge of the ages of stars hosting exoplanets allows us to obtain an overview on the evolution of exoplanets and understand the mechanisms affecting their life. The measurement of the ages of stars in the Galaxy is usually affected by large uncertainties. An exception are the stellar clusters: For their coeval members, born from the same molecular cloud, ages can be measured with extreme accuracy. In this context, the project PATHOS is providing candidate exoplanets orbiting members of stellar clusters and associations through the analysis of high-precision light curves obtained with cutting-edge tools. In this work, we exploited the data collected during the second year of the Transiting Exoplanet Survey Satellite mission. We extracted, analysed, and modelled the light curves of |$\sim 90\, 000$| stars in open clusters located in the Northern ecliptic hemisphere in order to find candidate exoplanets. We measured the frequencies of candidate exoplanets in open clusters for different orbital periods and planetary radii, taking into account the detection efficiency of our pipeline and the false positive probabilities of our candidates. We analysed the age–R_P distribution of candidate and confirmed exoplanets with periods <100 d and well constrained ages. While no peculiar trends are observed for Jupiter-size and (super-)Earth-size planets, we found that objects with |$4 \lesssim R_{\rm P} \lesssim 13R_{\rm Earth}$| are concentrated at ages ≲200 Myr; different scenarios (atmospheric losses, migration, etc.) are considered to explain the observed age–R_P distribution.

1 INTRODUCTION

In Summer 2020, the Transiting Exoplanet Survey Satellite (TESS; Ricker et al. 2015) concluded its main mission after about 2 yr of observations. In this period, the spacecraft has observed millions of stars in about ≳70 per cent of the sky with an unprecedented photometric precision and temporal coverage, and new data from the extended mission, characterized by (in part) a new observing strategy, are coming.

Stellar clusters and associations offer the unique opportunity to derive precise stellar parameters (like radius, mass, chemical content, and especially age) for their members simply using theoretical models. During the main mission TESS observed many hundreds stellar (open and globular) clusters and associations in sectors of ∼27 d. However, the low resolution of the four cameras (∼21 arcsec pixel⁻¹) makes the extraction of high-precision light curves difficult for the stars located in these dense regions.

The project ‘a PSF-based Approach to TESS High quality data Of Stellar clusters’ (PATHOS; Nardiello et al. 2019, hereafter Paper I) was born to exploit the TESS data in order to extract high-precision photometry for members of stellar clusters and associations adopting an innovative approach, based on the use of empirical point spread functions (PSFs) and neighbour subtraction. Scope of the project is the discovery and characterization of new candidate exoplanets in stellar clusters and the analysis of possible correlations between well-measured star properties and candidate exoplanet characteristics. Field stars’ age is usually affected by large uncertainties, but it is also an essential information to constrain the formation and understand the evolution of exoplanets, like for example how and on which temporal scales the mechanisms that bring to the atmosphere evaporation of low-mass close-in exoplanets happen (Lammer et al. 2003; Baraffe et al. 2005; Murray-Clay, Chiang & Murray 2009; Owen & Jackson 2012; Owen & Wu 2017; Wu 2019; Owen 2020). In this context, the PATHOS project is providing interesting candidate exoplanets orbiting stars with well constrained ages. Moreover, the high-precision light curves generated in our project, and publicly available to the astronomical community, allow us to obtain results not only in the research field of exoplanets, but also in other fields like asteroseismology (Mackereth et al. 2021), or the analysis of the spin axis orientations of cluster members (Healy & McCullough 2020).

We already successfully applied the PATHOS pipeline in Paper I, when we studied the stars in an extremely crowded region containing the globular cluster 47 Tuc and the Small Magellanic Cloud. In Nardiello et al. (2020, hereafter Paper II), we extracted and analysed the light curves of open cluster members located in the Southern ecliptic hemisphere, finding 33 objects of interest and deriving a first estimate of exoplanet frequency in open clusters. Nardiello (2020, hereafter Paper III) studied the light curves of the members of five young associations, having ages ≲10 Myr; in particular, the author performed a gyrochronological analysis of association members to constrain the age of the stars, analysed the dust in the circumstellar discs of the young members and identified and characterized six strong candidate exoplanets.

In this work, we exploited the TESS data collected during Cycle 2 (Sectors 14–26) to obtain high-precision light curves of cluster members in the Northern ecliptic hemisphere by using our cutting-edge tools (Section 2), find and characterize candidate exoplanets in stellar clusters (Section 3), and analyse their frequency and properties as a function of host stars’ characteristics (Sections 4). We summarized and discussed the joined results obtained in this work and in Paper II in Section 5.

2 OBSERVATION AND DATA REDUCTION

In this work, we extracted and analysed the light curves of the stars likely members of Northern ecliptic hemisphere open clusters observed by TESS during the second year of the mission. The observations used in this work were carried out between 2019 July 18 and 2020 July 4 (∼352 d), and are divided into 13 sectors (Sectors 14–26); in Sectors 21, 22, and 23, no open clusters fell in the TESS field of view and therefore the final number of analysed sectors is 10.

For the light-curve extraction and correction, we used the PATHOS pipeline described in detail in Papers I and II. We extracted the light curves of stars in a given catalogue from TESS full frame images (FFIs) by using the light-curve extractor IMG2LC. This software was developed by Nardiello et al. (2015a, 2016a) for ground-based observations, and it is a versatile tool that can be used with photometric time-series collected also with space-based observatories (see, e.g. Libralato et al. 2016a,b; Nardiello et al. 2016b).

The three main inputs of our PSF-based approach are (i) FFIs, (ii) PSFs, and (iii) input catalogue. For each star in the input catalogue, the light-curve extractor searches for the neighbours within a radius of 20 TESS pixels in the Gaia DR2 catalogue (Gaia Collaboration et al. 2018), transforms their positions and luminosities in the reference system of the FFI, models them by using a local PSF and then subtracts them from the FFI. Finally, it extracts PSF-fitting and aperture (1-, 2-, 3-, and 4-pixel radius) photometries of the target star from the neighbour-subtracted FFI. This approach has two advantages: (i) It minimizes the dilution effects due to the neighbour contaminants, and (ii) it allows us the extraction of high-precision photometry for stars in the TESS faint regime of magnitudes (T ≳ 15; see, e.g. Apai, Nardiello & Bedin 2021). We corrected the extracted raw light curves for systematic effects by fitting and applying the Cotrending Basis Vectors, as widely discussed in Papers I and II.

As in Paper II, we used as input list the catalogue of cluster members published by Cantat-Gaudin et al. (2018); this catalogue contains the positions, colours, magnitudes, proper motions, parallaxes, and membership probabilities of likely members in 1229 stellar clusters. From this catalogue, we selected all the stars that satisfy these two conditions: (i) magnitude G < 17.5, because stars with larger magnitude are too faint to be detected by TESS; and (ii) ecliptic latitude β > 4^○, which corresponds to the part of the Northern ecliptic hemisphere covered by TESS.¹

Fig. 1 shows the 126 372 stars (red points) in the input catalogue overlapped with the TESS fields of view (grey points): About 1/3 of them fall outside the TESS observations. We extracted 150 216 light curves of 89 858 stars in 411 clusters; about 50.3  per cent, i.e. 45 182 stars, were observed in only one sector, 30 957 stars (∼34.4  per cent) were observed in two sectors, 11 844 (∼13.2  per cent) in three sectors, and 1 875 (∼2.1  per cent) in four or more sectors.

Light curves are released on the Mikulski Archive for Space Telescopes (MAST) as a High Level Science Product (HLSP) under the project PATHOS² (DOI: 10.17909/t9-es7m-vw14). A detailed description of the light curves (that are both in ascii and fits format) is reported in Papers I and II and in the MAST web page of the PATHOS project.

2.1 Photometric precision

We explored two different quality parameters, already defined in Papers I, II, and III, to identify for each star the photometric method that gives the best light curve.

The first quality parameter is the simple rms, defined as the 68.27th percentile of the 3.5σ-clipped sorted residual from the median value. This parameter is sensitive to the (high) variability of some stars, and, for this reason, is not recommended to estimate the photometric precision of the light curve. The mean trends of the rms as a function of the TESS magnitude T for the five photometric methods are reported in Fig. 2 (top panel).

The second quality parameter is the P2P rms, defined as the 68.27th percentile of the 3.5σ-clipped sorted residual from the median value of the vector δF_j = F_j − F_{j + 1}, with F the flux at a given epoch j. This parameter is not sensitive to the intrinsic luminosity variations of the stars, and we used it to define the interval in which each photometric method works, on average, better than the others. From the mean trends shown in bottom panel of Fig. 2, we found (confirming the results obtained in the previous works) that for stars with 5.5 ≲ T ≲ 7.0, the aperture photometry with radius 4-pixel gives the best results; the best photometric methods for the intervals 7.0 ≲ T ≲ 9.0, 9.0 ≲ T ≲ 10.0, and 10.0 ≲ T ≲ 13.0 are 3-pixel, 2-pixel aperture, and PSF-fitting photometries, respectively. For faint stars with T ≳ 13.0, the 1-pixel aperture photometry gives the lower P2P rms.

In the following analysis, we used, for each star of magnitude T_⋆, the light curve associated with the photometric method that has the lower mean P2P rms in T_⋆. We excluded from the analysis all the stars whose mean light-curve instrumental magnitude (T_instr) is too different from the expected magnitude T_calib, following this procedure: We extracted the δT = T_instr − T_calib distribution, we calculated its mean (⁠|$\bar{\delta T}$|⁠) and the standard deviation (σ_δT), and we excluded the ith light curve if |$|\delta T_i-\bar{\delta T}|\gt 4\sigma _{\delta T}$|⁠. We also excluded all the light curves that have |$\lt 75\, {{\ \rm per\ cent}}$| of well-measured points (i.e. DQUALITY=0 and FLUX≠ 0). The final number of analysed light curves is 138 924 associated with 84 967 stars.

3 CANDIDATE EXOPLANETS: SEARCHING, VETTING, AND CHARACTERIZATION

We searched for signals of transiting objects among the selected light curves following the procedure described in Papers II and III. Briefly, we removed the intrinsic stellar variability interpolating the light curve with a fifth order spline defined on N_knots knots. In order to model short- and long-period variability, we considered two different grids of knots, with knots spaced by 6.5 and 13.0 h, respectively. The grids of knots are defined on continuous parts of light curves that does not present ‘breaks’ >0.5 d, in order to avoid the introduction of artifacts in the flattened light curve. We also removed from the light curves the photometric points associated with high values of the local sky (>5σ_SKY from the mean local background), with DQUALITY>0 and 4σ above the median normalized flux. We extracted the transit-fitting least-squares (TLS) periodograms of the flattened light curves adopting the Python package tls³ (Hippke & Heller 2019), and searched for transit signals with period 0.6 d ≤ P ≤ T_LC, with T_LC being the temporal length of the light curve. We used the output parameters for the first selection of candidate transiting objects, as follows: (i) We selected the stars associated with a depth of the transit |$\delta _{\rm t}\lt 10\, {{\ \rm per\ cent}}$| and to a significance between odd and even transits σ_odd-even < 2.5 ; (ii) we divided the signal detection efficiency (SDE) and the signal-to-noise ratio (SNR) distributions in bins of periods δP = 0.5 d, and calculated the 3.5σ-clipped mean and standard deviation of SDE (⁠|$\bar{\rm SDE}$| and σ_SDE) and SNR (⁠|$\bar{\rm SNR}$| and σ_SNR) within each bin. Then we interpolated the binned points |$\bar{\rm SDE} + 3.5\sigma _{\rm SDE}$| and |$\bar{\rm SNR} + 3.5\sigma _{\rm SNR}$| with a spline, and we selected as good candidates all the stars above these splines. Panels (a) of Fig. 3 show an example of selections based on the output parameters of TLS. We visually inspected the light curves that passed the above selections to exclude false positives due to the presence of artefacts. We applied the above procedure both to the light curves obtained in each single sector, and then to the stacked light curves of stars observed in more than one sector, in order to increase the efficiency of transit detection. The number of stars that passed the first-step selection is 279 (∼0.3 per cent).

These candidates were subjected to a series of vetting tests widely described in the previous papers of the PATHOS series. They are as follows: (i) Check for the transit depths in the light curves obtained with different photometric methods; (ii) check for the presence of secondary eclipses in the light curves phased with a period 0.5 × P, 1.0 × P, and 2.0 × P, with P the period found by the TLS routine; (iii) check for the transit depths by comparing binned odd/even transits (panels b of Fig. 3); and (iv) check for contamination through the analysis of the in/out-of-transit difference centroid (panel c of Fig. 3). After this second-step selection, 39 transiting objects of interest survived (∼0.05  per cent); one of these objects showed only one transit in its light curve.

3.1 TESS Objects of Interest

We cross-matched the TESS Objects of Interest (TOI) list⁴ with our input catalogue of cluster members. Four candidates, found by the Quick-Look Pipeline (QLP; Huang et al. 2020), are also in our input catalogue, but only one of them (TOI-1535, TIC 420288086) is in our final list of transiting objects of interest. Two of them (TOI-1497, TIC 371673488 and TOI-1321, TIC 195199644) were not detected by our pipeline because no transit signals are present in the light curves we analysed, even if the mean scatter of our light curves is lower than the scatter of the light curves shown in the QLP data validation report. We checked the notes about these two candidates on ExoFOP⁵: (i) the depth-aperture correlation for TOI-1497.01, reported in a note on ExoFOP, is confirmed by our analysis; (ii) for TOI-1321, a depth-aperture correlation is also reported; moreover, a note associated with a photometric follow-up with MuSCAT2 reports a deep transit signal from a nearby star at ∼1 arcmin from TOI-1321. Therefore, both these candidates are likely contaminated by neighbour sources. The fourth QLP candidate, TOI-1188 (TIC 15247229), was excluded from our final list after the centroid analysis. Its vetting tests are reported in Fig. 3. Our conclusion is also supported by the notes reported in the photometric follow-up section of the ExoFOP website.

3.2 Stellar parameters

We fitted theoretical models from the last release of BaSTI (‘a Bag of Stellar Tracks and Isochrones’) models (Hidalgo et al. 2018) to the CMDs of the 32 open clusters that host the stars associated with the 39 transiting objects of interest. In this way, we were able to extract primary information (stellar radius, mass, density, effective temperature) of the stars that host candidate transiting objects. Because metallicity measurements are not available for the large part of the clusters and because open clusters have, on average, metallicities similar to that of the Sun, in our fit, we used isochrones with [Fe/H] = 0.0 ± 0.3 as already done in Paper II, and we added the contribution of the uncertainties on the metallicity to the final errors on the stellar parameter estimates.

We transformed the isochrones from the theoretical to the observational plane using the distance modulus of the clusters obtained by Cantat-Gaudin et al. (2018), and the reddening and ages measured by Kharchenko et al. (2016), Röser et al. (2016), and Bossini et al. (2019); since some of the catalogues do not provides error estimates, for homogeneity, we used a conservative error of 10 per cent on the age and reddening values. Gulliver 49 is an open cluster discovered by Cantat-Gaudin et al. (2018), and no age estimate is provided in literature: We followed the technique by Nardiello et al. (2015b) based on the use of the χ²-minimization between isochrones and fiducial lines to derive an estimate of the cluster age. We found an age of 200 ± 20 Myr.

Clusters’ parameters adopted for the isochrone fitting are reported in Table 1. Stellar parameters of the transiting candidates’ hosts obtained from isochrone fitting were used as priors in transit modelling described in the next section, and are reported in Table A1.

3.3 Transit modelling

We modelled the transits of the objects of interest using the python package pyorbit⁶(Malavolta et al. 2016, 2018, see also Benatti et al. 2019; Carleo et al. 2021; Lacedelli et al. 2021), based on the combined use of the package batman (Kreidberg 2015), the global optimization algorithm pyde⁷ (Storn & Price 1997), and the affine invariant Markov chain Monte Carlo sampler emcee (Foreman-Mackey et al. 2013).

For the transit modelling, we included the central time of the first transit (T₀), the period (P), the impact parameter (b), the planetary-to-stellar-radius ratio (R_P/R_⋆), the stellar density (ρ_⋆), and the dilution factor (df). The latter quantity is included as a free parameter, with a Gaussian prior obtained considering all the stars in the Gaia DR2 catalogue that fall in the same pixel of the target,⁸ and transforming their Gaia magnitudes in TESS magnitudes adopting the equations by Stassun et al. (2019). Host star parameters, like the stellar radius (R_⋆), mass radius (M_⋆), gravity (log g), and effective temperature (T_eff) come from the isochrone fits described in the previous section. On the basis of log g and T_eff, we obtained information on the limb-darkening (LD) by using the grid of values published by Claret (2018); we adopted the LD parametrization described by Kipping (2013). In the modelling process, the routine takes into account the local variability of the star by fitting a second-degree polynomial to the out-of-transit part of the light curve. The routine modelled the transits with a fixed circular orbital eccentricity (e = 0), and taking into account the 30-min cadence of the TESS FFIs (Kipping 2010).

The adopted priors on stellar parameters are reported in Table A1. The package pyorbit explored all the parameters in linear space. In the emcee run, the number of walkers used is 10 times the number of free parameters. We ran, for each model, the sampler for 80 000 steps, excluding the first 15 000 steps as burn-in and using a thinning factor of 100. Fig. 4 shows an example of the modelling process in the case of PATHOS-53, a mono-transit object of interest.

In Table 2, we report the results of the transit fitting; Figs A1–A3 give an overview on the main properties of each transiting object of interest (position on the CMD, proper motions, in-/out-of-transit centroid analysis, transit modelling).

4 RESULTS

We improved the calculation of the frequencies of exoplanets reported in Paper II, taking into consideration the detection efficiency of our method and the false positive probabilities (FPPs) of our candidates. In this analysis, we considered candidates identified in this work and in Paper II, excluding all the PATHOS objects with R_P > 2.5R_J (because of their doubtful planet nature): 23 candidates survived after this selection (14 and 9 in the Southern and Northern ecliptic hemispheres, respectively).

4.1 Detection efficiency

We calculated the detection efficiency of our finding pipeline injecting transit signals of planets having radii R_P in a sample of light curves extracted randomly from the collection of light curves analysed in Paper II and in this work.

To inject the transit signals in each light curve, we roughly estimated the stellar radius and mass of each the 232 669 stars analysed in Paper II and in this work by using their absolute M_{G, 0} versus (G_BP − G_RP)₀ CMD. To obtain the absolute CMD, we transformed the apparent magnitudes into absolute magnitudes by using the Gaia DR2 distances obtained by Bailer-Jones et al. (2018). For each star, we corrected the effects of the extinction on the colour and the magnitude of the stars calculating the reddening value [E(B − V)] by using the python routine mwdust⁹ (Bovy et al. 2016) and the Combined19 dustmap (Drimmel, Cabrera-Lavers & López-Corredoira 2003; Marshall et al. 2006; Green et al. 2019), and the colour-dependent equation and the coefficients reported by Bossini et al. (2019). Fig. 5 shows the M_{G, 0} versus (G_BP − G_RP)₀ CMD for all the stars in our sample; red starred points indicate the 23 candidates. We selected the likely main-sequence stars as follows: We performed a first guess selection of the main-sequence stars by hand, excluding all the stars that clearly are evolved stars (subgiant and red-giant stars); in the second step, we calculated the fiducial line of the likely main-sequence stars by using the naive estimator (Silverman 1986; see also Nardiello et al. 2015b for the description of the method) and we selected the stars whose colours are within 2σ from the mean colour of the fiducial line. We calculated the radius and the mass of each star in our catalogue by using the main-sequence points of the PARSEC (Girardi et al. 2002; Bressan et al. 2012; Marigo et al. 2017) isochrones with ages between 10 and 1000 Myr. In Fig. 5, we reported as reference the mean M_{G, 0} absolute magnitude for stars having radius |$R_{\star }=0.5, 1.0, 1.5, \ {\rm and} \ 2.5\, R_{\odot }$|⁠.

Because our candidate exoplanets orbit stars with R_⋆ ≲ 2.5 R_⊙, we first selected the light curves associated with these stars; we divided our sample of light curves on the basis of the TESS magnitude of the associated stars, considering bins of size ΔT = 1.0 in the magnitude range 6.0 ≤ T ≤ 18.0. In each bin, we randomly caught 350 light curves¹⁰ associated with stars having the previously estimated radii and masses, and, by using the routine INJECTTRANSIT implemented in vartools v.1.39 (Hartman & Bakos 2016), we injected in each light curve a periodic transit signal of a planet having radius |$R^{\rm inj}_{\rm P}$|⁠, orbital period P^inj, random inclination i^inj (with the constraint that there must be a transit), and eccentricity e = 0. We considered 12 different cases in which |$R^{\rm inj}_{\rm P}$| and P^inj randomly vary between a given minimum and a maximum; (i) we first injected (super-)Earth-size planets having radii |$0.85 \le R^{\rm inj}_{\rm P} \le 3.9R_{\rm Earth}$|⁠; (ii) the second case we considered was for (super-)Neptune planets having radii |$3.9 \lt R^{\rm inj}_{\rm P} \le 11.2R_{\rm Earth} \sim 1 R_{\rm J}$|⁠; finally, we considered giant planets with radii |$1.0 \lt R^{\rm inj}_{\rm P} \le 2.5R_{\rm J}$|⁠. For each |$R^{\rm inj}_{\rm P}$| case, we considered random orbital periods in the following intervals: (a) short periods |$0.5 \le P^{\rm inj} \le 2.0\, {\rm d}$|⁠; (b) short-medium periods |$2.0 \lt P^{\rm inj} \le 10.0\, {\rm d}$|⁠; (c) long-medium periods |$10.0 \lt P^{\rm inj} \le 85.0\, {\rm d}$|⁠; and (d) long periods |$85.0 \lt P^{\rm inj} \le 365.0\, {\rm d}$|⁠. We have not made any selection on the time length of the light curves, so, on the basis of the associated star observability, the light curve can span randomly between ∼27 and ∼365 d.

We injected the transit signals in the raw light curves and we followed the same pipeline for the correction of the light curves and the detection of transit signals described in Sections 2 and 3. We considered the injected planet as recovered if it passes the selections illustrated in Section 3 and in Fig. 3 and if |$|P^{\rm inj}-P^{\rm inj}_{\rm TLS}|\lt 4.0\sigma (P^{\rm inj}_{\rm TLS})$|⁠, where |$P^{\rm inj}_{\rm TLS}$| is the period obtained by the extraction of the TLS periodogram, and |$\sigma (P^{\rm inj}_{\rm TLS})$| its error. We also considered the planet as recovered if |$P^{\rm inj}_{\rm TLS}$| is equal to 0.5 × P^inj or 2.0 × P^inj, within |$4.0\sigma (P^{\rm inj}_{\rm TLS})$|⁠. On a sample of 100 light curves selected randomly from those that passed the selections, we performed the vetting tests described in Section 3¹¹: all the selected objects passed the tests. We finally calculated the detection efficiency in each magnitude interval ΔT as the ratio N_rec/N_inj, where N_rec is the number of simulated planets recovered, and N_inj the number of light curves in which we injected transit signals. In Fig. 6 and Table 3, we reported the detection efficiency in the P^inj versus T grid for the three different planetary radius intervals: For Earth-size planets with R_P ≲ 3.9R_Earth, the detection efficiency is low (⁠|$\lesssim 40~{{\ \rm per\ cent}}$|⁠), mainly because the large part of low-mass stars, for which the detection of small size planets is easier, is concentrated at magnitudes T ≳ 11 (see fig. 8 of Paper II); for short period Neptune-size planets (P < 10 d) orbiting stars with magnitude T ≲ 12, the detection efficiency is ≳30–|$50{{\ \rm per\ cent}}$|⁠, and decreases at fainter magnitude to 20–30 per cent; we obtained a detection efficiency always |$\gtrsim 30~{{\ \rm per\ cent}}$| for giant planets with periods ≲10 d. For periods P ≳ 10 d, the detection efficiencies are always |$\lesssim 20~{{\ \rm per\ cent}}$|⁠: This is an effect caused by the fact that only |$\sim 15{{\ \rm per\ cent}}$| of the stars in our sample are observed in more than two sectors, the detection of long-period planets difficult. Finally, even if we used different grids to flat the light curves to take into account the different kinds of variability, the detection efficiency can be lower for very active stars in young stellar clusters.

4.2 FPP estimation

We used the tool vespa v. 0.6¹² (Morton 2012, 2015) to estimate the FPP in a Bayesian framework of each of the 23 candidate exoplanets previously selected. This software estimates the probability that the transit signal detected in a light curve is due to a real planet or to a different source. Different scenarios are taken into consideration to explain the signal: (i) simple eclipsing binary (EB, single/double period); (ii) a hierarchical triple system where two components eclipse (HEB, single/double period); (iii) a blended eclipsing binary inside the photometric aperture of the target (BEB, single/double period); and (iv) a transiting planet on the target star (pl). We want to highlight that vespa does not consider the scenario in which transit signals are due to non-astrophysical sources, even if in our case the probability it happens is low, because the analysed candidates passed a series of vetting tests that allow to exclude signals due to systematic effects. For each configuration, vespa simulates a representative stellar population, constrained by the information we have from the isochrone fitting and the light-curve modelling; in particular, we put constraints on the G, G_BP, G_RP, J_2MASS, H_2MASS, K_2MASS magnitudes of the target star, the (α, δ) coordinates of the target, the effective temperature T_eff, the density ρ_⋆, and the surface gravity log (g) of the star, the orbital period P of the candidate exoplanet, the planet-to-star radius ratio R_P/R_⋆. Moreover, we gave as an input the light curve of the stars hosting the candidates. The software takes into account different false positive scenarios for those populations, and uses them to define a prior likelihood that a specific configuration actually exists and the likelihood of transit for those configurations. On the basis of the results obtained for different scenarios, at the end, it calculates the FPP that the transit signals are due to a false positive: The lower the FPP, the larger is the probability that the signal is due to a planet. We reported in Table 4 the FPPs and the most likely scenario (with the respective probability) derived with vespa for the 23 candidate exoplanets: ∼40 per cent of them have a probability ≳50 per cent to be a real exoplanet. It means that ∼3/5 of the signals we detected are likely false positives.

To test the reliability of the FPP estimation performed by vespa, we used this tool to validate 10 simulated transiting exoplanets extracted randomly from the sample of the recovered candidates described in Section 4.1 (considering all the scenarios). We ran vespa on the simulated exoplanets following the same procedure we adopted for the real candidate exoplanets. In an ideal case, the sum of all the FPPs obtained for the simulated exoplanets is equal to 0, but we obtained that ∑FPP_i ∼ 0.5 ; it means that, on average, the measured FPP is overestimated of |$\sim 5~{{\ \rm per\ cent}}$|⁠. We used this result to correct the FPPs measured for the real candidate exoplanets, as reported in Table 4 (FPP_c).

We also reported in Table 4 the Renormalized Unit Weight Error (RUWE) index calculated as described by Belokurov et al. (2020); its value is ∼1 for sources that can be fitted with single-star models, and its value is larger (≳1.4) if the astrometric fit of the source by using single-star model is not the best solution. Among the candidates labelled as ‘pl’, only PATHOS-33 shows a high RUWE value (∼2.6), even if the FPP obtained with vespa is ∼0.

4.3 Candidate exoplanets’ frequency in open clusters

4.4 Age–planetary radius distribution

Fig. 7 shows the stellar age versus planetary radius distribution for the candidate transiting exoplanets in open clusters identified in this work and in Papers II, and reported as ‘pl’ in Table 4, and the candidates around stars in young associations identified in Paper III; the colours indicate the insolation flux S calculated on the basis of the stellar radius and effective temperature and of the planetary semimajor axis. Moreover, we included confirmed and candidate exoplanets from literature, orbiting stars in stellar clusters and associations, i.e. with a well-constrained age (coloured triangles). The list of literature objects plotted in Fig. 7 and the corresponding references are reported in Table 6. Finally, we also added objects (coloured squares, Nardiello et al. in preparation) that are under investigation in the context of the ‘GAPS Young Objects’ (GAPS-YO) project (Carleo et al. 2020), aimed at the monitoring of young and intermediate-age stars for the discovery and characterization of young planets. All the objects showed in Fig. 7 have orbital periods <100 d, i.e semimajor axis a ≲ 0.5 au.

We found that Jupiters with R_P ≳ 1R_J are distributed randomly between ∼10 Myr and ∼10 Gyr. Objects having Neptune-sizes or smaller (Earths/super-Earths) are concentrated at ages >200–300 Myr; anyway the lack of this kind of objects around young stars might be an observational bias due to the difficulty of detecting their transits in the highly variable light curves of active young stars. Objects having a planetary radius |$4 \lesssim R_{\rm P} \lesssim 10R_{\rm Earth}$| are concentrated at ages <100–200 Myr. Is it an observational bias, or the lack of super-Neptune/sub-Jovian-size planets orbiting (on short periods) stars with (well-measured) ages ≳200 Myr is due to an evolutionary effect of the planets?

The majority of low-mass close-in exoplanets (with R_P= 1–4R_Earth) detected until today orbit field stars with ages ≳1–3 Gyr. Some of these planets are almost totally rocky (e.g. Kepler-93; Dressing et al. 2015), others have low densities that can be explained by the presence of an extended H/He atmosphere. Among them, the Kepler-36 system is particular because formed by two planets with very similar semimajor axis, but totally different densities (Carter et al. 2012), with the inner planet less massive than the external planet. Lopez & Fortney (2013) proved that both the exoplanets in the Kepler-36 system were born with H/He atmospheres and were more massive in the early stages of their life, and that, given the lower core mass of the inner planet, the latter has lost large part of its atmosphere, despite the outer planet that was able to retain about half of its initial atmosphere. As explained in detail in the review by Owen (2019, see also references therein), there is strong evidence that atmospheric escape is the mechanism that prevails in the first stages of a low-mass close-in exoplanet’s evolution, and it depends (in a first approximation) on the characteristics of the host star (and its high-energy emissions), the distance of the planet from the star, and its core mass. During their formation, these planets accreted large amount of H/He into their expanded atmospheres, inflating their radius at 5–13R_Earth; in few hundreds Myr, if irradiated strongly enough, these planets lose the large part of their atmosphere, resulting in (super-)Earth/sub-Neptune-size planets (Owen et al. 2020). Even if Fig. 7 does not show any particular dependence on the insolation flux, the distribution of short-period (candidate) exoplanets with radii |$1 13 \lesssim R_{\rm P}\lesssim 13R_{\rm Earth}$| seems to confirm the idea of atmospheric escape on timescales of ≲100–200 Myr; anyway, as demonstrated by Owen (2020), even if challenging, for many of these candidates, mass measurements are mandatory in order to understand the mechanisms and constraint the time-scales of planetary atmospheric evolution.

Other possible explanations for the lack of planets older than ∼200 Myr in the interval 4 ≲ R_P ≲ 10R_Earth may be linked to the dynamical evolution of this kind of planets after their formation, such as migration due to interaction with the protoplanetary disc or planet-planet scattering (Fabrycky & Tremaine 2007; Chatterjee et al. 2008; Kley & Nelson 2012). These phenomena are expected to modify the orbital characteristics of the planets, such as inclination, eccentricity, and semimajor axis, making the detection of a possible transit of these planets more difficult. These phenomena occur on scales of a few tens of Myr or less, and may not yet have come into play for the objects reported in this analysis.

5 SUMMARY AND CONCLUSION

The aim of the PATHOS project is the discovery and first characterization of transiting objects around stars in stellar clusters and associations observed by TESS. Stellar clusters and associations offer the rare opportunity to obtain precise measurements of the ages of the stars (usually affected by large uncertainties), in addition to stars’ physical parameters such as radius, mass, effective temperature, etc., and analyse planet characteristics as a function of its host’s properties. Because of the low resolution of TESS cameras, stellar clusters appear as very crowded regions on FFIs, and appropriate tools are necessary to obtain high-precision light curves for cluster members. In order to obtain the best light curves, we developed a cutting-edge technique for the extraction of high-precision photometry of stars in crowded fields, based on the use of empirical PSFs and neighbour subtraction, that allows us to minimize neighbour contamination and extract light curves for very faint objects.

In this work, the fourth of the PATHOS series, we extracted and corrected 150 216 light curves of 89 858 open cluster members listed in the Cantat-Gaudin et al. (2018)’s catalogue and observed during the Cycle 2 (Sectors 14–26) of the TESS mission. By using the pipeline already tested in the previous papers of the series, we searched for transit signals among the light curves, finding 39 transiting objects of interest, which are added to the 33 objects orbiting open cluster members identified in Paper II. We modelled their light curves to extract planet parameters. From the two lists of objects of interest, we isolated 23 candidates with planetary radius |$R_{\rm P} \lesssim 2.5R_{\rm J}$|⁠, and we calculated their FPPs to be a planet, finding that about 3/5 of them are likely false positives.

Taking into account of the transit detection efficiency of our pipeline and of the rate of false positives, we calculated the frequencies of candidate exoplanets in open clusters (f_⋆) for different (Earth/Neptune/Jupiter) sizes planets in different orbital periods intervals (between 0.5 and 365 d). Because we did not detect any strong candidate with orbital period <2 d, we calculated the frequency upper limit, finding f_⋆ ≲ 0.06, ≲ 0.007, and ≲0.004 per cent for Earth-, Neptune-, and Jupiter-size candidate exoplanets, respectively. For candidate exoplanets with periods 2.0 < P < 10.0 d, we found for Earth-size objects, an upper limit f_⋆ ≲ 0.2 per cent, while for candidate planets with R_P ≳ 3.9R_Earth, the frequency is ∼0.06 per cent; in the range 10.0 < P < 85.0 d, we obtained that the frequency of Earth and Neptune-size candidates is |$\lesssim 2.~{{\ \rm per\ cent}}$|⁠, while the measured frequency for giant planets is |$\sim 1.6~{{\ \rm per\ cent}}$|⁠; finally, for long period planets (>85 d), we found that the frequency of candidate exoplanets is |$\lt 5\!-\!7~{{\ \rm per\ cent}}$|⁠. The large part of the measured and upper limit frequencies are lower than the values reported by Fressin et al. (2013) in the same period intervals for the corresponding size exoplanets around field stars (even if bias can be introduced by the poor statistic of our sample); an exception are the giant planets with 10–85 d period: We found a mean frequency of 1.6 ± 0.9 per cent, in agreement with the value obtained by Fressin et al. (2013, |$1.5 \pm 0.2~{{\ \rm per\ cent}}$|⁠), but also by Santerne et al. (2016, f_⋆ = 0.90 ± 0.26 per cent) and Deleuil et al. (2018, f_⋆ = 1.86 ± 0.68 per cent).

We investigated the stellar age versus planetary radius distribution by using the results obtained in the PATHOS project and the results obtained in other works for exoplanets orbiting stars with well constrained ages. We only used planets with orbital periods <100 d and we divided the distribution on the basis of the planetary radius R_P: (1) for (candidate) exoplanets with R_P ≳ 1R_J, we found no particular trends, and the planets are randomly distributed on over the range of ages; (2) objects with R_P ≲ 4R_Earth are concentrated at ages >100–200 Myr, but it might be an observational bias due to the difficulty to detect (super-)Earth transits in the light curve of young active stars; and (3) in the range |$4 \lesssim R_{\rm P} \lesssim 13R_{\rm Earth}$|⁠, there is a concentration of objects around young stars with ages <100 Myr. A possible explanation of such concentration is that these objects are young planets with rocky cores that, in the early stages of their formation, have accreted large amount of hydrogen and helium in their atmospheres, inflating their radius. On the time-scales of ∼100–200 Myr, they lose large part of their atmosphere mainly because of the strong irradiation of the host star (see Owen 2019 for a review of all the mechanisms that contribute to the atmosphere escape), and within few hundreds Myr, their radius decreases to that typical of (super-)Earth/Neptune planets. Other explanations can be related to the dynamical evolution of these exoplanets (planet-planet scattering, migration, etc.).

For this reason, the analysis of light curves of members of young associations (ages ≲200–300 Myr), the subject of the next PATHOS works, will be essential to shed light on the mechanisms of planet formation and evolution of close-in sub-Jovian planets.

SUPPORTING INFORMATION

Supplementary data are available at https://archive.stsci.edu/hlsp/pathos online.

Please note: Oxford University Press is not responsible for the content or functionality of any supporting materials supplied by the authors. Any queries (other than missing material) should be directed to the corresponding author for the article.

ACKNOWLEDGEMENTS

DN acknowledges the support from the French Centre National d’Etudes Spatiales (CNES). GL acknowledges support by CARIPARO Foundation, according to the agreement CARIPARO-Università degli Studi di Padova (Pratica n. 2018/0098). LB acknowledges the funding support from Italian Space Agency (ASI) regulated by ‘Accordo ASI-INAF n. 2013-016-R.0 del 9 luglio 2013 e integrazione del 9 luglio 2015 CHEOPS Fasi A/B/C’. The authors thank the anonymous referee for carefully reading this paper and for the useful suggestions that improved this work. This paper includes data collected by the TESS mission. Funding for the TESS mission is provided by the NASA Explorer Program. 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. Some tasks of the data analysis have been carried out using vartools v 1.39 (Hartman & Bakos 2016) and tlspython routine (Hippke & Heller 2019).

DATA AVAILABILITY

The data underlying this paper are available in MAST at doi:10.17909/t9-es7m-vw14 and https://archive.stsci.edu/hlsp/pathos.

Footnotes

REFERENCES

APPENDIX A: LIGHT-CURVE MODELLING