THE 5:1 NEPTUNE RESONANCE AS PROBED BY CFEPS: DYNAMICS AND POPULATION

All astrometric measurements have inherent uncertainties that come from a variety of sources. Errors in the initial orbit solution and location prediction can make target identification and astrometric fits difficult. Stellar crowding, plate boundary errors, stellar motion, and sparse catalogs can make astrometric plate solutions less accurate. Because of these issues, systematic offsets for all observations in a single observing run are common. When a target has clusters of astrometry separated by months without observations, a systematic offset from one observing run can significantly bias the orbital fit (Gladman et al. 2008). Accounting for all sources of uncertainty in astrometric measurements can result in a wider range of orbital fit parameters than a simple orbital fit would suggest.

During our orbital fitting procedure, we examined all outliers in the astrometric fits to determine if the discrepant astrometry should be kept. With years of astrometry, much of it at monthly cadences, outliers are immediately obvious in the fit residuals calculated with the routine by Bernstein & Khushalani (2000). We revisited the original data for all high residual astrometric points, and re-examined the astrometric plate solution. For most of these observations, we calculated a new astrometric position. For several others, we found that a good astrometric solution could not be derived, which occurred as a result of low signal on the TNO (centroiding was impossible) or an insufficient number of (unsaturated) reference stars. We replaced or removed the astrometry in the orbital fit and repeated the check for outliers. Our careful re-measurement of discrepant points resulted in significantly improved astrometry for all of our objects.

Systematic astrometric offsets result in poor orbital predictions for TNOs observed once or twice a year for a limited period of time (Gladman et al. 2008). The orbit parameter estimation by Gladman et al. (2008) allows systematic offsets for some of the astrometric points and calculates the largest and smallest a clones possible from the astrometry. The larger orbit uncertainties occasionally reveal the possibility of resonance behavior for TNOs near a resonance boundary. These targets may have high a or low a clones which show a resonant angle libration. Gladman et al. (2008) noted this effect for L3y02, which at the time would not have been flagged as a potentially resonant object based on the Bernstein & Khushalani (2000) uncertainty range. The Gladman et al. (2008) classification is useful for targeting possible resonant objects with small or sparsely sampled arcs for additional astrometric followup, because these objects may have systematic offsets in their astrometry that are not apparent from the fit residuals.

Our combination of lower than average astrometric uncertainty and better observation cadence over long arcs allows us to develop our own estimation of the uncertainty in the objects' orbital parameters. After discrepant astrometry is removed from our results, the only sources of uncertainty in our data are the uncertainty in the position of reference stars and the centroid calculation of sources. Our astrometry is accurate to ∼0.2 arcsec. Instead of using the method by Gladman et al. (2008), we resampled each of our astrometric points linearly within their uncertainty of 0.2 arcsec . This resamples the astrometry of the objects within 1–1.5σ uncertainty range in R.A. and decl. We calculated a new orbital fit to the resampled astrometry for the object, so we have new orbital parameters for each resampled clone. Repeating this resampling process ∼1500 times provides the range of orbital parameters that are consistent with the astrometric measurements. Using multiple orbital fits instead of resampling the uncertainty on each orbital parameter from Bernstein & Khushalani (2000) also allows us to preserve the interdependency of the orbital parameters. The result was a cloud of clones of each TNO, shown in Figure 1. The behavior of the clones over time provides an exploration of the resonance behavior for orbits that are consistent with our knowledge of the intrinsic orbit of these four TNOs.

Zoom In Zoom Out Reset image size

