A publishing partnership

The following article is Free article

The Gaia–Kepler Stellar Properties Catalog. II. Planet Radius Demographics as a Function of Stellar Mass and Age

, , , , and

Published 2020 August 12 © 2020. The American Astronomical Society. All rights reserved.
, , Citation Travis A. Berger et al 2020 AJ 160 108DOI 10.3847/1538-3881/aba18a

Download Article PDF
DownloadArticle ePub

You need an eReader or compatible software to experience the benefits of the ePub3 file format.

1538-3881/160/3/108

Abstract

Studies of exoplanet demographics require large samples and precise constraints on exoplanet host stars. Using the homogeneous Kepler stellar properties derived using the Gaia Data Release 2 by Berger et al., we recompute Kepler planet radii and incident fluxes and investigate their distributions with stellar mass and age. We measure the stellar mass dependence of the planet radius valley to be / = , consistent with the slope predicted by a planet mass dependence on stellar mass (0.24–0.35) and core-powered mass loss (0.33). We also find the first evidence of a stellar age dependence of the planet populations straddling the radius valley. Specifically, we determine that the fraction of super-Earths (1–1.8 ) to sub-Neptunes (1.8–3.5 ) increases from 0.61 ± 0.09 at young ages (<1 Gyr) to 1.00 ± 0.10 at old ages (>1 Gyr), consistent with the prediction by core-powered mass loss that the mechanism shaping the radius valley operates over Gyr timescales. Additionally, we find a tentative decrease in the radii of relatively cool (Fp < 150 ) sub-Neptunes over Gyr timescales, which suggests that these planets may possess H/He envelopes instead of higher mean molecular weight atmospheres. We confirm the existence of planets within the hot sub-Neptunian "desert" (2.2 R < Rp < 3.8 , Fp > 650 ) and show that these planets are preferentially orbiting more evolved stars compared to other planets at similar incident fluxes. In addition, we identify candidates for cool (Fp < 20 ) inflated Jupiters, present a revised list of habitable zone candidates, and find that the ages of single and multiple transiting planet systems are statistically indistinguishable.

Export citation and abstractBibTeXRIS

1. Introduction

One of the most impactful exoplanet discoveries in recent years has been the planet radius "valley," a dip in the occurrence of Kepler planets at ≈1.9 separating super-Earth- and sub-Neptune-sized exoplanets (Owen & Wu 2013; Fulton et al. 2017; Fulton & Petigura 2018). The discovery of the radius valley was enabled by precise stellar parameters for subsamples of Kepler host stars, such as those derived in the California-Kepler Survey (CKS; Johnson et al. 2017; Petigura et al. 2017) and from asteroseismic constraints (van Eylen et al. 2018). More recently, Gaia parallaxes (Gaia Collaboration et al. 2018; Lindegren et al. 2018) have better constrained the stellar radii of the vast majority of Kepler host stars, followed by more detailed investigations of the valley as a function of stellar mass (Fulton & Petigura 2018), metallicity (Owen & Murray-Clay 2018), planet orbital period, and stellar incident flux (Berger et al. 2018b; Fulton & Petigura 2018). Most recently, the radius valley has also been identified in the K2 sample (Cloutier & Menou 2020; Hardegree-Ullman et al. 2020).

Several models have been proposed to explain the planet radius valley, including planet formation in a gas-poor disk (Lee et al. 2014; Lee & Chiang 2016), extreme ultraviolet (EUV) photoevaporation of planet atmospheres (Owen & Wu 2013, 2017; Lopez & Rice 2018; Owen & Murray-Clay 2018; Wu 2019), and core-powered mass loss (Ginzburg et al. 2016, 2018; Gupta & Schlichting 2019, 2020). Currently, photoevaporation and core-powered mass loss are the two leading theories that can effectively explain the dependence of the gap on stellar mass, orbital period, and incident flux. However, observational studies have not yet been able to differentiate between these two theories. For example, to explain the radius valley as a function of stellar mass, photoevaporation requires a planet mass dependence on stellar mass (Wu 2019). Multitransiting systems hosting planets with mass measurements and radii both above and below the gap can distinguish between these models, but only a few examples are available (e.g., Cloutier et al. 2020; Nowak et al. 2020).

Stellar ages provide a new dimension to determine the physical mechanisms shaping exoplanet populations. Exoplanet properties are expected to change over time, such as a decrease in their radii from cooling and contraction (Lopez et al. 2012) and atmosphere loss (Ginzburg et al. 2016; Owen & Wu 2017) or an increase in orbital eccentricity due to dynamical interactions between planets (Weiss et al. 2018). However, ages are difficult to determine for stellar populations, such as Kepler host stars, because available methods differ considerably across the Hertzsprung–Russell (H-R) diagram. For instance, isochrone ages are effective on the upper main sequence (M ≳ 1 ) but uninformative on the lower main sequence, where stellar rotation, activity, and lithium abundances provide more discriminatory power (Pont & Eyer 2004; Epstein & Pinsonneault 2014). Asteroseismology provides precise stellar ages but is generally only available for a small subset of mostly evolved exoplanet host stars (Silva Aguirre et al. 2015).

So far, only a few studies have compared properties of exoplanets orbiting stars of different ages. Berger et al. (2018a) found tentative evidence for the shrinking of planetary radii with stellar age based on the lithium abundances of CKS planet hosts differentiated by the Hyades 650 Myr empirical lithium abundance (A(Li))– curve (Boesgaard et al. 2016). Additionally, the Zodiacal Exoplanets In Time survey yielded evidence for larger, younger planets in clusters where ages are already known (Mann et al. 2018). While both core-powered mass loss and photoevaporation predict that the ratio of super-Earths to sub-Neptunes should increase over time, their timescales are very different: photoevaporation acts on timescales of ∼100 Myr (Lopez et al. 2012; Owen & Wu 2017), and core-powered mass loss acts on timescales of ∼Gyr (Gupta & Schlichting 2019, 2020). Stellar ages are also critical to address other open questions in exoplanet radius demographics, such as the hot sub-Neptunian desert (Lundkvist et al. 2016; Berger et al. 2018b; Dong et al. 2018), hot Jupiter inflation (Guillot & Showman 2002; Baraffe et al. 2010, 2014; Fortney & Nettelmann 2010; Laughlin & Lissauer 2015; Laughlin 2018; Komacek et al. 2020), and the dynamical evolution of multiplanet systems (Armitage & Rice 2005; Spalding & Batygin 2016; Rizzuto 2017; Weiss et al. 2018).

Previous Kepler stellar properties catalogs have not estimated ages due to inhomogeneous input parameters and the lack of precise parallaxes (Huber et al. 2014; Mathur et al. 2017). Here we rederive and analyze planet parameters using the updated stellar parameters by Berger et al. (2020, hereafter B20), the first homogeneous catalog of stellar , , radii, masses, densities, luminosities, and ages of Kepler stars.

2. Sample Selection and Methodology

2.1. Host Star and Planet Sample

First, we downloaded the KOI table on 2019 October 13 from the NASA Exoplanet Archive, including 9564 planet candidates. Then, we cross-matched this table with our Table 2 in B20, leaving 8875 planets. To avoid using stars with likely binary companions, we eliminated all stars with Gaia DR2 renormalized unit-weight error (RUWE) >1.2 (Evans 2018; Rizzuto et al. 2018; Bryson et al. 2020; B20; and see A. Kraus et al. 2020, in preparation). In addition, we discarded stars with unreliable isochrone-derived parameters (iso_gof < 0.99; B20). We also removed all planets designated as false positives according to the koi_disposition flag and those without reported planet-to-star radius ratios. We did not remove adaptive optics (AO)–detected binaries (Furlan et al. 2017) to preserve number statistics, but we comment on their influence where relevant. Following these sample cuts, we retained 2956 stars hosting 3898 planets.

