The Bimodal Color Distribution of Small Kuiper Belt Objects^*

From the apparent g and i magnitudes obtained through our photometric calibration (generally two in each band, unless masked owing to a problematic detection in one exposure), we computed error-weighted mean g and i magnitudes for each object, from which we subsequently derived the g − i color. The uncertainty in color ${\sigma }_{g-i}$ is given by the quadrature sum of the individual g and i mean magnitude errors.

Asteroid rotation during the approximately 4 hr observational arc of each detected object can entail an additional contribution to the color uncertainty. We followed the methods of Wong & Brown (2015) and utilized the repeated observations in each band to empirically determine the rotational contribution to the color uncertainty. The standard deviation of differences between two consecutive g-band magnitude measurements is composed of a contribution from the photometric uncertainties of the individual magnitude measurements alone (given by the quadrature sum of the individual magnitude uncertainties) and any additional contribution from magnitude modulation of the object due to rotation. When comparing the standard deviation of g-band magnitude differences (binned by 0.5 mag intervals) with the corresponding photometric error contribution, we found that the values were comparable; in other words, the photometric error contribution does not systematically underestimate the standard deviation values. The same conclusion was reached for i band. Therefore, we find that asteroid rotation does not have an appreciable effect on the color measurements overall.

Only two KBOs detected by our survey and matched with known objects have previous color measurements: 2000 QN251 and 2005 SC278. Both of these KBOs were observed and analyzed by Sheppard (2012), who reported g − i colors of 1.21 ± 0.07 and 1.51 ± 0.04, respectively. The colors that we measured are 1.73 ± 0.05 and 1.40 ± 0.03, respectively. While the colors measured for 2005 SC278 are comparable, our measured color for 2000 QN251 is much redder than the earlier measurement.

2000 QN251 is one of the brightest objects detected by our survey, with calculated $g,g,i,i$ magnitudes of 24.3, 24.4, 22.6, and 22.6, respectively. All four images had comparable seeing, and our computed absolute magnitude (derived from interpolation of g- and i-band magnitudes) matches the value listed by the Minor Planet Center. Our reported magnitude values preclude any rotational modulation in apparent magnitude between consecutive images taken in the same band. We suggest that the discrepancy between the measured colors in our work and Sheppard (2012) may be due to inconsistent observing conditions during the previous observation or a different observing cadence, which may have led to a rotational contribution to the Sheppard (2012) color measurement. We note that neither 2000 QN251 nor 2005 SC278 is a hot KBO, so they are not included in the color–magnitude analysis in this paper.

The distribution of colors g − i for all the hot KBOs detected by our survey is shown in Figure 3 with the unfilled histogram. The bimodality of the distribution is clearly discernible, with two color modes centered at around $g-i=0.9$ and $g-i=1.3$. Moreover, the bimodality is present throughout the entire magnitude range covered by our survey ($H=6.5-10.4$). This color distribution is indicative of two subpopulations within the small hot KBO population whose members have characteristically different surface colors. Hereafter, we refer to these subpopulations as the red (R) and very red (VR) KBOs.

Zoom In Zoom Out Reset image size

To carry out a more detailed analysis of the individual subpopulations, we first filter our survey KBO data set by removing the low-S/N objects. The photometric errors of the detected objects become large with decreasing apparent brightness; the color uncertainties are also increased for objects located in images that were obtained when the seeing was poor. A sizable fraction of the objects in the survey KBO data set have color uncertainties that are comparable to the difference between the two color modes. These large uncertainties contribute additional scatter to the color distribution and prevent clear categorization of those objects as R or VR. For the results presented below, we constructed a filtered KBO data set by only considering the 118 objects with ${\sigma }_{g-i}\leqslant 0.15$.