We integrated the clones of the four potential 5:1 objects for approximately 4.5 Gyr using SWIFT (Levison & Duncan 1994) in order to determine resonance occupation for each initial condition. SWIFT integrates gravitationally interacting objects using a provided time-step, 0.5 years for this simulation. We used the Regularized Mixed Variable Symplectic method which handles test particle–planet close approaches. The position, velocity, and mass of the Sun, Jupiter, Saturn, Uranus, and Neptune were included as massive bodies, and their initial conditions were calculated for the epoch of the TNO's orbital fit. The clones of the possible 5:1 TNOs were added as massless test particles. We observed a variety of behaviors over the integration time, including short ($\sim {{10}^{5}}$ years) or long ($\gt {{10}^{7}}$ years) term resonance occupation in different 5:1 libration islands, scattering, capture into other Neptunian resonances, and stable non-resonant behavior. Because we are interested in the region near the 5:1 resonance, test particle information was recorded only for heliocentric distances between 20 and 150 AU. The positions of the clones that went outside these ranges were not recorded while they were beyond the limits. The behavior of all of the clones provides possible classifications of the TNOs found in our survey.

We determined whether an initial condition was resonant by examining the behavior of the resonant angle. Typically, secure resonant classification requires the resonant angle to librate for 10⁷ years from the start of integration for the best fit clone as well as the minimum and maximum a clones (Gladman et al. 2008). We classified objects as resonant if their resonant angle oscillated instead of exploring a full 360° circulation. Following the method of Alexandersen et al. (2013) for the 1:1 resonance, we used a moving window to examine the resonant angle behavior. Because many of these objects experienced multiple periods of resonant libration, the beginning and end of each resonant period was recorded. HL8k1 is never resonant, while L3y02, HL7j4, and HL7c1 all show periods of resonant behavior during the simulation.

Based on our simulations, HL7j4 is in the 5:1 resonance. The behavior of the best fit clone over the first 500,000 years of the simulation is shown in Figure 2. This symmetric oscillation was observed for all of the clones for a minimum of 2 × 10⁶ years, so we classify this object as a symmetric 5:1 librator. The median duration of the first period of resonance is $4.9\times {{10}^{6}}$ years. See Figure 1 for a plot of resonance duration. The higher occurrence of longer lived clones in the large a, e, and i corner suggest that this region may be more likely to contain the "true" orbit of the TNO. When the resonant periods for each clone are summed, more than 99% of the clones are resonant for at least a total of 10⁷ years of the entire simulation, and 80% are resonant for 10⁹ years. Of the clones that leave the 5:1 resonance, several are captured into other distant resonances, including the 6:1 and 7:1 resonance. Of our four objects, HL7j4 has the highest fraction of resonant clones that are also in the Kozai resonance (Kozai 1962; Lidov 1962), which can lead to larger variations in the test particles' eccentricities and inclinations. The clones of HL7j4 display the most stable resonant behavior of our test particles, and the "true" TNO may be on an orbit in a portion of the phase space stable for ∼Gyr.

Zoom In Zoom Out Reset image size

HL7c1 is also a 5:1 resonator. The best fit clone is resonant for the first $6\times {{10}^{7}}$ years, before exiting the resonance. The clone experiences two additional short resonance captures ($\sim {{10}^{5}}$ years) before scattering outward. The best fit clone passes beyond 150 AU at $3.2\times {{10}^{8}}$ years, then re-enters the region of interest at $7.1\times {{10}^{8}}$ years. It moves beyond 150 AU again briefly before scattering inward at $9.5\times {{10}^{8}}$ years where the simulation output ends. All clones of this object begin in symmetric 5:1 resonance, and remain in symmetric resonance for a minimum of 10⁶ years. The duration of the first resonance occupation for these clones is shown in Figure 1. The seemingly random distribution of resonance occupation times suggests that no particular aei region is favored. The median duration of the first period of resonance is $1.6\times {{10}^{7}}$ years. Approximately 90% of the clones of this object are resonant for a total of at least 10⁷ years. However, only 30% are resonant for a total of 10⁸ years and <10% for a total of 10⁹ years. After the clones leave the resonance, the majority move into scattering orbits and some become detached objects. Some of these scattering clones are captured into other resonances and some are ejected from the region of our SWIFT output (20–150 AU). HL7c1 is currently a 5:1 resonator.