2.2. Updated Planet Parameters

We computed the updated planet radii utilizing the planet-to-star radius ratios provided in the KOI table from the NASA Exoplanet Archive and the stellar radii computed in B20. In addition, we updated semimajor axes using the stellar masses in B20 and the orbital periods in Thompson et al. (2018). Finally, we updated the incident fluxes for each planet by using the semimajor axes and stellar luminosities from B20. We tested the effect of our new stellar parameters on the planet-to-star radius ratios by computing new planet-to-star radius ratios using the transit period, duration, and depth values from Thompson et al. (2018), the quadratic limb-darkening coefficients from Claret & Bloemen (2011), Equation (9) from Seager & Mallén-Ornelas (2003), and the small-planet-limit version of Equation (8) from Mandel & Agol (2002). We found that the differences between the planet-to-star radius ratios derived from Mathur et al. (2017) and B20 were on the order of 3% and within the 8% median uncertainty in Thompson et al. (2018). Therefore, these systematic effects are small, and we neglect them here. We provide all of our planet parameters in Table 1.

Table 1.  Planet Parameters

KIC ID KOI ID Planet Disposition Planet Radius () Semimajor Axis (au) Incident Flux () ZAMS Flux () Interesting Object Flag
11446443 1.01 Confirmed 524.75 AO
10666592 2.01 Confirmed 2394.68 AO
10748390 3.01 Confirmed 58.22  
3861595 4.01 Confirmed 3647.63 AO
11853905 7.01 Confirmed 719.80  
6922244 10.01 Confirmed 965.36 AO
5812701 12.01 Confirmed 141.77 YoungAO
10874614 17.01 Confirmed 450.14  
8191672 18.01 Confirmed 1079.02 AO
11804465 20.01 Confirmed 425.39  
9631995 22.01 Confirmed 181.79  
6521045 41.01 Confirmed 91.63 AO
6521045 41.02 Confirmed 209.73 AO
6521045 41.03 Confirmed 23.70 AO
8866102 42.01 Confirmed 91.23 AO
10905239 46.01 Confirmed 385.33  

Note. The interesting object flag denotes whether the host star is included in the older ("Old") or younger ("Young") than 1 Gyr samples selected in Section 4.2, the planets are located within the valley ("Gap"; see Section 4.3), the planet is located within the hot sub-Neptunian desert ("Desert"; see Section 5), planets are in the habitable zone ("HZ"; see Section 6.2), and/or the host has an AO-detected companion (Furlan et al. 2017). A subset of our planet parameters is provided here to illustrate the form and format.

Only a portion of this table is shown here to demonstrate its form and content. A machine-readable version of the full table is available.

Download table as:  DataTypeset image

3. Kepler Planet Host Stars

3.1. Histograms

We plot histograms of the physical parameters of the host star sample in Figure 1. Figure 1(a) shows the distribution of masses, which peak at solar mass, with a lower 1σ bound of 0.78 and an upper 1σ bound of 1.19 . The lowest- and highest-mass hosts are 0.14 and 2.8 , respectively.

Figure 1. Refer to the following caption and surrounding text.

Figure 1. Histograms of host star properties. The black dashed vertical lines illustrate the median value for each parameter. In panel (b), we plot both the overall host sample (black) and those with reliable ages (red).

Standard image High-resolution image

Figure 1(b) displays the distributions of stellar ages, separated into the entire host sample (black) and an "informative" age sample (red). The informative age distribution ignores any hosts that have a terminal age of the main sequence (TAMS) >20 Gyr, which are stars that evolve too slowly in their main-sequence lifetimes for isochrone fitting to constrain their ages (B20). The median ages are close to the age of the Sun, as expected from the Kepler target selection (Batalha et al. 2010). Both distributions peak at 3 Gyr with a tail toward old ages. Ages older than 14 Gyr occur because the B20 model grid of MESA (Paxton et al. 2011, 2013, 2015) Isochrones and Stellar Tracks (MIST v1.2 with rotation; Choi et al. 2016; Dotter 2016) purposefully used an upper limit of 20 Gyr to minimize grid edge effects, which can bias the parameter posteriors. The two distributions in Figure 1(b) differ the most at ≈10 Gyr, corresponding to M dwarfs with uninformative ages producing flat posteriors with medians at half the age of the grid. The small number of hosts older than 14 Gyr are stars whose input parameters place them on or close to the edge of the grid. Some of these old hosts are probably cool main-sequence binaries similar to those identified in Berger et al. (2018b), as stellar models are unable to reproduce their cool and large radii at the age of the universe. Finally, we note that the overall age distribution is consistent with asteroseismic ages provided by the APOKASC2 catalog of Kepler red giants (Pinsonneault et al. 2018).

Figures 1(c) and (d) show the stellar surface gravity () and mean stellar density distributions of Kepler host stars. These distributions are significantly different from the entire Kepler stellar sample, which has another, smaller peak at ≈ 2.5 dex and ρ ≈ 10−3 ρ. These peaks do not appear here because the percentage of giants with detected planets is much lower, given observational biases. Our median and ρ are slightly smaller than solar, and the tails to smaller values are comprised of subgiants (Verner et al. 2011; Everett et al. 2013; Gaidos & Mann 2013; Huber et al. 2014). There are a few giant hosts at the lowest and density values, but many of these stars host unconfirmed, potential false-positive planets (Sliski & Kipping 2014).

3.2. The Host Star H-R Diagram

Figure 2 shows stellar radii versus temperature for the Kepler host star sample, color-coded by age and maximum absolute age uncertainty. For the hottest stars, the age values and uncertainties are less than 2 Gyr. Both ages and age uncertainties increase smoothly toward cooler effective temperatures. Subgiants, compared to their main-sequence counterparts, have smaller error bars due to rapid evolution on the subgiant branch, while their ages can vary significantly based on their . On the giant branch, we also see that stars with the youngest ages and smallest age uncertainties are also the hottest and most massive, as more massive stars evolve more quickly than their lower-mass counterparts. We note that giant ages are potentially unreliable because of the limitations of isochrone fitting with photometric colors and a solar neighborhood metallicity prior for stars that do not have spectroscopic measurements (see B20 for details).

Figure 2. Refer to the following caption and surrounding text.

Figure 2. An H-R diagram of Kepler planet host stars, colored by their isochrone age (top; colors capped at 12 Gyr) and maximum absolute age uncertainties (bottom; colors capped at 6 Gyr). The gray points have uninformative ages (TAMS > 20 Gyr) and/or low goodness-of-fit values. Nine stars hotter than 8000 K are omitted from this plot.

Standard image High-resolution image

The ages reach a maximum of >12 Gyr at the main-sequence turnoff for the least massive stars, while the maximum age uncertainties increase toward the lower main sequence until they reach >6 Gyr. Zero-age main-sequence (ZAMS) F stars have low age uncertainties between zero and 2 Gyr, while ZAMS G dwarfs have moderate age uncertainties between 1 and 5 Gyr. For K dwarfs, the ZAMS uncertainties are typically larger than 6 Gyr. Some cool dwarfs with large/small radii have underestimated uncertainties due to grid edge effects. Finally, we see that all of the late K–M dwarfs have uninformative ages (TAMS > 20 Gyr; B20). In particular, their observables provide limited information with which we can distinguish between the ages of these stars, which evolve slowly in the H-R diagram over 14 Gyr.

