OSSOS. II. A SHARP TRANSITION IN THE ABSOLUTE MAGNITUDE DISTRIBUTION OF THE KUIPER BELT'S SCATTERING POPULATION

The Survey Simulator determines whether a given object would have been detected and tracked by one of our surveys. The simulator is given a list of survey pointings and the detection efficiency for each pointing in order to perform a simulated survey. A randomly selected model object, with an assigned H-magnitude and color, is tested for detection by the survey simulator. The simulator first checks that the object is bright enough to have been seen in any of our surveys' fields, then checks that the object was in a particular field, and that it was bright enough to be detected in that field. Based on a model object's simulator-observed magnitude and the field's detection and tracking efficiencies, model objects are assessed for "observability." The survey simulator reports the orbital parameters, the specified H-magnitude, and color conversion for all orbital model objects determined to have been "detected" and "tracked." The object's "observed" magnitude, and corresponding H-magnitude, which includes accounting for measurement uncertainties, are also reported. Using the Survey Simulator, we produce a statistically large model sample that has been biased in a way that matches the biases present in our observed survey sample. This large Survey Simulator produced model sample is then compared directly to the detected sample.

In order to carry out a simulated survey, one requires an orbital model for the objects being "observed." For our model we select out the scattering TNOs from a modified version of the Kaib et al. (2011, KRQ11) orbital model of the TNO population. The KRQ11 model (see Figure 1) is the end-state of a dynamical simulation of the evolution of the solar system that includes the gravitational effects of the giant planets, stellar passages, and galactic tides. The simulation begins with an initial disk of massless test particles between semimajor axis a = 4 AU and a = 40 AU following a surface density proportional to a^−3/2, eccentricities, e < 0.01, and inclinations, i, drawn from sin(i) times a Gaussian, as introduced by Brown (2001). The giant planets are placed on their present-day orbits (see Section 5.1 for a discussion on the effects of the planet configuration and how it does not affect our result), the stellar neighborhood is modeled assuming the local stellar density, and the effects of torques from galactic tides are added (for more detail, see Kaib et al. 2011). This system is then integrated forward in time for 4.5 Gyr, resulting in a model for the current state of bodies in the outer solar system. The resulting orbital distribution is then joined with a candidate H-distribution, and a TNO color distribution derived from Petit et al. (2011). This joint orbit, H-distribution and color distribution model forms the input for the Survey Simulator.

Zoom In Zoom Out Reset image size

Shankman et al. (2013) demonstrate that the inclinations in KRQ11 are too low to match the observed scattering TNOs; the model's assumed initial inclination distribution is too "cold." A model with a "hotter" initial i distribution (with Gaussian width σ = 12°) was created and dynamically evolved forward for the age of the solar system. We continue to use this modified KRQ11 as our orbital model in this analysis. In Section 5.1 we discuss the effect that the choice of model has on our analysis.

We use the modified KRQ11 orbital model as a representation of the current-day population of scattering objects in the solar system in order to perform a simulated survey of scattering TNOs. We select out the scattering TNOs from the orbital model (i.e., those with Δa ≥ 1.5 AU in 10 Myr) in plausibly observable ranges (a < 1000 AU, pericenter q < 200 AU), which results in a relatively small number of objects. To account for the finite size of the model, we draw objects from the model and add a small random offset to some of the orbital parameters. The scattering TNOs do not have specific orbital phase space constraints (unlike the resonances which are constrained in a, e, and resonant angles) and thus we can better sample the space the model occupies by slightly adding this small random offset. This extends the model beyond the set of TNOs produced in the original run (∼29k for above a, q slices) that was necessarily limited by computation times. We resample a, q, and i by randomly adding up to ±10%, ±10%, and ±1°, respectively, to the model-drawn values. We also randomize the longitude of the ascending node, the argument of pericenter, and the mean anomaly of each model object. This modified KRQ11 orbital model, resampled to increase its utility, is used as the input orbital model for the scattering TNOs.

A variety of forms have been used to try to match the observed magnitude or H-magnitude distributions of various TNO populations. Single slopes (e.g., Jewitt et al. 1998; Gladman et al. 2001; Fraser & Kavelaars 2008; Gladman et al. 2012; Adams et al. 2014), knees (Fuentes & Holman 2008; Fraser & Kavelaars 2009; Fraser et al. 2014), knees with smooth rollovers (Bernstein et al. 2004), and divots (Shankman et al. 2013; Alexandersen et al. 2014) have all been proposed. Here we present a generalized form of the H-magnitude distribution for testing, which in limiting cases becomes either a knee or a single-slope.