The color distribution of objects in the filtered KBO data set is shown by the solid histogram in Figure 3. To quantitatively determine the mean colors of the R and VR subpopulations, we followed the technique employed in our past studies of bimodal asteroid populations (Wong et al. 2014; Wong & Brown 2015) and fit the filtered KBO color distribution to a two-peaked Gaussian probability function. A standard log-probability minimization yielded the mean R and VR colors: ${\overline{(g-i)}}_{{\rm{R}}}=0.91$ and ${\overline{(g-i)}}_{\mathrm{VR}}=1.42$. The standard deviations of the two color modes are ${\sigma }_{g-i,{\rm{R}}}=0.12$ and ${\sigma }_{g-i,\mathrm{VR}}=0.09$, respectively.

This normal mixture model fitting also allows us to determine the significance of the bimodality. We computed the Bayesian Information Criterion (BIC) for both unimodal and bimodal distribution fits to the KBO color distribution. The BIC is defined as BIC $\equiv -2\mathrm{log}(L)+k\mathrm{log}(n)$, where L is the likelihood of the best-fit solution, k is the number of free parameters, and n is the number of data points. The difference in BIC between the unimodal and bimodal fits is ΔBIC = 23.0, which shows that the bimodal interpretation is very strongly favored (ΔBIC > 5 indicates a strong preference for the model with smaller BIC).

Having established the mean colors, we can now sort individual objects into the R and VR subpopulations based on their measured colors. We classified all objects in the filtered KBO data set with $g-i\leqslant {\overline{(g-i)}}_{{\rm{R}}}+{\sigma }_{g-i,{\rm{R}}}=1.03$ as R and all objects with $g-i\geqslant {\overline{(g-i)}}_{\mathrm{VR}}-{\sigma }_{g-i,\mathrm{VR}}=1.31$ as VR. While our previous studies used the mean colors as the cutoff locations for classifying objects into subpopulations, limiting the filtered KBO data set to just objects with ${\sigma }_{g-i}\leqslant 0.15$ enabled us to extend our categorization into the region between the two mean colors. As a result, we categorized 63 objects as R and 30 objects as VR.

Our moving object detection procedure returned the right ascension and declination of all objects in each of the four exposures taken at the corresponding fields. From here, we converted to geocentric ecliptic coordinates to obtain the ecliptic longitude λ and latitude β of every object. Having computed the heliocentric distance d through our orbit fitting routine (Section 2), we calculated the geocentric distance Δ of each object using standard coordinate transformations:

$$\begin{eqnarray}{\rm{\Delta }} & = & {r}_{S}\cos (\beta )\cos (\lambda -{\lambda }_{S})\\ & & +\sqrt{{r}_{S}^{2}[{\cos }^{2}(\beta ){\cos }^{2}(\lambda -{\lambda }_{S})-1]+{d}^{2}},\end{eqnarray} \tag{ 1 }$$

where r_S and ${\lambda }_{S}$ are the geocentric distance and geocentric ecliptic longitude of the Sun at the time of each observation. We first interpolated the fluxes in g band and i band to derive the V-band apparent magnitude, v, for each object, following the spectral slope computation described in Section 4.1. Next, we converted the apparent v magnitudes to absolute H magnitudes via

$$\begin{eqnarray}&&H=v-2.5{\mathrm{log}}_{10}({d}^{2}{{\rm{\Delta }}}^{2}).\end{eqnarray} \tag{ 2 }$$

For the identified KBOs in our data set with previously obtained visible-band or multiband photometry, our calculated absolute magnitudes agree with the Minor Planet Center magnitudes to within 0.3 mag.

Figure 4 shows the cumulative absolute magnitude distribution N(H) of all objects in our survey KBO data set (black points). The gentle rollover starting at $H\sim 8.0$ indicates the onset of detection incompleteness in our survey. The magnitude distribution of KBOs has been the subject of intensive study over the past decade, with numerous dedicated surveys having been undertaken to explore the shape of the distribution through $H=11\mbox{--}12$. Many recent published works in the literature present incompleteness-corrected KBO magnitude distributions extending well beyond the completeness limit of our photometric survey (e.g., Fraser et al. 2014; Shankman & Kavelaars 2016).