Figure 2 illustrates variations in the effectiveness of isochrone placement for different types of stars, from subgiants (very effective) to K and M dwarfs (not effective at all). In addition, it demonstrates that the ages determined here can be used for the majority of Kepler planet host stars, and hence Kepler exoplanets, enabling one of the first investigations of how exoplanet properties change with stellar age. In the following investigations of stellar age, we ignore all gray points in Figure 2.

4. The Planet Radius Valley

4.1. Dependence on Stellar Mass

4.1.1. Results

We will now use our revised planet parameters to address whether features in the planet radius distribution are dependent on host star properties, particularly stellar mass. Figure 3 shows the distribution of planet radii as a function of incident flux, with individual planets colored according to host star mass. Planets at higher incident fluxes tend to orbit more massive hosts because these main-sequence stars have higher luminosities than their less massive counterparts. Additionally, we observe a separation between super-Earths (Rp  ≲ 2 ) and sub-Neptunes (Rp ≳ 2 ) in a narrow region from (≈10 , ≈1.5 ) to (≈1000 , ≈2.2 ). This is the planet radius valley, a 2D manifestation of the planet radius gap. The center of the gap is the planet radius at which the number density of planets is at a local minimum. Consistent with previous results (van Eylen et al. 2018; Martinez et al. 2019), we observe that the gap location increases in planet radius with increasing incident flux. Therefore, we also expect to find a dependence on stellar mass due to the stellar mass–luminosity relation.

Figure 3. Refer to the following caption and surrounding text.

Figure 3. Planet radius vs. incident flux for Kepler exoplanets. Points are colored according to the host star mass, as indicated by the color bar on the right. Planet candidates are shown as translucent points. The dashed box shows the sub-Neptunian desert identified in Lundkvist et al. (2016), and the green bar indicates the approximate "optimistic" habitable zone defined by Kane et al. (2016).

Standard image High-resolution image

Figure 4 displays the planet radius distributions in five equally populated stellar mass bins ranging from 0.14 to 2.8 . We see a clear dependence of the planet radius gap's location in the planet radius with stellar mass, increasing from 1.67 at M < 0.81 to 2.05 at M > 1.18 . To quantify this dependence, we computed the gap location in planet radius and its uncertainty by (1) drawing planet radii from normal distributions centered on the expected values and with standard deviations equal to their uncertainties and (2) computing the location of the gap by finding the relative minimum between 1 and 4 in the kernel density estimate (KDE) distribution. If no gap was found, we repeated steps (1) and (2). We then computed the standard deviation for 100 successful gap samples to determine the typical uncertainties in the gap location.

Figure 4. Refer to the following caption and surrounding text.

Figure 4. Distribution of Kepler exoplanet radii, binned by stellar mass. Each panel is labeled by stellar mass and includes 767 planets. The teal and green histograms represent confirmed and candidate planets, respectively, while the purple histogram represents all planets. The solid purple lines show the 0.12 bandwidth KDE of the combined planet population for that panel. The vertical dashed purple lines and the shaded regions show the gap locations and their uncertainties from our Monte Carlo simulations.

Standard image High-resolution image

Figure 5 directly compares planet radii with stellar mass. To quantify the stellar mass dependence of the planet radius gap's location in the planet radius, we used gapfit (Loyd et al. 2020) to fit a line of the form

where the average gap depth is deepest in the 2D KDE distribution (red line). To ensure gapfit found the line corresponding to the deepest valley, we constrained the range of the gapfit search to 1.5–2.4 .

Figure 5. Refer to the following caption and surrounding text.

Figure 5. Planet radius vs. stellar mass for Kepler exoplanets. The contours represent the 2D KDE distribution of the individual planets (small circles). Higher planet densities are darker colors. The blue shaded region illustrates the 1σ bounds of our Monte Carlo and bootstrap simulations using gapfit (Loyd et al. 2020). The red line represents our best fit to the data, with a slope of / = . The black lines illustrate the slope assuming a relation between the planet mass and stellar mass of Mp ∝ M (0.24; dashed line) and (0.35, dotted line; Wu 2019).

Standard image High-resolution image

To determine the uncertainty of the best-fit line while accounting for finite sampling, we drew each of the planet radii and stellar masses from normal distributions centered at each planet radius/stellar mass value with standard deviations given by the maximum of the upper and lower uncertainties for each planet radius and stellar mass. Next, we computed a 2D KDE from these simulated observations within the bounds of Figure 5 and constructed a sample of artificial planets drawn from the 2D KDE distribution. Then, we ran gapfit with M0 = 0.95 , R0,init = 1.89 , minit = 0.27, and sig = 0.15 (bandwidth in units of for the y-axis and for the x-axis with no covariance) with 100 bootstraps to determine the best-fit slope for the valley in the newly sampled distribution of planet radii and stellar mass. We repeated this process for 100 redraws of the planet radii and stellar masses and then computed uncertainties from these 100 × 100 determinations of the best-fit line using the uncertainty_from_boots routine provided within gapfit. The blue shaded region in the slope in Figure 5 represents the 1σ uncertainty region of the best-fit line (/ ).

Because core-powered mass loss predicts that the positive slope of the radius valley is due to the stellar mass–luminosity relation (Gupta & Schlichting 2020), we also compared the slopes of the valley across different stellar mass ranges. The slope of the mass–luminosity relation in log–log space, α (Eker et al. 2018), can be computed directly from the masses and luminosities of the host stars. We used gapfit with different but appropriate M0 and to measure the slopes of the radius valley for stars of α ≈ 2.9 (M = 0.14–0.72 ) and α ≈ 5.5 (M = 0.81–1.05 ). We also compared the slopes of the radius valley for stellar mass ranges including M < 0.81 , M > 1.18 , 0.8 MM <1.0 , and 1.0 MM < 1.2 . However, we found that all of these 2D KDE distributions and their radius valley best-fit lines were very sensitive to the choice of the KDE bandwidth and initial parameters; hence, we do not report their results here.

4.1.2. Discussion

Our results confirm the stellar mass dependence of the radius gap (Fulton & Petigura 2018; Cloutier & Menou 2020) with a best-fitting slope of / = . This result is consistent with the range of slopes predicted by a dependence of planet mass on stellar mass (which is required by photoevaporation models to explain the stellar mass dependence of the radius gap; Wu 2019) and core-powered mass-loss models (≈0.33; Gupta & Schlichting 2020). Therefore, we are unable to differentiate between these scenarios using the valley's slope in the planet radius–stellar mass diagram (see also Loyd et al. 2020).

To evaluate whether our inability to differentiate between core-powered mass loss and a planet mass dependence on stellar mass in planet radius–stellar mass space is simply a problem of sample size that might be ameliorated by future discoveries, we ran Monte Carlo simulations similar to those used in Figure 5 to test various sample sizes and measurement precisions. We found that a sample size of 20,000 planets (≈4 times the current number) assuming typical planet radius and stellar mass errors of 1% is needed to reduce the uncertainties in / to 0.04. This uncertainty would allow a 2σ–3σ separation between the lower and upper bounds of a planet mass dependence on stellar mass (0.24–0.35; Wu 2019) that contain our measured value (0.26) and the value predicted by core-powered mass loss (0.33; Gupta & Schlichting 2020). Even after combining all planets discovered by Kepler, K2 (Howell et al. 2014), and those already and yet to be discovered by the Transiting Exoplanet Survey Satellite (TESS; Ricker et al. 2014), both the sample size of ∼10,000 planets and the precision of our planet and stellar parameters would be insufficient to differentiate between photoevaporation, which requires a planet–stellar mass dependence to explain the stellar mass dependence of the radius gap, and core-powered mass loss in planet radius–stellar mass space.

