Kepler Multi-planet Systems Exhibit Unexpected Intra-system Uniformity in Mass and Radius

NASA's Kepler mission (Borucki et al. 2010) discovered an abundance of close-in ($P\lesssim 100$ days), tightly spaced multiple-planet systems (Latham et al. 2011; Lissauer et al. 2011, 2014; Rowe et al. 2014). These systems are nearly coplanar (e.g., Fang & Margot 2012), have low eccentricities (Van Eylen & Albrecht 2015), and are rarely in low-order mean-motion resonances (e.g., Fabrycky et al. 2014; Steffen & Hwang 2015). Recent work by Weiss et al. (2017) has further revealed that these multi-transiting systems tend to comprise planets that are regularly spaced and more similarly sized than planets drawn randomly from the aggregate population. Simultaneous to staying similarly sized, the radii also tend to increase slightly with orbital separation (Kipping 2017; Weiss et al. 2017).

A key question is whether this statistically significant "peas-in-a-pod" trend in planet radii extends to masses as well. Most planets comprising the Kepler multis are so-called super-Earth/sub-Neptune-type planets with ${R}_{p}=1\mbox{--}4{R}_{\oplus }$. These planets are common, occurring around roughly half of Sun-like stars (Fressin et al. 2013; Petigura et al. 2013). Because they bridge the gap between rocky and volatile-rich, super-Earths/sub-Neptunes exhibit a startlingly wide range of densities and plausible compositions (e.g., Rogers & Seager 2010). The variance in possible gas fractions yields planet radii that can be largely independent of mass (Lopez & Fortney 2014) and produces significant intrinsic scatter in the mass–radius relation (Weiss & Marcy 2014; Wolfgang & Lopez 2015; Wolfgang et al. 2016).

Given Weiss et al.'s (2017) observation of radius uniformity in multiple-planet systems, there are two options for the masses, each of which is unusual and consequential. Either the masses are dissimilar and the planets in a system are somehow capable of coordinating their radii, or the masses—like the radii—are similar, so that for a given system, the scatter in the mass–radius relation is significantly reduced. Our goal here is to determine which one of these scenarios dominates.

Most Kepler multis consist of stars that are too dim and/or planets that are too small and distant for reliable radial velocity (RV) mass measurements. However, if the planets in a system exhibit significant interplanetary gravitational perturbations—most often by being near a mean-motion resonance—the perturbations can add constructively and produce substantial transit-timing variations (TTVs; Agol et al. 2005; Holman & Murray 2005). Measuring TTVs for a multiple-planet system permits the estimation of the planet masses and eccentricities (see Agol & Fabrycky 2017 for a review).

Recently, Hadden & Lithwick (2017) conducted a comprehensive analysis of Kepler multiple-planet systems exhibiting significant TTVs (55 systems in total) in order to derive mass estimates. This data set offers a unique opportunity to assess the degree of intra-system similarity of planet masses within Kepler multis.³

We first examine the arrangement of Kepler multiple-planet systems in the mass–radius plane. Weiss et al. (2017) already demonstrated that planets in the same system tend to be similar in radius. We ask: does this hold for mass too?

Our sample selection begins with the 145 Kepler planets (55 systems) with masses derived by the TTV analysis conducted by Hadden & Lithwick (2017). They used two prior distributions for the masses: "default" and "high mass." We restricted the sample to the 37 systems for which all of the planets in a system had default and high-mass estimates that agreed within $2\sigma$ uncertainties, and we adopted the default estimates. The 37 multi-planet systems contain 89 planets; there are 26 systems with two planets, 8 with three, 2 with four, and 1 with five. Planet radii and host-star properties were taken from the California-Kepler survey (CKS) catalog⁴ (Johnson et al. 2017; Petigura et al. 2017); the CKS parameters were available for 73 of the 89 planets. When they were not available, or when the CKS and Hadden & Lithwick (2017) planet radii deviated by more than 100% (as they did in four cases), we used the Hadden & Lithwick (2017) radii for the discrepant planet and all of the others in the system. We also updated the planet mass estimates from Hadden & Lithwick (2017) by using the CKS catalog stellar masses.

Figure 1 displays the orbital architectures of 25 of the 37 systems, with points scaled to planet radii and coloration according to mass. The systems were selected based on a metric to be defined in the next section. Planets in a given system tend to be remarkably similar, not only in radius, but also in mass. In the following section, we will rigorously quantify this tendency for all 37 systems as an aggregate population.