Zoom In Zoom Out Reset image size

For the sake of comparison with previous studies of the hot KBO population, we chose a conservative threshold magnitude of ${H}_{\max }=7.5$ and fit the differential magnitude distribution, ${\rm{\Sigma }}(H)={dN}(H)/{dH}$, of our survey KBOs through ${H}_{\max }$ with a single power-law distribution of the form ${\rm{\Sigma }}(\alpha ,{H}_{0}| H)\,={10}^{\alpha (H-{H}_{0})}$, where α is the power-law slope and H₀ is the normalization magnitude for which ${\rm{\Sigma }}({H}_{0})=1$. Following the methods described in Wong & Brown (2015), we used a Markov chain Monte Carlo algorithm to estimate the best-fit parameters and corresponding $1\sigma$ uncertainties, which yielded $\alpha ={1.45}_{-0.37}^{+0.40}$. This value is steeper than the slope of ${0.87}_{-0.2}^{+0.07}$ reported by Fraser et al. (2014) for the same magnitude range, though still consistent at the $1.5\sigma$ level.

Our analysis of the g − i color distribution and the classification of individual hot KBOs as R or VR allow us to carry out the first study of the individual magnitude distributions of the R and VR color subpopulations. The cumulative absolute magnitude distributions of the categorized R and VR KBOs are shown in Figure 4. It is immediately apparent that the R and VR magnitude distributions have nearly identical shapes throughout the entire magnitude range covered by our survey. The similarity between the R and VR magnitude distributions has important implications regarding the origin of the color bimodality, which we discuss in the next section.

Although not all KBOs in our filtered sample were classified as R or VR owing to their intermediate measured colors, this incompleteness in categorization is magnitude independent and therefore has no effect on the relative shapes of the color–magnitude distributions. Meanwhile, nonuniform seeing and differing g- and i-band S/N can cause survey biases. Among the objects in our filtered sample, the majority had the worst seeing during g-band observations (see supplemental table). Overall, when the seeing during g-band observations and that during i-band observations were comparable, the S/N was slightly lower in g band. Both of these trends establish g as the limiting band for detection. For a KBO with a given absolute magnitude near the detection limit, R objects are roughly 0.1 mag brighter in g band than VR objects.

Thus, we conclude that our survey has a bias toward R objects. We note, however, that this bias only affects the faintest objects in our sample. Therefore, the relative shape of the R and VR magnitude distributions and the R-to-VR number ratio among the brightest objects should not be influenced by the bias.

We used mean R and VR geometric albedos from the literature (${p}_{v,{\rm{R}}}=0.06$ and ${p}_{v,\mathrm{VR}}=0.12$, respectively, as computed in Fraser et al. 2014) to convert absolute magnitudes to approximate diameters via the relation $D=1329\,\times {10}^{-H/5}/\sqrt{{p}_{v}}$ km. To obtain our first estimate of the R-to-VR number ratio r, we considered only categorized KBOs larger than 100 km in diameter, corresponding to threshold absolute magnitudes of H = 8.7 and H = 7.9 for R and VR objects, respectively. In this range, we detected 36 R and 10 VR objects, corresponding to r = 3.6. Since our method for classifying objects into subpopulations relied on the results of a double-Gaussian fit (see Section 3.1), the number ratio we calculate can be affected by deviations in the measured color distribution from a true two-peaked Gaussian shape (e.g., slight asymmetries relative to mean R and VR colors).