4.2. Dependence on Stellar Age

4.2.1. Sample Selection

The photoevaporation and core-powered mass-loss mechanisms operate on different timescales and consequently predict different dependences of the planet radius distribution with age. Thus, we grouped planets according to the age of their host star. To ensure we only included planets with reliable properties, we removed 64 planets with radii greater than 30 . Additionally, we removed 97 planets with grazing transits defined as + b > 1, where b, the impact parameter, is how far away from the center of the star's disk the center of the planet transits at mid-transit (zero being at the center and 1 being at the limb). We used (similar to Equation (1) in Fulton et al. 2017) to remove giant hosts, which have a high false-positive rate (Sliski & Kipping 2014). We also removed 1081 planets orbiting hosts without spectroscopic metallicity constraints because isochrone age estimates are metallicity-sensitive (Howes et al. 2019). After this metallicity cut, 2545 planets remained.

We selected 1 Gyr as our separator between young and old systems because core-powered mass loss acts on ∼Gyr timescales (Gupta & Schlichting 2020). We could not make an age cut at 100 Myr, the timescale relevant for photoevaporation (Owen & Wu 2017), because we have no stars with isochrone ages <100 Myr. Likewise, a cut at the age of the Hyades (≈650 Myr; Boesgaard et al. 2016) includes too few young planets for a robust statistical comparison. We removed seven young planets with host effective temperatures hotter than 7900 K because their and masses were outliers compared to the rest of the age sample. Eighty-five planets orbit hosts with median posterior ages younger than 1 Gyr, while the remaining 2453 planets orbit hosts with median posterior ages older than 1 Gyr.

To remove degeneracies between stellar mass, age, and metallicity in the old and young planet samples, we used the NearestNeighbors function in scikit-learn (Pedregosa et al. 2011) to choose the two nearest old neighbors for every young host in stellar mass, radius, and metallicity. Because our matching function occasionally chose old neighbors such that multiple young stars are matched to the same old star, we removed any duplicate old hosts to avoid counting the same old host twice. We chose to use the two nearest neighbors instead of either one or three because the former selected fewer old planets than young planets after dropping duplicates, while the latter produced inferior stellar property-matched samples, especially in stellar mass. Our resulting property-matched sample included 90 old planets. In addition to removing stellar population biases, this careful sample selection also reduced potential detection biases for small planets between the old and young samples, all while retaining the full young planet sample. We used Kolmogorov–Smirnov (K-S) tests to compare the stellar mass, radius, and metallicity distributions of the old and young samples, producing p-values of 0.24, 0.58, and 0.997, respectively, confirming that the distributions are statistically similar. Figure 6 shows an H-R diagram of the host star sample after making the above cuts, separated by the relevant age bins. We flag these planets as "Old" and "Young" in Table 1.

Figure 6. Refer to the following caption and surrounding text.

Figure 6. An H-R diagram showing host stars colored according to their ages, as well as the marginalized distributions of stellar radii, with shaded areas representing the 16th–84th percentile ranges. Purple host stars have ages greater than 1 Gyr, and green host stars have ages younger than 1 Gyr. The gray points are host stars that we do not include in our old and young samples: all evolved stars (above the teal line), stars hotter than 7900 K (left of the dashed green line), and old stars with dissimilar mass, radius, and/or metallicity to the young stellar sample.

Standard image High-resolution image

4.2.2. Results

Figure 7 shows the planet radius distributions of the old and young planet samples. We observe that the gap occurs at approximately the same planet radius in both the old and young planet distributions. Remarkably, we also observe that the ratio of super-Earths compared to sub-Neptunes significantly increases from young to old stellar ages.

Figure 7. Refer to the following caption and surrounding text.

Figure 7. Radius distributions in KDE and histogram of Kepler exoplanets with ages younger (green) or older (purple) than 1 Gyr. The individual planet radii are plotted at the bottom as vertical ticks. The black dashed vertical line at 1.8 is the Fulton et al. (2017) gap radius, which separates super-Earths and sub-Neptunes. The ratio of super-Earths (1–1.8 ) to sub-Neptunes (1.8–3.5 ) significantly increases from young (0.61 ± 0.09) to old (1.00 ± 0.10) ages.

Standard image High-resolution image

To quantify this age dependence, we computed the uncertainties in the ratio of super-Earths to sub-Neptunes using Monte Carlo simulations to draw each old and young planet radius from a normal distribution given its measured value and uncertainty. We define super-Earths as planets between 1 and 1.8 and sub-Neptunes as planets between 1.8 and 3.5 (Fulton et al. 2017). We then counted the number of super-Earths and sub-Neptunes, repeated this process 1000 times for both the old and young distributions, and computed the standard deviation of the ratios of super-Earths to sub-Neptunes.

Our data show a significant increase in the fraction of super-Earths as a function of age, with the fraction of super-Earths to sub-Neptunes increasing from 0.61 ± 0.09 at young ages (<1 Gyr) to 1.00 ± 0.10 at old ages (>1 Gyr). This result is insensitive to the choice of impact parameter cut (≳2σ for b < 0.7–0.9), gap location (≳3σ for 1.9–2.0 ), and radius range used to define super-Earths and sub-Neptunes (≳4σ for 0.8–1.8 and 1.8–5 , respectively). Similarly, if we instead used 1533 old planets hosted by old stars larger than 0.9 rather than the property-matched sample described in the previous section, we computed a ≳3σ difference in the ratios of young and old super-Earths to sub-Neptunes. Reassuringly, Fulton et al. (2017) computes the occurrence ratio of super-Earths to sub-Neptunes for the entire planet sample to be 0.8 ± 0.2, approximately the average of our old and young ratios.

We also compared the low (Fp < 150 ) and high (Fp > 150 ) flux planet radius distributions for old and young exoplanets (Figure 8). We chose 150 because it splits the young sample of planets almost in half: 41 young and 41 old planets receive more than 150 , and 44 young and 49 old planets receive less than 150 . Interestingly, we observe stark differences between old and young planets. At high incident flux (Figure 8(a)), we observe a large difference between the young and old planet radius distributions. We compute the ratios of super-Earths to sub-Neptunes to be 1.00 and 2.67 for the young and old planets, respectively. However, at low incident flux (Figure 8(b)), the overall distributions do not show a strong difference as a function of age, with tentative evidence that old sub-Neptunes are smaller than young sub-Neptunes. We ran K-S tests to quantitatively compare the old and young distributions in both panels and found that they were statistically distinguishable in the top panel (p-value = 0.02) and indistinguishable in the bottom panel (p-value = 0.11).

Figure 8. Refer to the following caption and surrounding text.

Figure 8. Same as Figure 7 but separating the samples into high (>150 ) and low (<150 ) incident flux.

Standard image High-resolution image

A potential bias in our results is the sensitivity of planet detection to stellar age. This is because young stars are typically noisier due to their increased activity (Skumanich 1972); thus, smaller planets might not be detected around them. To test this, we evaluated the CDPP3 values (Christiansen et al. 2012) to determine the single-transit signal-to-noise ratio (S/N) for an Earth-sized planet with a 3 hr transit duration (as in Petigura et al. 2018) for each host star in Figure 9. This comparison also accounts for differences in the stellar radius distributions, although we control for these differences in our sample selection.

Figure 9. Refer to the following caption and surrounding text.