We characterize the H-distribution by four parameters: a bright (large-size) slope, α_b, a faint (small-size) slope, α_f, a break absolute magnitude, H_b, and a contrast, c, that is the ratio of the differential frequency of objects just bright of the break to those just faint of the break. Depending on the parameters, our H-distribution takes one of three forms: single-slope, knee, or divot (schematic shown in Figure 2). We model the transition in the H-distribution as an instantaneous break, as the sample size of our observation does not merit constraining the form of a potential rollover. All H-distributions and values given in this work are presented in H_r unless otherwise specified. Our formulation of the H-distribution allows for the testing of the proposed H-distribution from a single framework.

Zoom In Zoom Out Reset image size

In the literature, the single slope distribution has been referred to as a single power-law distribution because it corresponds directly to the theorized distributions of diameters, D, which is a power law to an exponent, q : $\frac{{dN}}{{dD}}\propto {D}^{-q}$. This distribution is convertible to absolute magnitude and parameterized by a logarithmic "slope" α : $\frac{{dN}}{{dH}}\propto {10}^{\alpha H}$, with q = 5α + 1. H-distributions can be mapped to D-distributions, with an albedo, providing an observable way to probe size-distributions.

In order to create a synthetic H-distribution sample with a transition, one randomly samples from the two single-slope distributions (bright and faint ends), choosing which to sample from according to the fraction of the total distribution each section comprises. If the bright end of the distribution accounts for 60% of the whole distribution, then when randomly drawing objects, 60% of the time they should be drawn from a single-slope H-distribution corresponding to the bright distribution, and thus 40% of the time they should be drawn from the faint end distribution. The ratio of the number of objects in the bright end of the distribution to the total distribution depends on the two slopes, the contrast, the break magnitude, and the faintest magnitude, H_faintest; no normalization constants or knowledge of population size are required. The ratio can be calculated with (see Shankman 2012 for a derivation)

$$\begin{eqnarray}&&\displaystyle \frac{{N}_{\mathrm{bright}}}{{N}_{\mathrm{total}}}={\left(1+\displaystyle \frac{{\alpha }_{b}}{{\alpha }_{f}}\displaystyle \frac{1}{c}({10}^{{\alpha }_{f}({H}_{\mathrm{faintest}}-{H}_{b})}-1)\right)}^{-1}.\end{eqnarray} \tag{ 1 }$$

The H-distribution has four parameters: α_b, α_f, H_break, and c; with arguments from other TNO populations, we fix two of these parameters. Our sample of large-size objects is not large enough to measure the large-end slope, α_b. Motivated by other hot TNO populations, we fix α_b = 0.9, which matches our bright-end sample and is consistent with the slopes found for the hot Classical belt (Fuentes & Holman 2008; Fraser & Kavelaars 2009; Petit et al. 2011), the aggregated hot population (Adams et al. 2014; Fraser et al. 2014), the 3:2 resonators (Gladman et al. 2012; Alexandersen et al. 2014), and the pre-transition scattering plus Hot Classical TNOs (Adams et al. 2014). As in Shankman et al. (2013), we fix the break location to H_r = 8.3 (corresponding to H_g = 9.0 and D ≈ 100 km for 5% albedo). A single-slope H-distribution of α = 0.9 to the break magnitude of H_r = 8.3 is not rejectable at even the 1σ level by our sample, and thus provides good agreement to our observations. A steep slope and break at H_r = 8.3 is consistent with the Adams et al. (2014) scattering and Hot Classical TNO H-distribution with α = 1.05 measured down to H_r ∼ 6.7 that did not find a break, and the hot population breaks of H_r = 8.4 and H_r = ${7.7}_{-0.5}^{+1.0}$ found by Fraser et al. (2014). Our detected sample near the transition is not large enough to constrain the exact position of the break; moving the break location by several tenths of magnitude does not affect the conclusions of this work. Having fixed two parameters, we are left with α_f and c as free parameters to test and constrain.