As a second estimate of the R-to-VR number ratio, we calculated the error-weighted mean color of all 46 objects in our filtered KBO sample larger than 100 km: $\overline{g-i}=1.04$. We modeled the mean color as a mixture of R and VR objects with typical colors ${\overline{(g-i)}}_{{\rm{R}}}=0.91$ and ${\overline{(g-i)}}_{\mathrm{VR}}=1.42$, respectively: $\overline{g-i}=q\times {\overline{(g-i)}}_{{\rm{R}}}+(1-q)\times {\overline{(g-i)}}_{\mathrm{VR}}$, where q is the R object fraction and is related to r via $r=q/(1-q)$. Through this procedure, we obtained r = 3.0. Acknowledging the sensitivity of the computed number ratio to the particular method of calculation, we give the general conclusion that the R hot KBO subpopulation is roughly three to four times the size of the VR subpopulation.

Constraining for relatively high S/N objects via the condition ${\sigma }_{g-i}\leqslant 0.15$, we obtained a filtered KBO set of 118 objects (Section 3.1). Applying the same constraint to the 16 Centaurs detected by our survey yields 15 objects.

In order to combine our data with previously published color measurements of Centaurs and hot KBOs compiled in the literature, we must convert our g − i colors to spectral slopes S. We followed the prescription described in detail in Wong et al. (2014) and first converted the measured apparent magnitudes to flux, normalizing to unity in g band: ${F}_{g}=1$, ${F}_{i}={10}^{0.4((g-i)-0.55)}$, where we have removed the solar g − i color (Ivezić et al. 2001). The associated relative flux uncertainties were computed using the second-order method of Roig et al. (2008): ${\rm{\Delta }}{F}_{g}/{F}_{g}\,=\sqrt{2}{\sigma }_{g}$, ${\rm{\Delta }}{F}_{i}/{F}_{i}\,=0.9210{\sigma }_{g-i}(1+0.4605{\sigma }_{g-i})$. We then calculated the spectral slope S using standard linear regression.

We queried the database of objects analyzed in Peixinho et al. (2015) for hot KBOs (inclinations greater than or equal to 5°; Section 2) and Centaurs, both with the condition $H\geqslant 7$, finding 82 and 38 objects, respectively. Combining these objects with our survey data, we constructed a total spectral slope data set for 200 hot KBOs and 53 Centaurs. The corresponding spectral slope distributions are shown in Figure 5.

Zoom In Zoom Out Reset image size

We note that the transformation from g − i color to S is not linear, with the distance between spectral slope values corresponding to evenly spaced g − i color values being amplified as the color becomes redder. As a result, the shape of the VR mode in the spectral slope distribution is broad and noticeably asymmetric, as opposed to the shape of the analogous peak in the g − i color distribution (Figure 3). This serves to reduce the saliency of the bimodality in the spectral slope distribution. However, even with this effect, the bimodality in spectral slope among the hot KBOs is still discernible.

The bimodality in spectral slope among the Centaurs is well attested by numerous previous studies (e.g., Peixinho et al. 2003; Barucci et al. 2005; Perna et al. 2010). The comparison of Centaurs with hot KBOs demonstrates not only that the hot KBO spectral slope distribution is also bimodal, but also that the two color modes closely correspond in location to the two color modes in the Centaur spectral slope distribution. In other words, the R and VR Centaurs (often referred to as gray/neutral and red in the literature) have the same colors as R and VR hot KBOs.

From this observation, we can conclusively resolve the issue of the color bimodality of Centaurs in this magnitude range (which includes all known Centaurs except Chariklo and Chiron): before being scattered inward onto their current comet-like orbits, R and VR Centaurs were members of the R and VR hot KBO subpopulations, respectively, as previously suggested by Fraser & Brown (2012) and Peixinho et al. (2012).

We also place the color data of small KBOs from our survey in the context of the wider hot KBO population as a whole. Figure 6 shows the color–magnitude distribution of all hot KBO and Centaurs with spectral slope errors less than $10\times {10}^{-5}\,{\mathring{\rm A} }^{-1}$. The plot illustrates how our survey data extend the known body of high-quality KBO color measurements into the size region that had previously only been adequately studied using Centaur colors as proxy.