Figure 9. Single-transit S/N for an Earth-sized planet with a 3 hr transit duration orbiting each individual planet host used in Figures 68 (large circles) and Kepler target star (dots). Larger values on this plot indicate that planets are easier to detect. We use the CDPP3 values from Christiansen et al. (2012) and Equation (B2) in Petigura et al. (2018) to compute the S/N.

Standard image High-resolution image

We find that the young planets typically have similar S/Ns compared to the old planets and hence are just as difficult to detect. The large scatter of points in Figure 9 dominates what would be small differences in the old and young median single-transit S/N values. Similarly, when we compare the median single-transit S/N values, the overall Kepler sample with cuts similar to those in Figure 6, we find no significant difference between the old (≈0.43) and young (≈0.40) median S/N, given the large point-by-point scatter. Therefore, we conclude that our results are robust against a planet detection bias with stellar age.

Another bias that could affect our results is the presence of undetected binary companions, which would cause stars to be incorrectly assigned old ages. While we have removed most wide binaries using Gaia and debiased stellar photometry for binary contamination (B20), it is inevitable that some close binaries will contaminate our old and young samples. After removing all AO-detected companions (Furlan et al. 2017), the young and old distributions in Figure 7 remain mostly unchanged. However, the distributions in Figure 8(a) are sensitive to the removal of AO-detected binaries, as binaries make stars appear more luminous and increase the apparent irradiance of the planets.

To test the potential for binaries to bias our results, we estimated the number of planets that may have been mistakenly placed in the older sample because of binary contamination. We adopted a typical binary mass ratio ( ≈ 0.5) and binary fraction for solar-type stars (Fbin ≈ 40%) determined in Raghavan et al. (2010) and Moe & di Stefano (2017). We also used the +50% isochrone age bias determined from the 1.15 star in Figure 4 of B20. Using these assumptions, we predict that 55 planets × Fbin ≈ 22 planets are orbiting stars with isochrone ages between 1 and 1.5 Gyr that should be younger than 1 Gyr due to binary contamination. We then shifted these 22 planets into the young distribution, assuming that (1) they mimic the old distribution's ratio of super-Earths to sub-Neptunes and (2) their measured radii are not affected significantly by binary contamination. Consequently, we computed a 0.76 ± 0.09 ratio of super-Earths to sub-Neptunes for the young planet distribution. Comparing this to the 1.00 ± 0.10 ratio of super-Earths to sub-Neptunes for the old planet distribution, we still arrive at a >2σ difference between the old and young distributions. Given our conservative assumptions, we conclude that undetected stellar companions will not significantly affect our results on the age dependence of the radius valley.

4.2.3. Discussion

Figure 7 suggests that sub-Neptunes evolve to become super-Earths over Gyr timescales. This is consistent with the core-powered mass-loss mechanism (Gupta & Schlichting 2020), which predicts that the transition of sub-Neptunes to super-Earths occurs over Gyr timescales as the planets gradually lose their atmospheric envelopes due to their heated cores. While our result supports the core-powered mass-loss mechanism, it does not rule out the possibility that photoevaporation is also acting on the observed planet population. Given that photoevaporation is expected to occur in the first ∼100 Myr (Owen & Wu 2017; Wu 2019) of a planet's lifetime and our 100 Myr minimum grid age (B20), we are mostly insensitive to any planetary evolution before 100 Myr ages, barring extreme planet radius evolution during that time. Our Monte Carlo simulations do not account for the uncertainties in the stellar ages; hence, it is likely that there is contamination between the old and young planets. However, given the magnitude of the errors and the difference in the number of old and young stars in our sample, it is more likely that older stars are contaminating the younger bins (Eddington bias). If this bias truly exists in our sample, any observed differences between the old and young samples will be reduced by this bias; thus, the true difference between young and old planets is greater.

Our result agrees with Berger et al. (2018a), who used lithium abundances relative to the Hyades to separate young (A(Li) > A(Li)Hyades) and old (A(Li) < A(Li)Hyades) planets. In particular, Berger et al. (2018a) found that there is a significant difference in the sizes of planets in the old and young planet radius histograms and that young planets are larger. This result was tentative (≈2σ–3σ) due to the small sample size of the young planets as compared to the old planets. Although Berger et al. (2018a) did not explicitly quantify the number of super-Earths and sub-Neptunes, the largest difference between their old and young samples occurs at ≳2 , where the number of sub-Neptunes in the young sample is significantly greater than the number of sub-Neptunes in the old sample.

Our observation of a similar gap radius for young and old planets is expected, as we carefully chose stellar samples with similar stellar mass, which should produce a gap at a similar location in the distribution of planet radii (Section 4.1). While core-powered mass loss predicts that the gap's location should increase slightly to larger radii in the first ∼Gyr (Gupta & Schlichting 2020), such behavior will be hard to definitively establish given typical uncertainties on stellar age, mass, and planet radius. Meanwhile, photoevaporation does not predict significant gap movement after the first ∼100 Myr (Owen & Wu 2016, 2017).

We interpret that the differences between the old and young distributions as a function of incident flux in Figure 8 are consistent with core-powered mass loss, which predicts that planets that receive higher incident fluxes may experience increased and potential runaway mass loss over timescales of ∼Gyr due to increased equilibrium temperatures at the Bondi radius (Ginzburg et al. 2016, 2018; Gupta & Schlichting 2019, 2020). Conversely, sub-Neptunes that receive low incident fluxes (Figure 8(b)) will simply cool and contract, shifting from larger radii (≈3 ) at young ages to smaller radii (≈2.5 ) at old ages. This is predicted by photoevaporation as well (Lopez et al. 2012; Owen & Wu 2017), but it cannot describe the difference in high incident flux distributions over Gyr timescales. In addition, the marginally significant difference (p-value = 0.11) in the sizes of old and young sub-Neptunes at low incident fluxes may suggest that these planets have significant H/He envelopes instead of higher mean molecular weight envelopes that cannot produce the same magnitude of contraction (Nettelmann et al. 2011; Lopez et al. 2012; Lopez & Fortney 2014; Lopez 2017; but see also Howe & Burrows 2015). Alternatively, the shrinking of these planets could also be caused by H2/He ingassing at orbital periods <100 days (Kite et al. 2020). More theoretical and observational work is required to evaluate the composition of these atmospheres (Owen 2019), and we caution that these inferences are tentative, especially because of the small number of planets contained in the old and young low-flux (49 and 44 planets, respectively) and old and young high-flux (41 and 41 planets, respectively) distributions. Additionally, both the young and old low-flux planet radius distributions are likely biased more significantly by selection effects than the high-flux distributions, as small planets at low incident fluxes are more difficult to detect.

4.3. Planets within the Gap

According to photoevaporation, which produces radius changes within the first 100 Myr of a host's lifetime (Owen & Wu 2017), we should not see any planets within the gap (van Eylen et al. 2018) at old ages. Conversely, the ∼Gyr timescales of core-powered mass loss (Gupta & Schlichting 2020) suggest that we should find a few planets in the gap region as they transition from sub-Neptunes to super-Earths. It is currently unclear whether any planets firmly exist within the radius valley. Large population studies have revealed that there are at least a few planets that fall within the planet radius gap (Berger et al. 2018b; Fulton & Petigura 2018), although smaller, more precise planet samples have revealed a complete lack of planets within the gap (van Eylen et al. 2018).

To investigate this, we used gapfit (Loyd et al. 2020) to determine the best-fit parameters in the planet radius–incident flux diagram using