Zoom In Zoom Out Reset image size

2.1. Clustering in the ${M}_{p}\mbox{--}{R}_{p}$ Plane

One way to display intra-system clustering in the M_p–R_p plane is with a set of directed graphs. The left panel of Figure 2 uses a directed graph for each of the 37 systems. All of the planets in a given system have the same color, and have lines connecting them in orbital period order. For clarity, the brightly colored systems are the 25 with the smallest normalized total distance of the graph, where the normalization is by the number of planets in the system. The remaining systems are plotted in light gray. Inspection of the left panel of Figure 2 makes it clear that planets in a given system tend to be closer together than expected if the systems were randomly assembled from planets in the population. That is, if we know the location of one planet on the M_p–R_p plane, we can guess with better than random accuracy where the other planets in the system will fall.

Zoom In Zoom Out Reset image size

Although this phenomenon of low M_p and R_p system dispersion is identifiable by eye, a quantitative metric is necessary for a rigorous analysis. The lengths of the lines in Figure 2 are a good measure of the dispersion of an individual system. To measure all of the systems at once, we can consider the total distance, ${ \mathcal D }$, of all of the systems' directed graphs in $\mathrm{log}{M}_{p}$–$\mathrm{log}{R}_{p}$ space. This may be expressed as

$$\begin{eqnarray}&&{ \mathcal D }=\displaystyle \sum _{i=1}^{{N}_{\mathrm{sys}}}\displaystyle \sum _{\begin{array}{c}j=1\\ {P}_{j}\lt {P}_{j+1}\end{array}}^{{N}_{\mathrm{pl}}-1}{\left[{\left(\mathrm{log}\displaystyle \frac{{R}_{j+1}}{{R}_{j}}\right)}^{2}+{\left(\mathrm{log}\displaystyle \frac{{M}_{j+1}}{{M}_{j}}\right)}^{2}\right]}^{1/2}.\end{eqnarray} \tag{ 1 }$$

Calculating the distance in log-space is beneficial because the resulting metric involves intra-system ratios of radii and masses, which are known to better precision than absolute values. In part, this is due to the elimination of the stellar radius and mass uncertainties. An additional boost comes from the fact that TTV-derived planet mass ratios are often more tightly constrained than individual masses (e.g., Jontof-Hutter et al. 2016; Hadden & Lithwick 2017).

This metric for system clustering on the M_p–R_p plane is only meaningful if provided with a control population for comparison. To this end, we constructed 50,000 control sets of shuffled systems. In each realization, the total number of systems and the number of planets in each system were conserved, but the planets themselves were randomly shuffled. To denote this, let ${ \mathcal M }=({M}_{1},{M}_{2},...,{M}_{N})$ and ${ \mathcal R }=({R}_{1},{R}_{2},...,{R}_{N})$ be the set of all planet masses and radii in the population. We define the symmetric group, S_N, as the set of all permutations of the set, $\{1,2,...,N\}$. If $X\in {S}_{N}$, then ${ \mathcal M }(X)$ and ${ \mathcal R }(X)$ represent a random shuffling of the planets, where mass/radius pairs are maintained. The set of distance metrics for the control population is therefore given by $\{{ \mathcal D }(X)={ \mathcal D }[{ \mathcal M }(X),{ \mathcal R }(X)]\ | \ X\in {S}_{N}\}$. The right panel of Figure 2 shows the directed graphs for one realization of randomly shuffled systems. Just like the left panel, the 25 systems with the smallest normalized total distances are plotted in color, with the others in gray. It is clear that the real systems are much less dispersed than the systems in this random realization. For instance, most of the gray systems in the left panel are significantly more clustered than the average among the colored, best 25 in the right panel.

The panels in Figure 2 display the histogram of $\{{ \mathcal D }(X)| \ X\in {S}_{N}\}$ for the 50,000 sets of shuffled systems. The value of ${ \mathcal D }$ for the real systems is indicated with a vertical green line. The distance (or equivalently, dispersion) metric for the real systems is quite remarkably small in comparison to the distribution of $\{{ \mathcal D }(X)| \ X\in {S}_{N}\}$ for the control population. It is smaller than the mean of the control set by $8.3\sigma$. We can also consider the individual mass and radius components of the distance, ${{ \mathcal D }}_{R}$ and ${{ \mathcal D }}_{M}$. Specifically we define