L3y02 is also classified as resonant. The best fit clone begins with $5\times {{10}^{4}}$ years of resonance, exits the resonance, then re-enters the resonance at $2\times {{10}^{6}}$ years. Similar behavior is observed for the other clones; there are multiple periods of resonance broken by non-resonant behavior. All of these clones begin their first resonance as symmetric librators. The median duration of the first period of resonance is $6.0\times {{10}^{5}}$ years. Most of the resonance periods are of short duration; the length of resonance for each clone is shown in Figure 1. Over 95% of the clones are resonant for a total of at least 10⁷ years during the simulation, 80% for 10⁸ years, and 50% for 10⁹ years during the simulations. When leaving the resonance, these clones preferentially exit to low a, just sunward of the 5:1 resonance. Many clones remain there, becoming detached objects while the others enter the scattering population, some of which are captured into other resonances.

HL8k1 is not in the 5:1 Neptune resonance, in spite of being extremely close to the resonance boundary. All of its clones experience small periodic oscillations in a and e, shown in Figure 2, but not with the characteristic periodic pattern seen in the resonant object integrations. The resonant angle shows only a slow circulation instead of a bounded oscillation, so HL8k1 is classified as a detached object according to the Gladman et al. (2008) classification scheme. HL8k1 is not resonant, but it may have been associated with the 5:1 resonance at some point in the past.

The resonant angle behavior in Figure 2 for HL7j4 is representative of a typical symmetric libration behavior. Most of the clones of our three 5:1 objects are primarily symmetric librators, although a small number of clones transition into periods of asymmetric libration. The resonance angle ϕ librates with a high amplitude centered around 180°. This large amplitude symmetric libration was typical for ∼90% of the resonant clones; the smaller amplitude asymmetric libration was far less common and typically occurred between periods of large amplitude symmetric libration. The intrinsic asymmetric fraction of the resonance could be as high as ∼30% if the fraction is similar to the n:1 resonance symmetric fractions from Gladman et al. (2012). Our smaller observed fraction may be the result of the longitudinal biases of the discovery surveys. Understanding the distribution of the libration amplitudes and libration centers will be critical to understanding the 5:1 resonance capture process.

The detected 5:1 TNOs all have large eccentricity ($\gt 0.54$), but ∼5% of their resonant clones were found to undergo a decrease in eccentricity as a result of resonant dynamics. The eccentricity of resonant clones ranged from approximately 0.65 to 0.30, with the majority of clones having values between 0.5 and 0.65. After 10⁹ years 3% of resonant clones had cycled to an e below 0.5, and after 2 × 10⁹ years 5% of the remaining resonant clones had cycled to an $e\lt 0.5$. This cycling is visible in Figures 3 and 4. This suggests that a Kozai mechanism (Kozai 1962; Lidov 1962), where there is an exchange between eccentricity and inclination, is acting on the clones. The evolution of e and i results in objects with inclinations up to 45°. These objects remained primarily resonant, with low e and larger q than the detected objects in the 5:1 resonance. The greater stability of low e objects in the 5:1 resonance, combined with the supply of this population from large e objects, implies that there is likely a population in the 5:1 resonance stored at low e. Many of the 5:1 resonators that are also in Kozai resonance will cycle back to a larger e over time, so the fraction of objects at low e remains relatively constant after 10⁹ years. Current survey depths make it nearly impossible to detect all but the brightest low e 5:1 objects; for CFEPS with limiting magnitude of $g\;\sim$ 23.5–24.4 (Petit et al. 2011), a 5:1 TNO with $e=0.3$ (q ∼ 61 AU) would have to have an absolute magnitude ${{H}_{g}}\;\sim$ 5.6–6.5 in order to be detected at perihelion. We did not detect any low e objects, so to avoid complication from this possibly interesting component we ignore this ∼5% component of low e objects that may have evolved from larger e objects when we calculate our population estimate in Section 4.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The clones of the resonant objects that left the resonance with e values such that they do not immediately scatter off Neptune were preferentially deposited or retained at a and e values similar to the current orbit of HL8k1, the non-resonant object. This overabundance of clones is shown in Figure 4. These objects are found along the resonance borders and at a slightly smaller a and e than the 5:1 resonance. HL8k1 exhibits stable behavior over the age of the Solar System. Figure 2 shows the first $5\times {{10}^{5}}$ years of HL8k1's evolution, and the behavior is representative of the entire 4.5 Gyr history of all clones of HL8k1. The clones experience regular, but small, changes in a and e as well as slow evolution of the resonance angle. None of the clones of HL8k1 are captured by the resonance during their 4.5 Gyr integrations, which reflects the low probability of reentering the resonance. This object could have been emplaced in its current location during Neptune's migration (Gomes et al. 2008) or diffused out of the resonance without the aid of planetary migration, which was not included in our simulations. We use "resonant diffusion object" to refer to a TNO that appears to have leaked out of a resonance and entered a pseudo-stable region of phase space. The near-resonant object, HL8k1, was likely dynamically associated with the 5:1 resonance and may be representative of a large number of resonant diffusion objects.