assuming a pivot point F0 = 100 and the optimal Gaussian kernel width sig = 0.15 (gapfit; Loyd et al. 2020). We found m = / and R0 = 1.86 . Next, we computed parallel lines by varying the R0 parameter by typical uncertainties (±0.09 ) determined from a combination of Monte Carlo and bootstrap simulations, ignoring uncertainties in the slope. We then isolated all confirmed planets that were within the log–log lines with slopes, m = /, and pivot point central gap radii, R0 = 1.77 and 1.95 . We also removed all planets with 1σ errors in planet radius that would place them outside the log–log lines representing bounds of the planet radius gap as a function of incident flux.

Following these cuts, five planets remain in the gap: Kepler-11 b, Kepler-110 b, Kepler-114 c, Kepler-634 b, and Kepler-887 b. We conclude that these five planets may currently undergo core-powered mass loss but caution that they are only ∼1σ removed from the gap boundaries. Additional follow-up observations will be required to definitively identify planets inside the radius valley.

5. The Hot Sub-Neptunian Desert

The hot sub-Neptunian desert, defined as Rp = 2.2–3.8 and Fp > 650 , is another region of parameter space believed to be devoid of planets (Lundkvist et al. 2016). The lack of planets in the hot sub-Neptunian desert can be explained by photoevaporation (Lopez et al. 2012; Owen & Wu 2016, 2017), with hydrogen and helium atmospheres being completely lost to the large EUV flux at small orbital separations (Lopez 2017). At lower incident EUV fluxes, planets still lose significant portions of their atmospheres but stabilize at an envelope mass fraction of 1%–2% (Owen & Wu 2017). Berger et al. (2018b) found 74 confirmed planets/planet candidates in the "desert" and suggested that some of these planets (1) may be the remnants of the photoevaporation of a giant planet's envelope (Baraffe et al. 2005), (2) did not receive enough EUV flux to lose their low molecular weight atmospheres (Owen & Wu 2017), or (3) may have high molecular weight atmospheres (Lopez 2017). Our newly derived stellar masses and ages can shed additional light on the properties and formation of this intriguing class of planets.

Figure 10(a) shows the age distribution of planets within (red) and outside (blue) the hot sub-Neptunian desert at high incident fluxes (>650 ; Lundkvist et al. 2016). We also plot the overall Kepler sample for comparison. We observe similar stellar age distributions for planets within the hot sub-Neptunian desert and those at high incident fluxes with radii outside 2.2–3.8 . The median ages and distributions are almost identical, which is also confirmed with a K-S test (p-value = 0.14). We therefore conclude that most of these "desert dwellers" are not young planets that are currently losing mass.

Figure 10. Refer to the following caption and surrounding text.

Figure 10. Stellar age (panel (a)) and mass (panel (b)) distributions for planets with high incident flux (Fp > 650 ) within the hot sub-Neptunian desert Rp = 2.2–3.8 (red), at high incident flux outside the hot sub-Neptunian desert (blue), and the overall Kepler target sample with reliable ages (RUWE < 1.2, iso_gof > 0.99, and TAMS < 20 Gyr; black). The blue and red ticks represent the individual ages/masses used to calculate the KDEs using bandwidths following Scott's rule (Scott 1992). The dashed vertical lines and shaded areas are the median and 1σ bounds.

Standard image High-resolution image

Unlike stellar ages, Figure 10(b) indicates a difference between the stellar mass distributions of high incident flux planets inside and outside the desert. A K-S test yielded a p-value of 0.02, indicating a difference at ≈2σ significance. Therefore, we tentatively conclude that desert planets tend to be around more massive stars. This conclusion is also supported by a K-S test using host star radii and uncertainties, yielding a p-value of 0.007. Because hot sub-Neptune planet hosts appear to have higher stellar masses and larger stellar radii, we hypothesize that these planets could have shifted into the desert through stellar evolution.

To investigate this, Figure 11(a) shows the incident flux history of planets within and outside the desert. We count 35 confirmed and 15 planet candidates within the hot sub-Neptunian desert (as defined by Lundkvist et al. 2016), and we denote them with the "Desert" flag in Table 1. Figure 11(b) shows the ratios of the current flux compared to the ZAMS flux (see Table 1 for planets in the different regimes of panel (a)). Together, Figures 11(a) and (b) suggest that 60% of desert planets have moved into the desert as a result of stellar luminosity evolution. Therefore, we infer that in the first ∼100 Myr, the majority of desert planets were not exposed to enough EUV flux to completely strip their atmospheres. However, this incident flux evolution does not explain all desert planets.

Figure 11. Refer to the following caption and surrounding text.

Figure 11. Panel (a): planet radius vs. incident flux for Kepler exoplanets in the hot sub-Neptunian desert (red) defined by Lundkvist et al. (2016), at incident fluxes >650 and outside the desert (blue), and at incident fluxes <650 (gray). Solid and translucent points are confirmed and candidate planets, respectively. The black bars show the incident flux history of each planet, starting at the incident flux the planet received at the ZAMS. Panel (b): current divided by ZAMS flux ratio for desert planets (red), other high incident flux planets (blue), and all other planets (gray).

Standard image High-resolution image

We find nine confirmed planets within the desert by >1σ even after accounting for the effects of stellar evolution on the incident flux: Kepler-234 b, Kepler-541 b, Kepler-611 b, Kepler-644 b, Kepler-645 b, Kepler-656 b, Kepler-1016 b, Kepler-1171 b, and Kepler-1518 b. If we assume a ZAMS luminosity that corresponds to the ZAMS luminosity of the lower-mass uncertainty bound for each system, we still find three planets within the desert: Kepler-644 b, Kepler-645 b, and Kepler-1171 b. Kepler-644 and Kepler-645 have AO-detected companions; thus, their stellar properties are likely inaccurate, in addition to the incident fluxes their planets receive. Using ZAMS luminosities corresponding to the lower-mass uncertainty bound is a pessimistic assumption that overestimates the increase in the star's luminosity since the ZAMS and accounts for uncertainties. Given 100 Myr H/He atmosphere loss timescales (Owen & Wu 2017), it is unlikely that these desert-dwelling planets are typical sub-Neptune-mass planets with H/He envelopes, unless they migrated to their current orbital separations (Dong et al. 2018) and/or have higher molecular weight atmospheres (Lopez 2017; Gaidos et al. 2020). It is also possible that these planets are the bare cores of 2–3× more massive planets (Armstrong et al. 2020). While photoevaporation is not expected to strip enough atmospheric mass off of these massive cores, tidal disruption (Vick et al. 2019) or giant planet collisions (Mordasini 2018) can produce the required mass loss. Alternatively, these bare cores may form in situ by opening up a gasless gap in the protoplanetary disk and avoiding runaway accretion (Lee 2019). Ultimately, these desert dwellers represent interesting tests of planet formation and evolution theories and warrant additional scrutiny.

6. Cool Planets

6.1. Cool Inflated Jupiters

The mechanism producing the inflated radii of hot Jupiters is still a major unsolved problem (Guillot & Showman 2002; Baraffe et al. 2010, 2014; Fortney & Nettelmann 2010; Fortney et al. 2011; Laughlin & Lissauer 2015; Laughlin 2018; Komacek et al. 2020). Most theories are linked to the observation that all inflated Jupiters experience high incident flux (Demory & Seager 2011; Laughlin et al. 2011; Miller & Fortney 2011), whether they are orbiting main-sequence stars or have been reinflated from the effects of post-main-sequence stellar evolution (Grunblatt et al. 2016, 2017, 2019; Lopez & Fortney 2016). Thus, finding examples of cool, inflated Jupiters may present an interesting challenge to these theories. Berger et al. (2018b) identified three Jupiters with anomalously large radii at <150 . If these planets are inflated at low incident fluxes, there must be some other mechanism causing their inflation. For example, these planets may be young and hot from the gravitational energy of accretion, and this additional energy could produce an inflated atmosphere (Lopez et al. 2012). These cool inflated Jupiters could also be heated from recent tidal interactions with other planets and their host star (Jackson et al. 2008; Fortney et al. 2010).