$$\begin{eqnarray}&&{{ \mathcal D }}_{R}=\displaystyle \sum _{i=1}^{{N}_{\mathrm{sys}}}\displaystyle \sum _{\begin{array}{c}j=1\\ {P}_{j}\lt {P}_{j+1}\end{array}}^{{N}_{\mathrm{pl}}-1}\left|\mathrm{log}\displaystyle \frac{{R}_{j+1}}{{R}_{j}}\right|,\end{eqnarray} \tag{ 2 }$$

with a similar expression for ${{ \mathcal D }}_{M}$. The values for the real systems for radius and mass are, respectively, $5.3\sigma$ and $7.6\sigma$ smaller than the control population means.

Compared to the random expectation, the extremely statistically significant reduction in the system dispersion metric, ${ \mathcal D }$, definitively indicates that planets in systems tend to be comparable to one another in mass and radius. Among various important implications, this is of prime relevance for understanding the arrangement of planets in the mass–radius diagram.

Many authors have noted that the mass–radius relation exhibits astrophysical or intrinsic scatter (Weiss & Marcy 2014; Wolfgang & Lopez 2015; Wolfgang et al. 2016). That is, at a particular radius, there is a wide distribution of possible masses, reflecting a diversity in compositions. This scatter is particularly pronounced for super-Earths/sub-Neptunes (Lopez & Fortney 2014). Our conclusion here is that astrophysical scatter in the mass–radius relation is dominated by—or at least largely linked to—system-to-system compositional dispersion rather than within-system dispersion.

2.2. Ordering in ${M}_{p}$ and ${R}_{p}$

Another property of interest is intra-system radius and mass ordering when the planets are ranked by orbital period. In a recent paper, Kipping (2017) examined the radius size ordering of observed Kepler multiple-planet systems and found that systems tend to be statistically significantly size-ordered (thereby occupying low-entropy configurations) when compared to random expectation. A schematic way of investigating this for the data set under consideration is to plot all of the M_p–R_p directed graphs in Figure 2 end-to-end. Drawing an analogy to a random walk, the first planet in a given system is co-located with the last planet in the preceding system.

Figure 3 presents this schematic. The colored line represents the random walk of the set of 37 real systems, ordered from low to high stellar mass. The black lines represent the random walks of 100 sets of shuffled systems. These sets were randomized slightly differently than in Section 2.1 because here we are isolating intra-system mass and radius ordering. Rather than assembling the systems from permutations of all planets in the population, we simply randomized the order of the planets in each system for each realization.

Zoom In Zoom Out Reset image size

The trajectory for the real systems tends upwards and to the right, confirming that, on average, the planets in these systems tend to be size-ordered in both mass and radius. At its end, the trajectory is also at the edge of the "cloud" of the control population. We can quantify ordering metrics, ${{ \mathcal O }}_{R}$ and ${{ \mathcal O }}_{M}$, as the total displacement of the trajectories in the radius and mass directions, respectively. Specifically, we define ${{ \mathcal O }}_{R}$ as

$$\begin{eqnarray}&&{{ \mathcal O }}_{R}=\displaystyle \sum _{i=1}^{{N}_{\mathrm{sys}}}\displaystyle \sum _{\begin{array}{c}j=1\\ {P}_{j}\lt {P}_{j+1}\end{array}}^{{N}_{\mathrm{pl}}-1}\mathrm{log}\displaystyle \frac{{R}_{j+1}}{{R}_{j}},\end{eqnarray} \tag{ 3 }$$

with a similar expression for ${{ \mathcal O }}_{M}$. Considering 50,000 sets in the control population, we find that ${{ \mathcal O }}_{R}$ for the real systems is $2.4\sigma$ larger than the control mean, and ${{ \mathcal O }}_{M}$ is $1.8\sigma$ larger.⁵

Although this is not the first time that planet size ordering in the Kepler multis has been shown (Ciardi et al. 2013; Kipping 2017; Weiss et al. 2017), here we have not only confirmed the observation of radius ordering but also extended it to mass.

We have shown that Kepler multi-transiting systems tend to host "the same" planets, that is, planets that are self-similar in their radii and masses. Consequently, if provided with all of the parameters of a host star and its natal protoplanetary disk, it should be possible (hypothetically speaking) to predict the most probable outcome of planet formation. What sets this default outcome? Here we use regression analyses to investigate which characteristics of host stars might be important in controlling the properties of their planetary systems.