Zoom In Zoom Out Reset image size

For both hot KBOs and Centaurs, a bimodality in color is evident at all magnitudes among objects fainter than $H\sim 7$, as was discussed in Section 4.1. The Subaru hot KBO data nearly triple the number of high-S/N KBO color measurements between H = 7 and H = 9 contained in the Peixinho et al. (2015) database and are instrumental in solidifying the observation of bimodality in that size range, which was previously hinted at using far fewer data points (Fraser & Brown 2012; Peixinho et al. 2012, 2015).

Notably, however, the colors of larger, brighter hot KBOs show a more uniform distribution, with many objects having spectral slope values intermediate between the two color types. The apparent clustering of near-neutral objects between H = 2 and H = 6 is primarily due to the presence of the Haumea collisional family. Previous studies have shown that, omitting the Haumea family members, no significant bimodality is present among objects brighter than $H\sim 6$ (Fraser & Brown 2012; Peixinho et al. 2012, 2015). In other words, only small hot KBOs are bimodal in color.

The threshold of $H\sim 6$ corresponds to a diameter of roughly 300–400 km. Observations of these large and intermediate-sized objects have revealed somewhat higher albedos than their smaller counterparts (Stansberry et al. 2008). This suggests that the larger KBOs have fundamentally different surface properties than the objects covered by our survey. The discrepancy may be due to the effects of secondary processes since formation, such as the onset of cryovolcanism or differentiation, which would not have occurred in smaller bodies (see, e.g., discussion in Fraser & Brown 2012).

Recent theories of solar system evolution posit that a mean motion resonance crossing between Jupiter and Saturn roughly 600 Myr after the formation of the solar system led to significant alterations in the orbital architecture of the giant planets (Gomes et al. 2005; Morbidelli et al. 2005; Tsiganis et al. 2005). An important consequence of this dynamical restructuring is the scattering of the primordial body of planetesimals that were situated beyond the orbits of the ice giants. While the majority of the mass within this region was expelled from the solar system as the orbits of the ice giants expanded, a fraction of the planetesimals were scattered outward to become the present-day Kuiper Belt. The remainder was scattered inward and eventually captured into resonance with Jupiter to become the Trojan asteroids. Therefore, the current paradigm of solar system evolution predicts a common origin for Jupiter Trojans and KBOs.

Several authors have previously noted that the overall magnitude distributions of Trojans and hot KBOs brighter than $H\sim 8-9$ are similar (e.g., Fraser et al. 2014). In Wong et al. (2014), we reported a power-law slope¹ of 1.11 ± 0.20 for Trojans brighter than $H\sim 8$. In our present study, we obtained a best-fit magnitude distribution slope of ${1.45}_{-0.37}^{+0.40}$ for hot KBOs brighter than H = 7.5. The comparison demonstrates that the bright hot KBOs and Trojans of the same size have statistically indistinguishable magnitude distributions.

Similar to the KBOs, the Trojans display a robust optical color bimodality, albeit with the peaks of the modes centered at bluer colors ($g-i=0.73$ and 0.86; Wong & Brown 2015). In Wong et al. (2014), we categorized Trojans into the two respective color subpopulations and analyzed the individual color–magnitude distributions. We showed that the two Trojan subpopulations have statistically indistinguishable magnitude distributions in the bright end ($H\lt 9$), just as the R and VR hot KBO subpopulations have nearly identical magnitude distribution shapes throughout the brightness range covered by our survey (Section 3.2).

Taken together, the analyses of Trojans and KBOs provide strong evidence to support one of the key predictions of current models of solar system evolution. Not only are the overall magnitude distributions of the hot KBOs and Trojans similar, but the individual constituent color subpopulations have comparable magnitude distributions. Therefore, we can posit that each of the two subpopulations in the Trojans and hot KBOs is sourced from a single progenitor body of planetesimals.