Similar to Berger et al. (2018b), we find a small sample of cool, confirmed (or validated), >1 RJ Jupiters (Figure 12). In addition to isochrone ages, we mark points according to whether they exhibit UV excess (red circles) or rapid rotation (orange squares) as additional indicators for youth. We flag a star as having UV excess when it meets two criteria: (1) (to avoid magnitude-limited cases) and (2) , where the condition is set by the Hyades relation evaluated at the of the Kepler star. The Hyades cluster (∼650 Myr; Boesgaard et al. 2016) has a well-defined trend in near-UV (NUV)–Ks versus , making it an effective separator for young/old stars in the Kepler field (Berger et al. 2018a). We obtained the NUV fluxes for the Kepler and Hyades stars from the Galaxy Evolution Explorer (GALEX; Martin et al. 2005; Olmedo et al. 2015). Similarly, we flag a star as having rapid rotation if its rotation period is more rapid than the Hyades gyrochrone with an initial rotation period of 3.4 days (Kundert et al. 2012). We use the rotation periods derived in McQuillan et al. (2013, 2014) or Mazeh et al. (2015).

Figure 12. Refer to the following caption and surrounding text.

Figure 12. Radius vs. incident flux for Kepler's cool Jupiters. Planets are colored by isochrone age. Orange squares and red circles indicate stars that are rotating more rapidly and exhibit more UV excess compared to the Hyades, respectively. Gray points are candidate planets. The red curve represents the maximum radius for a 4.5 Gyr, Jupiter-mass, pure hydrogen and helium object (Thorngren et al. 2016). The inset shows the position of the host stars on the H-R diagram compared to the Kepler target sample (gray).

Standard image High-resolution image

We find two planets that are significantly above the maximum radius for a 4.5 Gyr, Jupiter-mass, pure hydrogen and helium object (Thorngren et al. 2016). Kepler-468 b appears to be young, according to its host's isochrone age ( Gyr), excess UV flux (Skumanich 1972; Soderblom 2010), and rapid rotation period. It does not have a rotation period detection according to McQuillan et al. (2013, 2014) or Mazeh et al. (2015), but an inspection of the Kepler light curve reveals a rotation period of ≈5.7 days, consistent with a young age (Barnes & Kim 2010). Kepler-468 has an RUWE = 1.15, which is below the threshold for being a likely binary. Similarly, Kepler-468 does not appear to have any companions, according to high-resolution imaging (Law et al. 2014; Furlan et al. 2017). Even at 2.4 Gyr, a Jupiter-mass planet may still be cooling and contracting from its heat from formation, although the majority of this contraction occurs within the first ∼Gyr (Fortney et al. 2007; Linder et al. 2019). Hence, it is possible that Kepler-468 b is young and still cooling and contracting from its heat from formation.

Unlike Kepler-468 b, Kepler-706 b does not appear to orbit a young star. Kepler-706 has an isochrone age of ∼17 Gyr, a rotation period of ≈38 days measured from the Kepler light curve, and an NUV magnitude that is beyond the GALEX limiting magnitude of 22.6 (Olmedo et al. 2015). While the isochrone age is unreliable, the rotation period supports an old age for this star (Barnes & Kim 2010). Furthermore, neither the RUWE value (1.13) nor high-resolution imaging from Law et al. (2014) indicate the presence of a stellar companion.

If confirmed, the radius and age of Kepler-706 b would be highly interesting. Given its old age, it is very unlikely that it is inflated from its residual heat from formation based on the expected cooling timescales (Lopez et al. 2012). Kepler-706 b orbits every 41 days, so it is also unlikely that its radius is inflated by strong tidal or magnetic interactions with its host star, barring a highly elliptical orbit. Planet–planet interactions typically produce only a fraction of the tidal heating that is caused by the host star (Hay & Matsuyama 2019), but this heating is strongly dependent on the proximity and mass of the other planet. Rings are another potential explanation for Kepler-706 b's inflated radius (Schlichting & Chang 2011), but we cannot confirm the presence of rings with Kepler's light-curve cadence/precision. Additional follow-up observations will be required to confirm whether Kepler-706 b is indeed a bona fide cool inflated Jupiter.

6.2. Habitable Zone Planets

Our work also yields a revised list of planets in the habitable zone. Following the definition of Kane et al. (2016), we find 133 planet candidates and 32 confirmed planets within the habitable zone (0.25–1.5 ), all of which are flagged in Table 1. Compared to Berger et al. (2018b), we count two fewer confirmed and 24 more candidate planets. We report eight confirmed planets with radii <2 : Kepler-62 e, Kepler-62 f, Kepler-186 f, Kepler-283 c, Kepler-440 b, Kepler-442 b, Kepler-452 b (but see also Mullally et al. 2018), and Kepler-1544 b, although Kepler-62, Kepler-186, and Kepler-283 have AO-detected companions (Furlan et al. 2017). This increase of habitable zone candidate planets is mostly caused by the slightly cooler derived in B20, which produces smaller incident fluxes for planets at the same orbital periods.

7. Single- and Multiplanet Systems

Several observational results suggest that the Kepler multiplanet systems (multis) have smaller eccentricities and inclinations than the Kepler systems with just one small transiting planet (singles). The Kepler multis are nearly coplanar (Fang & Margot 2012; Fabrycky et al. 2014), with a typical mutual inclination of 1fdg0–2fdg2. These typical mutual inclinations are smaller than what is required for all of the planets to transit, indicating that detection bias alone cannot explain the observed coplanarity of the Kepler multis. Furthermore, the Kepler singles have higher eccentricities than the multis, as determined from their more varied transit durations (Xie et al. 2016; Mills et al. 2019; van Eylen et al. 2019). In multis with an ultrashort-period planet, the mutual inclinations between the planets are higher than in multis without an ultrashort-period planet (Dai et al. 2019), suggesting that the evolutionary pathway that creates ultrashort-period planets also excites planetary inclinations and eccentricities (e.g., Petrovich et al. 2019).

What physical processes affect the eccentricities and inclinations of the Kepler planets? The host star masses, metallicities, and values of the Kepler singles and multis are statistically indistinguishable (Weiss et al. 2018), suggesting that the host stars are not the main source of the observed differences in the eccentricity and inclination distributions. One possible source of dynamical excitation (and also instability) is the planets themselves: dynamically packed planetary systems become unstable on timescales of 107–1010 orbits (Obertas et al. 2017), with the closest-packed systems typically becoming unstable earliest. Systems with larger eccentricities are more likely to become unstable, excited in inclination, and tidally circularized (Petrovich et al. 2019; Pu & Lai 2019). Long-distance processes, such as a gravitational perturbation from a passing star, can also affect the coplanarity and stability of the Kepler multis (Spalding & Batygin 2016). The likelihood that any of these disruptive mechanisms has already occurred increases with time. Therefore, significant differences in stellar age between the single and multiplanet systems, in particular a preference for single planets to be around older stars, may point to an evolutionary pathway from the multis to the singles.