We first expand our sample to a much larger collection of Kepler multis rather than focusing on the 37 from Hadden & Lithwick (2017). The trends in mass uniformity have already been examined, so we withhold planet masses from the analysis and focus instead on system characteristics involving orbital and planetary radii. The Weiss et al. (2017) "CKS multis" form our base data set. For consistency, we applied the same cuts as Weiss et al. (2017) and discarded known false positives, systems in Furlan et al. (2017) with dilution from nearby stellar companions greater than 5%, grazing transits with impact parameters, $b\gt 0.9$, and planets with $\mathrm{SNR}\lt 10$. The total sample contains 909 planets in 355 multi-planet systems. Before proceeding, we report the values of the CKS multis intra-system radius dispersion and ordering metrics. We find that ${{ \mathcal D }}_{R}$ and ${{ \mathcal O }}_{R}$ are $14.9\sigma$ smaller/$10.5\sigma$ larger than their respective control population means.

Multiple regression is well suited to the problem of correlating planet system characteristics with stellar properties. We used the ordinary least squares, general linear regression from the scikit-learn Python machine learning package (Pedregosa et al. 2011).⁶ The stellar properties (independent variables) that are potentially of interest and provided in the CKS catalog are M_⋆, R_⋆, ${T}_{\mathrm{eff}}$, metallicity, $\mathrm{log}g$, stellar age, and $v\sin i$. Due to their strong correlations with M_⋆, we ignored R_⋆, $\mathrm{log}g$, and ${T}_{\mathrm{eff}}$ to avoid complications from multicollinearity.

For summary planetary system characteristics, we considered the semimajor axis and equilibrium temperature of the innermost planet (quantifying the "inner edge" of the system), the median planet radius, and a measure of dispersion of the planet radii, ${\mathrm{CV}}_{{R}_{p}}/{\mathrm{CV}}_{{R}_{p},\mathrm{rand}}$. This dispersion measure is the coefficient of variation ($\mathrm{CV}=\sigma /\mu$) of the planet radii within a system normalized by that of an equal number of planets pulled randomly from the population.

The regressions reveal several weak but statistically significant correlations of planet parameters with M_⋆ and $\mathrm{Fe}/{\rm{H}};$ $v\sin i$ and stellar age had no significant effects. Interestingly, the correlations we identified strengthen significantly when considering systems with progressively higher multiplicity; before discussing that phenomenon in more detail, we will first summarize the results for the 128 systems with three or more planets. For ${a}_{\mathrm{inner}}$, M_⋆ and $\mathrm{Fe}/{\rm{H}}$ explain $\sim 22 \%$ of the variation (${R}^{2}=0.22$), with a positive correlation of ${a}_{\mathrm{inner}}$ with M_⋆ ($\beta =0.069\,\mathrm{au}/{M}_{\odot }$, $p\lt {10}^{-5}$) and a negative correlation with $\mathrm{Fe}/{\rm{H}}$ ($\beta =-0.057\,\mathrm{au}/\mathrm{dex}$, $p\lt {10}^{-5}$).⁷ The correlations are slightly stronger for ${T}_{\mathrm{eq},\mathrm{inner}};$ M_⋆ and $\mathrm{Fe}/{\rm{H}}$ together describe $\sim 30 \%$ of the variation in ${T}_{\mathrm{eq},\mathrm{inner}}$, and both have positive correlations.

For the variation in $\mathrm{median}({R}_{p})$, $\mathrm{Fe}/{\rm{H}}$ can alone explain ∼21% with a positive correlation ($\beta =2.2\,{R}_{\oplus }/\mathrm{dex}$, $p\lt {10}^{-7}$, ${R}^{2}=0.21$). M_⋆ produces no significant impact.⁸ Likewise, the metallicity is positively correlated with ${\mathrm{CV}}_{{R}_{p}}/{\mathrm{CV}}_{{R}_{p},\mathrm{rand}}$ ($\beta \,=0.7\ {\mathrm{dex}}^{-1}$, $p\lt 0.0005$, ${R}^{2}=0.1$).

We stress that this data set exhibits several mild violations of the assumptions of a multiple linear regression, specifically multicollinearity ($\mathrm{Fe}/{\rm{H}}$ and M_⋆ are correlated with ${R}^{2}=0.19$), nonlinearity, and heteroscedasticity. To address these, we performed non-parametric tests using Spearman's and Kendall's rank correlations. While there are (expected) differences in the correlation coefficients, all of the conclusions are consistent with the multiple regression analyses.