Because CFEPS was primarily performed in g-band and Alexandersen et al. (2014) and OSSOS were performed in the r band, we must account for conversion between these two filters. The majority of our detections (13 versus 9) are r-band detections and so we choose to convert everything to r, requiring us to adopt a g − r color conversion for both the survey simulator and our detections. We do not have measurements of g and r for all of our objects, and so we adopt the CFEPS reported color conversion of $g-r=0.7$ (Petit et al. 2011) which is consistent with the color measurements we do have for our scattering TNOs. This conversion allows us to combine the samples from several surveys.

We convert all of our observed scattering TNOs into the r band with the above conversion. We select colors for our modeled objects from a Gaussian g − r distribution centered at 0.7 with a standard deviation of 0.2 (the range seen in the CFEPS object sample with g − r available). We test our analysis by choosing extreme values for our color conversions (i.e., 0.5 and 0.9). The effects of these color choices do not alter the conclusions of the analysis (see Section 5.1).

While the survey simulator allows for the modeling of light curve effects, we do not model the unknown light curves of our scattering TNO sample. Each of our detections is measured at a random phase of its light curve and the magnitude of light curve variation is small in comparison to the uncertainty in converting between detection band passes, which we show has no meaningful impact on our analysis (see Section 5.1). We conclude that light curves do not affect our analysis.

We want to test whether our model (orbital distribution, H-magnitude distribution, and observation biases) is a good representation for our observations. Having generated a set of simulated detections for a candidate H-distribution, we then test the hypothesis that our detected sample can be drawn from the simulated detections. If this hypothesis is rejected, the model used to generate our simulated detections is rejected; in particular, we conclude that this is a rejection of the candidate H-magnitude distribution used. We use a variant of the AD statistical test (Anderson & Darling 1954) to assess if our observations can be drawn from our simulated detections for a candidate H-magnitude distribution.

The AD test is a variant of the Kolmogorov–Smirnov two-sample test, weighted such that there is greater sensitivity to the tails of the distribution (Anderson & Darling 1954). The test metric can be thought of as the distance between two cumulative distributions. Typically, the test is used to determine if data are consistent with being drawn from a well known distribution (e.g., normal, uniform, lognormal) for which there is a look-up table of critical values for the AD metric. In our case, our simulated detections take the place of the well known distribution, and we have no look-up table of critical values for the metric. We therefore bootstrap our sample, repeating draws of 22 objects (the same as our observed sample), calculating their AD distance from their parent distribution in order to build a distribution of critical values for the simulated sample distribution. If the AD metric for our observed sample is in the tails (3σ of the bootstrap-built distribution for our simulated detections), we reject the hypothesis that our observations could have been drawn from the candidate simulated detection distribution and we thus reject the candidate H-distribution used.

We apply the AD test across several variables in our data set. We have a set of observed scattering TNOs which can be characterized by their orbital parameters and their observable characteristics. We test our model on a combination of three orbital parameters (semimajor axis, inclination, pericenter) and three observable parameters (magnitude at detection, distance at detection, H-magnitude). We sum the AD metrics for all of our tested distributions, and use this summed metric to test for the rejection of our model. While bootstrapping to build up our distribution of critical values as described above, we sample out 22 objects and simultaneously compute the AD metric for each distribution, as opposed to testing the distributions independently. This preserves the relationships between the parameters; in essence, we calculate the dimensionless AD metric for a set of 22 orbits across our six observational and dynamical characteristics to build a distribution of AD values. We use the sum of the AD metrics, as the sum will be small if all of the distributions are in good agreement, large if one is in poor agreement or several are in moderately poor agreement, and larger still if all are in poor agreement (Parker 2015). This approach has less power than rejecting based on the worst failing tested distribution, as has been commonly done in the literature, but is a more robust and multivariate approach.