Figure 13 shows cumulative distribution functions (CDFs) for ages of multis (red; 539 systems) and singles (blue; 1768 systems). We observe that the single transiting planet systems indeed appear to be, on average, older than the multiple transiting planet systems. To quantify the difference of these two distributions, we conducted Monte Carlo K-S simulations to account for individual errors in the ages of our stars. We also performed Monte Carlo Anderson–Darling (A-D, scipy's anderson_ksamp; Scholz & Stephens 1987; Virtanen et al. 2020) test simulations in a similar manner. The resulting p-values are 0.12 and 0.10 for the K-S and A-D tests, respectively, indicating that the single- and multiplanet ages are statistically indistinguishable. Requiring spectroscopic metallicities for single and multi systems, which improve isochrone ages (Howes et al. 2019), produced p-values of 0.59 and 0.25 (capped) for the K-S and A-D tests, respectively.

Figure 13. Refer to the following caption and surrounding text.

Figure 13. The CDFs of Kepler systems with multiple transiting planets (red; 539 systems) and one transiting planet (blue; 1768 systems).

Standard image High-resolution image

We therefore conclude that the ages of single and multiple transiting planet systems are statistically indistinguishable within the precision of our age measurements. This implies that dynamical interactions that affect planet eccentricities and inclinations, such as planet–planet and planet–passing star interactions, operate significantly faster than the ∼Gyr timescales that can be resolved by our sample (e.g., Spalding & Batygin 2016; Obertas et al. 2017) and may happen at almost any age. For instance, young singles could be the multis that were disrupted at young ages, while old multis have not been disrupted yet.

8. Summary and Conclusions

We presented a reanalysis of Kepler exoplanet properties using the uniform stellar properties provided by B20. In particular, we performed a systematic analysis of planet properties as a function of stellar mass and age for the entire Kepler sample. Our main conclusions are as follows.

  • 1.  
    We observe that the location of the planet radius gap increases in planet radius with increasing stellar mass with a slope of / = , consistent with previous results based on smaller samples by Fulton & Petigura (2018) and Cloutier & Menou (2020). The uncertainty on the slope, even for the full Kepler sample, is too large to discern between a planet mass dependence on stellar mass (as required by photoevaporation; Wu 2019) and core-powered mass loss (Gupta & Schlichting 2020). We estimate that differentiating between these theories would require ≳20,000 planets with 1% fractional precisions in planet radius and stellar mass.
  • 2.  
    We find the first evidence for a stellar age dependence of the fraction of super-Earths to sub-Neptunes, increasing from young (<1 Gyr; 0.61 ± 0.09) to old (>1 Gyr; 1.00 ± 0.10) planets. This is consistent with predictions of the core-powered mass-loss mechanism that sub-Neptunes evolve to become super-Earths over Gyr timescales (Gupta & Schlichting 2020), but we caution that our sample does not constrain evolution at times ≲100 Myr due to photoevaporation (Owen & Wu 2017). We also observe that the age dependence of the fraction of super-Earths to sub-Neptunes becomes stronger with higher incident flux. Additionally, there is a marginally significant (p-value = 0.11) decrease in the radii of cool (Fp < 150 ) sub-Neptunes, which may suggest that these sub-Neptunes have H/He envelopes as opposed to higher mean molecular weight atmospheres. However, we caution that the latter results may be influenced by small number statistics and possible effects of undetected binary companions.
  • 3.  
    We identify 35 confirmed planets and 15 planet candidates that occupy the hot sub-Neptunian desert (Lundkvist et al. 2016). We determine that most of the desert planets orbit evolved stars, and thus it is unlikely that they are young planets currently undergoing mass loss. In addition, we identify nine planets that are within the desert even after accounting for stellar evolution: Kepler-234 b, Kepler-541 b, Kepler-611 b, Kepler-644 b, Kepler-645 b, Kepler-656 b, Kepler-1016 b, Kepler-1171 b, and Kepler-1518 b.
  • 4.  
    We investigate Kepler's population of cool (Fp < 20 ) Jupiters and identify Kepler-706 b and Kepler-468 b as candidates for cool inflated Jupiters. While Kepler-468 b orbits a young star and is potentially inflated due to its heat from formation, Kepler-706 b apparently orbits an old star based on its rotation period. Future observations will be required to rule out binary companions that could bias the radius measurement.
  • 5.  
    We find that 32 confirmed and 133 planet candidates have incident fluxes from 0.25 to 1.50 and thus occupy a nominal habitable zone (Kane et al. 2016). Of these, 37 planet candidates and eight confirmed planets have radii smaller than 2 (but see also Mullally et al. 2018; Burke et al. 2019); the confirmed planets are Kepler-62 e, Kepler-62 f, Kepler-186 f, Kepler-283 c, Kepler-440 b, Kepler-442 b, Kepler-452 b, and Kepler-1544 b.
  • 6.  
    We find that stars hosting multiple transiting planets have ages that are, on average, higher but statistically indistinguishable from those of single-planet host stars. This implies that if dynamical interactions (planet–planet and/or planet–passing star interactions) frequently scatter planets out of our line of sight to produce single transiting planet systems, these interactions are quick (e.g., Spalding & Batygin 2016; Obertas et al. 2017) and may happen over a wide range of ages.

Our results demonstrate the importance of precise, homogeneous parameters and the power of stellar ages and masses to allow a more comprehensive investigation of exoplanet populations and their evolution over time. An extension to lower-mass stars will require the use of alternative age indicators, such as lithium abundances, rotation, and UV-excess measurements, to identify a more robust sample of young Kepler hosts. Additionally, a homogeneous stellar classification of planet hosts observed by K2 and TESS will offer new insights into their planet populations and provide additional clues about planet formation and evolution.

We gratefully acknowledge everyone involved in the Gaia and Kepler missions for their tireless efforts that have made this paper possible. We thank Yanqin Wu, James Owen, and Edwin Kite for helpful comments that improved the paper. We thank Erik Petigura, Sam Grunblatt, Jamie Tayar, Ashley Chontos, Erica Bufanda, Maryum Sayeed, Connor Auge, Vanshree Bhalotia, Nicholas Saunders, Michael Liu, Benjamin Boe, and Ehsan Kourkchi for helpful discussions in addition to feedback on the figures. T.A.B. and D.H. acknowledge support from a NASA FINESST award (80NSSC19K1424) and the National Science Foundation (AST-1717000). D.H. also acknowledges support from the Alfred P. Sloan Foundation. E.G. acknowledges support from NSF award AST-187215 and, as a visiting professor at the University of Göttingen, the German Science Foundation through DFG Research 644 Unit FOR2544 "Blue Planets around Red Stars." This research was partially conducted during the Exostar19 program at the Kavli Institute for Theoretical Physics at UC Santa Barbara, which was supported in part by the National Science Foundation under grant No. NSF PHY-1748958. 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 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. This research has made use of NASA's Astrophysics Data System. This research was made possible through the use of the AAVSO Photometric All-Sky Survey (APASS), funded by the Robert Martin Ayers Sciences Fund. This research made use of the cross-match service provided by CDS, Strasbourg. 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 Exoplanet Exploration Program.

Software: astropy (Astropy Collaboration et al. 2013), dustmaps (Green et al. 2018, 2019), gapfit (Loyd et al. 2020), GNU Parallel (Tange 2018), isoclassify (Huber et al. 2017), KS2D (Taillon 2018), Matplotlib (Hunter 2007), mwdust (Bovy et al. 2016), Pandas (McKinney 2010), scikit-learn (Pedregosa et al. 2011), SciPy (Virtanen et al. 2020), skimage (van der Walt et al. 2014).

10.3847/1538-3881/aba18a
undefined