Our model of the 5:1 resonance is based on the regions of dynamical stability. The structure of the n:1 resonances is complex, due to the presence of multiple libration islands. The symmetric resonance objects come to pericenter 180° from Neptune (plus or minus their libration amplitude which can approach 180°), and the objects in the two asymmetric resonant islands come to pericenter on each side, 120° from Neptune (plus or minus a smaller libration amplitude). These complexities make visualization difficult, but a simplified model of the 5:1 resonance is presented in Gladman et al. (2012). Here we use a similar parameterization of the resonance structure.

The semimajor axis range of a resonance is restricted by the resonance width. The exact width depends on the other orbital parameters (such as inclination and eccentricity), but for the purposes of detectability it is sufficient to establish a model with a single range of a values. While the density of points at the a resonance edges in the left plot in Figure 2 make it clear that the clones spend the majority of their time near the resonance edges, the resonance is narrow enough that this complication does not have a significant impact on the detectability of model objects. We model the a range as a simple uniform distribution between 87.9 and 88.9 AU.

The objects in the 5:1 resonance appear to have a very excited inclination distribution. The detected objects all have high inclinations, even L3y02, which was discovered in the ecliptic block of CFEPS. The ecliptic blocks of CFEPS were sufficiently deep that 5:1 objects with low i would have been detectable if they were present. The detection of a high inclination object near the ecliptic, where high i objects spend only a small fraction of their orbits (Gladman et al. 2012), as well as two object discoveries in the high latitude blocks indicates that this population is dynamically excited.

Based on our model analysis, we find that the 5:1 resonance probably has a wider inclination distribution than the Plutinos. We implemented an inclination distribution parameterization from Brown (2001) which is acceptable for other studied resonances (Gladman et al. 2012):

$$\begin{eqnarray}&&P(i)\propto {\rm sin} (i)\times {\rm exp} \left( \frac{-{{i}^{2}}}{2\sigma _{i}^{2}} \right).\end{eqnarray} \tag{ 1 }$$

We drew inclinations randomly from this inclination distribution for our 5:1 population model, using ${{5}^{\circ }}\lt {{\sigma }_{i}}\lt 65{}^\circ$ in increments of 1°. 10,000 simulated detections were produced using the CFEPS survey simulator, and then the inclination of these simulated detections was compared to the three measured object detections. We used the Anderson–Darling statistical test with a bootstrapped sample from the simulations to determine the likelihood of reproducing the three measured inclinations, shown in Figure 5. HL8k1 has the largest inclination of the discovered objects, so if the fourth object were included, the inclination width would be higher. The most acceptable inclination width, ${{\sigma }_{i}}$, for the three 5:1 objects was 22°. We reject inclination widths less than 15° and more than 43° at 95% confidence. The published inclination widths for the Plutinos are $14{^\circ} _{-4}^{+8}$ (Alexandersen et al. 2015, 95% confidence limits), and $16{^\circ} _{-4}^{+8}$ (Gladman et al. 2012, 95% confidence limits), and $10\buildrel{\circ}\over{.} 7_{-2.3}^{+2.0}$ (Gulbis et al. 2010, 68% confidence limits). We find that the 5:1 resonance likely has a wider inclination distribution than the Plutinos and other hot TNO populations, however we cannot rule out some similar inclination widths at 95% confidence.