The bimodal color distribution of hot KBOs is indicative of two different surface types, with distinct photometric properties. There are two classes of possible explanations for the origin of the observed color bimodality: (1) the R and VR KBOs formed in different regions of the protoplanetary disk and were sourced from material with different bulk compositions, and (2) the R and VR KBOs formed in the same region of the disk and out of the same material, but later developed distinct surface properties via some secondary evolutionary process.

The size distribution, or by proxy, the magnitude distribution, of a population contains information about both the accretion environment and the subsequent collisional evolution of the population. Therefore, the color–magnitude distributions we derived in our analysis are a crucial tool in helping us understand how the color bimodality developed.

In Figure 4, we showed that the R and VR hot KBO magnitude distributions are nearly identical in shape, which suggests that the two subpopulations have undergone a similar formation and evolutionary history and therefore likely formed out of the same body of material within the protoplanetary disk. This similarity in magnitude distribution shape contrasts with, for example, the strong discrepancy between the hot and cold KBO magnitude distributions (e.g., Fraser et al. 2014).

In Wong & Brown (2016), we proposed a simple hypothesis that accounts for the development of a color bimodality within the primordial planetesimal disk beyond the orbits of the ice giants. Assuming the dynamical framework of the Nice model, the model describes a primordial body of material situated between 15 and roughly 30 au, from which the modern hot KBO population is sourced (e.g., Morbidelli et al. 2005; Levison et al. 2008). After a residence time of roughly 600 Myr (Gomes et al. 2005), a period of chaotic dynamical evolution throughout the outer solar system scattered a fraction of these objects into the current Kuiper Belt region (see discussion in previous subsection).

The hypothesis posits that all the bodies within the primordial trans-Neptunian disk accreted with a roughly cometary composition, i.e., rocky material and water ice, with a significant volume of other volatile ices such as carbon dioxide, methanol, ammonia, hydrogen sulfide, etc. During the accretion process, large objects migrated differentially through the disk (e.g., Kenyon et al. 2008), so that the overall disk was well mixed compositionally, with a single characteristic size distribution. Once these bodies were exposed to insolation after the dispersal of the gas disk, sublimation-driven loss of volatiles from the outermost layers led to the development of distinct surface compositions throughout this region.

The loss rate of each individual volatile ice species is correlated with its vapor pressure, which in turn has a strong exponential dependence on temperature. As a result, each volatile ice species would have a sublimation line at a different location within the primordial disk: objects situated inward of the line became depleted in that species on their surfaces, while objects situated beyond the line retained that species on their surfaces.

From our modeling, we found that the sublimation line of hydrogen sulfide ice (H₂S) passed through the region of the primordial trans-Neptunian planetesimal disk. Therefore, objects within this single primordial population would have been divided into two groups—those that retained H₂S on their surfaces, and those that did not. Previous experimental work has shown that irradiation of volatile ices leads to a general reddening of the optical spectral slope (e.g., Brunetto et al. 2006). We posited that the presence of H₂S on the surface of an object would contribute additional reddening relative to the case where H₂S was absent. H₂S is known to induce a strong reddening of the surface upon irradiation, as has been observed in the H₂S frost deposits on Io (Carlson et al. 2007).

After the Nice model dynamical instability, some objects in the trans-Neptunian region were scattered outward, sharing a common size distribution by virtue of their common formation environment while inheriting the surface color bimodality that developed in the primordial planetesimal disk owing to differential loss/retention of H₂S. Likewise, some objects were scattered inward to eventually settle into resonance with Jupiter, thereby explaining the observed color bimodality in the Trojans and the similarity in magnitude distributions between the color subpopulations of the Trojans and hot KBOs. In Wong & Brown (2016), we speculate that the bluer colors of the two Trojan subpopulations may be due to secondary processing of the irradiated surfaces following their placement in a higher-temperature environment at 5.2 au.