Figure 4 displays the correlations graphically. The most interesting features are the following.

• 1.  
  System inner edges trend closer to the star for higher metallicities. This is true despite the positive correlations in ${a}_{\mathrm{inner}}$ with M_⋆ and in $\mathrm{Fe}/{\rm{H}}$ with M_⋆.
• 2.  
  The median planet radii tend to increase with metallicity.
• 3.  
  The within-system radius dispersion is generally smaller than random expectation (as already discussed in Section 2.1). When size dispersion is large, it is typically associated with high metallicities, $\mathrm{Fe}/{\rm{H}}\gtrsim 0$.

Zoom In Zoom Out Reset image size

As mentioned previously, the correlations strengthen significantly with planet multiplicity. For example, the R² of the correlation between M_⋆ and $\mathrm{Fe}/{\rm{H}}$ (as independent variables) and ${a}_{\mathrm{inner}}$ (as the dependent variable) is 0.1 for systems with 2+ planets, 0.22 for $3+$, and 0.39 for $4+$. By sampling subsets of the $2+$ systems, we confirmed that this observation is not simply a consequence of the decrease in sample size with multiplicity.

We expanded on the findings of Weiss et al. (2017) to investigate the within-system dispersion in the planet masses and radii of Kepler multi-transiting systems. The conclusions are striking: the masses and radii in a given system are significantly more uniform as compared to random expectation (to $\sim 8.3\sigma$). Moreover, both masses and radii tend to be ordered when the planets are ranked by orbital period. In light of these new observations, the large astrophysical scatter of the mass–radius relation has suddenly been endowed with much more structure.

It is important to point out the selection biases that may be affecting our results. We considered a sample of 37 Kepler multis with TTV masses measured by Hadden & Lithwick (2017). The covariance between planet masses in the TTV posteriors (which we mentioned following Equation (1)) could work to slightly boost the observed phenomenon of intra-system mass uniformity. However, the correlations are not always positive, particularly for systems with more than two planets. Furthermore, our results are strictly only relevant for near-resonant systems with sizable TTVs. It is conceivable that these systems have planet masses that are significantly more uniform than the broader population of Kepler close-in, coplanar systems. However, this seems unlikely given that the tendency for systems to exhibit radius uniformity exists beyond this TTV sample; it would be unusual for the planets in a system to coordinate their radii but not their masses.

The correlations that we have found allow us to sketch a plausible framework of explanation, while keeping in mind that planet formation is likely a highly nonlinear process, and that any such argument cannot be construed as strictly deterministic.

Lee & Chiang (2015) outlined an analytic framework that links the planetary gas-to-core mass ratio (GCR) of a planet to protoplanetary nebular variables. For disk lifetime, τ, and gas metallicity, Z, their theory predicts GCR $\propto {\tau }^{0.4}{M}_{\mathrm{core}}^{1.7}{Z}^{-0.4}{(1-Z)}^{-3.4}$. In other words, enhanced mean molecular weight associated with a high metallicity envelope can permit higher GCRs, and hence larger planets. This insight, together with the observational evidence that disk lifetimes are shorter in low-metallicity environments (Yasui et al. 2010), is consistent with the positive correlation between the median system radii and stellar metallicities shown in Figure 4.

For disks with surface density power laws concordant with the Minimum Mass Extrasolar Nebula (Chiang & Laughlin 2013), solid body accretion simulations (e.g., Hansen & Murray 2013) show the rapid development of multiple solid cores with similar masses. A set of similar cores will accrete gas at comparable rates, thereby finishing the formation process with roughly equal masses and radii. The faster accretion at high metallicity, however, makes it more likely for some planets to approach runaway accretion, resulting in greater intra-system size dispersion. Furthermore, the central role of metallicity as a coolant in optically thin disk environments permits planet formation at high equilibrium temperature, which may contribute to the correlations that we have noted.

We thank the anonymous referee for insightful comments and suggestions. S.M. is supported by the NSF Graduate Research Fellowship Program under Grant DGE-1122492. S.W. acknowledges the Heising-Simons Foundation for their support. G.L. acknowledges NASA Astrobiology Institute support under Agreement #NNH13ZDA017C.