Zoom In Zoom Out Reset image size

It is possible that a significant amount of objects in the n:1 resonances are resonance sticking from the scattering object population, either currently or in a primordial scattering event. Resonance sticking objects may not include a distribution that extends to low inclinations, in which case the Brown (2001) distribution may overestimate the size of the resonance population. Our sample of three objects is too small to rule out a low i tail, but simulations of planetary migration resulting in resonance sticking show an inclination distribution not centered at zero (Gomes et al. 2005). Our clones show a similar evolution to the particles Gomes et al. (2005), including the influence of the Kozai mechanism on the inclination distribution. We explore a toy model of the inclination distribution based on the i distribution presented in Figure 3 of Gomes et al. (2005); a Gaussian peak center, μ, at 35° instead of 0° with a width ${{\sigma }_{i}}$ of 7°. This conservative parameterization removes the undetected low i tail present in the Brown (2001) model.

$$\begin{eqnarray}&&P(i)\propto {\rm sin} (i)\times {\rm exp} \left( \frac{-{{(i-\mu )}^{2}}}{2\sigma _{i}^{2}} \right).\end{eqnarray} \tag{ 2 }$$

The apparent eccentricity distribution of a population at large semimajor axis is strongly biased. Objects from the 5:1 resonance are typically only sufficiently bright for detection near pericenter, and if their pericenter is too large they may not be detectable even at pericenter. For this reason, we did not assume an eccentricity range for our model population based on another resonant TNO population but instead used an eccentricity range based on the range in our real detections and the bulk of the integrated clones (see Figure 4). We used a uniform eccentricity distribution from 0.50 to 0.62, or pericenters of ∼30 to 45 AU. The known objects have eccentricities in the range 0.55–0.62. While there may be a component of the 5:1 resonance with lower eccentricity due to the Kozai mechanism (suggested at the 5% level by our simulation results) or some other mechanism not present in our simulations, we present a minimum population estimate by not extending the model to eccentricities significantly lower than our detected TNOs.

Our cloned objects explore much of the resonance space, so we used the resonant libration characteristics of the clones to constrain the libration centers and amplitudes in our parametric model. The population model assumes 90% symmetric and 10% asymmetric resonators. We did not detect any asymmetric resonators and do not place strong constraints on their population fraction. We also calculated the range of libration amplitudes explored by the 5:1 resonators and included that in our orbital distribution model. The symmetric objects had their expected libration centers at 180° from Neptune, and an amplitude of $175{^\circ} _{-5}^{+3}$. The trailing island was centered at 240° with an amplitude of $70{^\circ} _{-40}^{+20}$, and the leading island was centered at 120° with an amplitude of $70{^\circ} _{-40}^{+20}$. We calculated positions for each model object by drawing from a triangular distribution of the libration amplitude encompassing the libration amplitude range reported. An angle was randomly selected from within the object's libration range, instead of sinusoidally weighting the objects position. This gives a closer approximation to the sawtooth pattern of the resonant angle seen in Figure 2. This angle was added to its libration center to give the current resonant angle, ${{\phi }_{51}}$, for the model object. As described above, the resonance characteristics of the integrated clones were used to generate a representative sample of resonance angles.

The three additional orbital parameters (argument of pericenter, mean anomaly, and longitude of the ascending node) must fulfill the resonance condition. A detailed discussion of the resonance angles can be found in Murray-Clay & Chiang (2005). We selected the mean anomaly ($\mathcal{M}$) and the node (Ω) randomly, and calculated the necessary argument of pericenter (ω) using the resonant angle, ${{\phi }_{51}}$. This condition required Neptune's position, ${{\lambda }_{N}}$ = ${{\omega }_{N}}+{{{\Omega }}_{N}}+{{\mathcal{M}}_{N}}$, at the selected epoch. The resonance condition is:

$$\begin{eqnarray}&&{{\phi }_{51}}=5\lambda -{{\lambda }_{N}}-4\varpi \end{eqnarray} \tag{ 3 }$$

where ϖ = $\omega +{\Omega }$. This means that ω was calculated via:

$$\begin{eqnarray}&&\omega ={{\phi }_{51}}-{\Omega }-5M+{{\lambda }_{N}}.\end{eqnarray} \tag{ 4 }$$

An absolute H magnitude distribution is needed in order to estimate the number of objects in the resonance. This distribution determines the number of objects in the population of each size. Our three 5:1 resonant objects have H_g of 6.0, 7.3, and 7.7 (see Table 2). Our population estimate only extends to ${{H}_{g}}=8$, providing a more accurate population estimate for this sparsely sampled resonance than extrapolating to smaller sizes. ${{H}_{g}}=8$ corresponds to an approximate diameter of 170 km (Petit et al. 2011). A single power law with a slope of $\alpha =0.8$ has been shown to be appropriate for the scattering objects to ${{H}_{g}}\sim 9$ (Shankman et al. 2013). A similar result for ${{H}_{B}}\lt 7.7$ was reported by Fraser et al. (2014), who found the "hot" population had a bright end slope of $\alpha =0.87$. We expect our objects to be similar to scattering objects, so we used a single power law with a slope of $\alpha =0.8$ for our size distribution, and only estimate the population for ${{H}_{g}}\lt 8$. Because we only extended our population estimate over the range of our detected objects, a slightly steeper slope, as in Fraser et al. (2014), does not significantly impact our population estimate.

With the assumptions described in the previous section, we made a population model for the 5:1 resonant objects. We used inclination width σ = 22°, the eccentricity range of 0.50–0.62, and symmetric fractions and libration amplitudes consistent with the clones, discussed above. This results in a minimum model-based population estimate for the 5:1 resonance.

We estimated the size of the underlying population of the resonance by using a survey simulator (Jones et al. 2006; Kavelaars et al. 2009). The pointings and depths of each survey field are provided to the simulator from the CFEPS ecliptic and high latitude blocks and the Alexandersen et al. (2015) survey. Analysis similar to Kavelaars et al. (2009) showed consistency between the real object detection and the model objects, once biased by the survey simulator. The survey simulator works by testing the detectability of each model object against the survey fields. Synthetic objects were drawn from the model until three tracked TNOs were "discovered" by the survey simulator; the total number of objects drawn gives a possible population size. This process was repeated 6000 times, and a range of population estimates were computed, shown in Figure 6.

Zoom In Zoom Out Reset image size

We find that the 5:1 resonance contains $1900_{-1400}^{+3300}$ TNOs with ${{H}_{g}}\;\lt$ 8 and $e\gt 0.5$. This estimate is approximately the size of the Plutino population (Gladman et al. 2012) and possibly larger, shown in Figure 6. Our result is consistent with the 5:1 population from Gladman et al. (2012) based on a single detection; their population model assumed a narrower inclination distribution but extended the model to $0.35\lt e\lt 0.65$. For our toy model with the inclination distribution drawn from a Gaussian centered at 35°, we find $2400_{-1800}^{+4000}$ objects in the 5:1 resonance with ${{H}_{g}}\lt 8$ and $e\gt 0.5.$ We tested some additional model parameters to determine the extent of the model dependence of the population estimate. Changing the symmetric fraction from 10% to 30% resulted in a population estimate within 1% of our nominal model. Using the inclination distribution extremes (14° and 44°) gave population estimates 5% larger and 17% smaller, respectively. This outer resonance represents a significant reservoir of TNOs, regardless of the inclination distribution of the population.