As this work is model dependent, it is important to examine the effects of model choice. The scattering TNOs are insensitive to the history of the exact number and configuration of planets in the outer solar system (Shankman et al. 2013) as long as the end-state is the current planet configuration; the scattering TNO interactions are so disruptive and chaotic that the exact history is "erased." We have shown above that the KRQ11 i-distribution is discrepant with our observations, and we thus test the effects of the i-distribution on this work. We perform the same analysis as above on the original version of the KRQ11 orbital model, which has a colder initial i-distribution that is then evolved forward for the age of the solar system. Figure 6 shows the simulator-observed cumulative distributions for both the initially hot (green) and initially cold (blue) models with the same H and color distributions. While the i and a distributions vary dramatically between the two orbital models, the m, q, d, and H distributions only show small variations that do not statistically differentiate the two models. Figure 7 shows the confidence contours for the hot model (solid contour as in Figure 5) and cold model (blue line overlays) when testing the AD statistic on the m, q, d, and H distributions. There is only minor tension between these two models which have good agreement at the 3σ level. The observed objects Drac (peculiar i) and MS9 (peculiar a) have distinct orbital elements, but are classified as scattering TNOs by the Gladman et al. (2008) criterion. Given their unusual orbital elements, we consider that they may not be members of the scattering TNO population and we examine the effect of their inclusion in the analysis. Drac and MS9 have no effect on the conclusions of the analysis for the a and i distributions. Removing Drac and MS9 from the analysis also removes the tension between the hot and cold orbital models when testing on m, q, d, and H and has no effect on the conclusions of this work regarding the a and i distributions. Our analysis is not strongly dependent on the i-distribution of the orbital model when testing on the m, q, d, and H distributions.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

We also test the effects of the g − r color distribution used. In our analysis we draw model object g − r colors from a Gaussian centered at 0.7 with a standard deviation of 0.2, based on the measured CFEPS g − r colors (Petit et al. 2011). We test the effects of our adopted g − r distribution by performing the same analysis with fixed g − r values at the 1σ extremes of our adopted distribution. Figure 8 shows the results of performing the same analysis with color conversions of g − r = 0.5 and g − r = 0.9. There is no tension between rejectable (α_f, c) parameters for these two very different g − r color values. This analysis is not sensitive to exact knowledge of the intrinsic color of the TNOs in our sample, within reasonable bounds on those colors.

Zoom In Zoom Out Reset image size

We test a combination of the orbital model, the color distribution and the H-distribution. We have shown that our analysis is not strongly dependent on the choice of color distribution nor on orbital model and thus when a candidate model is rejected, it is a rejection of the modeled H-distribution.

Motivated by the measured JFC slopes (see Solontoi et al. 2012, Table 6), the Centaur slopes (Adams et al. 2014), and the theoretical slope of collisionally ground populations (O'Brien & Greenberg 2005), all of which suggest α_f ∼ 0.5, we favor α_f = 0.5. Our analysis has constrained the scattering TNO H-distribution to contrast values ranging from 1 (knees) to ∼32 (divots) for α_f = 0.5. Here we consider the implications of a shared H-distribution among the hot populations of the outer solar system. We use the observed Neptune Trojan and Plutino samples and the JFC source problem to inform the scattering TNO c value.

Requiring the scattering TNOs to be the source of the JFCs sets strong lower limits on the size of the scattering TNO population (Volk & Malhotra 2008), which can be used to constrain the form of the H-distribution. Volk & Malhotra (2008) find (0.8–1.7) × 10⁸ D > 1 km scattering TNOs in the 30–50 AU heliocentric distance range are required in order for the scattering TNOs to be the JFC source region. To apply this constraint to our measured H-distribution values, we extend our knee and divot population estimates down to D ∼ 1 km, using albedos of 5% and 15%, the range of measured albedos for the scattering TNOs and the Centaurs (Santos-Sanz et al. 2012; Duffard et al. 2014). H-distributions with c ≤ 10 and α_f = 0.5 produce a scattering TNO population that is large enough to act as the source of the JFCs for both 5% and 15% albedos. For c = 1 (knees), our sample requires α_f ≥ 0.4 for 5% albedo and α_f ≥ 0.5 for 15% albedo. If the H distribution transitions to a shallow slope of α_f ∼ 0.2 as seen in Fraser et al. (2014) and not rejected by this analysis, it must then transition to a third, steeper slope, α_f2, at a smaller H in order for the scattering TNOs to be numerous enough to act as the source population of the JFCs. A transition to α_f ∼ 0.2 in the scattering TNO H-distribution requires two transitions, the faint one at a break that has yet to be observed, and three slopes in order to explain the source of the JFCs. The H-distribution of the scattering TNOs must transition to α_f ≥ 0.4 slope to provide a large enough scattering TNO population to act as the source for the JFCs.

Under the hypothesis that all of the hot TNO populations were implanted from the same source, a contrast value of c > 1 (a divot) is preferred as this would simultaneously explain the Neptune Trojan and Plutino H-distributions. The Neptune Trojans must have been captured from the scattering TNOs and thus share the same H-distribution. The observed lack of small-sized Neptune Trojans (Sheppard & Trujillo 2010; Parker et al. 2013) requires that there be a sharp decrease in the number of Neptune Trojans below a certain size. Figure 3 of Sheppard & Trujillo (2010) shows that this decrease must be in the form of a drop in the differential number distribution of Neptune Trojans in order to explain the lack of detected Neptune Trojans within their magnitude limit. This drop can either be explained by a differential H-distribution of a transition to a steep negative slope, disfavored by our sample of scattering TNOs, or an H-distribution with c > 1. The preference of c > 1 from jointly considering the scattering TNOs and Neptune Trojans must also be consistent with the Plutino population if these three hot populations share a common source. Alexandersen et al. (2014) find a decrease in observed Plutinos at H magnitudes fainter (smaller) than our reported transitions in the scattering and Neptune Trojan H-distributions. Exploring a variety of H-distribution parameters, Alexandersen et al. (2014) find that the divot parameters favored by Shankman et al. (2013), with no tuning, provide a good representation of the observed H-distribution for 5 < H < 11 Plutinos. A divot (c > 1) solution simultaneously matches the observed scattering TNO, Neptune Trojan, Centaur, and Plutino H-distributions.

The Shankman et al. (2013) preferred H-distribution, with parameters of α_f = 0.5 and c ∼ 5.6, continues to provide an excellent match to our observations, even with the inclusion of two new surveys that add survey depth and double the observed sample. This H-distribution remains our preferred choice for its ability to match our observed Scattering TNO sample. Figure 9 shows the population estimates for our preferred divot in both cumulative and differential representations. We provide the cumulative representation to highlight how such a visually obvious feature is completely suppressed when viewed cumulatively. The flattening of the cumulative distribution seen in Figure 2 of Sheppard & Trujillo (2010) that then gradually transitions to a positive slope is characteristic of a divot. Figure 9 demonstrates that when you view a divot cumulatively, it produces exactly this signature. If there is a divot in the H-distribution of the TNO populations, prior studies of the luminosity function, or even H, could be suppressing the signature of this feature if they combine surveys which were not characterized in the same way.

Zoom In Zoom Out Reset image size

A divot can also explain the rollover to negative slopes that some groups have measured (Bernstein et al. 2004; Fraser & Kavelaars 2009). The transition in the size-frequency distribution for a divot would likely take the form of a steep negative slope over a short size range, rather than a discontinuous drop. If a survey terminates just past the break point, one could measure a signature of this transition as a negative slope without measuring the recovery.

Two different formation mechanisms have been shown to produce divot size-frequency distributions. Fraser (2009) shows that a kink in the size-frequency distribution collisionally propagates to larger sizes and produces a divot. Campo Bagatin & Benavidez (2012) explore the "born big" scenario (Morbidelli et al. 2009) in the Kuiper Belt, collisionally evolving a population with no initial objects smaller than 100 km. As the population evolves under collisions, a divot is produced and the small-sized objects have a collisional equilibrium slope and a contrast of the same scale (∼5; see Figure 6 of Campo Bagatin & Benavidez 2012) as our preferred model parameters (α_f = 0.5 and c ∼ 5.6). Divot size-frequency distributions are consistent with our understanding of small-body formation physics, and have been shown to arise from two different formation mechanisms.

Divot H-distributions are plausible, explain our observations and when generalized to other hot TNO populations can explain several outstanding questions. Our preferred H-distribution parameter set is one of the least rejectable H-distributions we tested, provides a good match with the observed H-distribution (see Figure 4), is consistent with the observed JFC slopes, can solve the problem of JFC supply, is consistent with the observed Centaur H-distribution, is consistent with the observed Neptune Trojan and Plutino H-distributions and is consistent with divots produced in collisional models. We find a divot of contrast c ∼ 5.6 and faint slope α_f to be compelling and argue that it provides a better, versus knees, solution to the observed lack of small-sized hot TNOs while simultaneously explaining the supply of JFCs.