# A Search of the Full Six Years of the Dark Energy Survey for Outer Solar System Objects

Pedro H. Bernardinelli, Gary M. Bernstein, Masao Sako, Brian Yanny, M.

Aguena, S. Allam, F. Andrade-Oliveira, E. Bertin, D. Brooks, E.

Buckley-Geer, et al.

► To cite this version:

Pedro H. Bernardinelli, Gary M. Bernstein, Masao Sako, Brian Yanny, M. Aguena, et al. A Search of the Full Six Years of the Dark Energy Survey for Outer Solar System Objects. The Astrophysical Journal Supplement Series, 2022, 258, 10.3847/1538-4365/ac3914. insu-03748297

HAL Id: insu-03748297

https://insu.hal.science/insu-03748297

Submitted on 9 Aug 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# A Search of the Full Six Years of the Dark Energy Survey for Outer Solar System Objects

Pedro H. Bernardinelli¹, Gary M. Bernstein¹, Masao Sako¹, Brian Yanny², M. Aguena³, S. Allam², F. Andrade-Oliveira³,⁴, E. Bertin⁵,⁶, D. Brooks⁷, E. Buckley-Geer²,⁸, D. L. Burke⁹,¹⁰, A. Carnero Rosell³, M. Carrasco Kind¹¹,¹², J. Carretero¹³, C. Conselice¹⁴,¹⁵, M. Costanzi¹⁶,¹⁷,¹⁸, L. N. da Costa³,¹⁹, J. De Vicente²⁰, S. Desai²¹, H. T. Diehl², J. P. Dietrich²², P. Doel⁷, K. Eckert¹, S. Everett²³, I. Ferrero²⁴, B. Flaugher², P. Fosalba²⁵,²⁶, J. Frieman²,²⁷, J. García-Bellido²⁸, D. W. Gerdes²⁹,³⁰, D. Gruen²², R. A. Gruendl¹¹,¹², J. Gschwend³,¹⁹, S. R. Hinton³¹, D. L. Hollowood²³, K. Honscheid³²,³³, D. J. James³⁴, S. Kent²,²⁷, K. Kuehn³⁵,³⁶, N. Kuropatkin², O. Lahav⁷, M. A. G. Maia³,¹⁹, M. March¹, F. Menanteau¹¹,¹², R. Miquel¹³,³⁷, R. Morgan³⁸, J. Myles⁹,¹⁰,³⁹, R. L. C. Ogando³,¹⁹, A. Palmese²,²⁷, F. Paz-Chinchón¹¹,⁴⁰, A. Pieres³,¹⁹, A. A. Plazas Malagón⁴¹, A. K. Romer⁴², A. Roodman⁹,¹⁰, E. Sanchez²⁰, V. Scarpine², M. Schubnell³⁰, S. Serrano²⁵,²⁶, I. Sevilla-Noarbe²⁰, M. Smith⁴³, M. Soares-Santos³⁰, E. Suchyta⁴⁴, M. E. C. Swanson¹¹, G. Tarle³⁰, C. To⁹,¹⁰,³⁹, T. N. Varga⁴⁵,⁴⁶, and A. R. Walker⁴⁷

(The DES Collaboration)

¹ Department of Physics and Astronomy, University of Pennsylvania, Philadelphia, PA 19104, USA; pedrobed@sas.upenn.edu

² Fermi National Accelerator Laboratory, P.O. Box 500, Batavia, IL 60510, USA

³ Laboratório Interinstitucional de e-Astronomia—LineA, Rua Gal. José Cristino 77, Rio de Janeiro, RJ—20921-400, Brazil

⁴ Instituto de Física Teórica, Universidade Estadual Paulista, São Paulo, Brazil

⁵ CNRS, UMR 7095, Institut d'Astrophysique de Paris, F-75014, Paris, France

⁶ Sorbonne Universités, UPMC Univ Paris 06, UMR 7095, Institut d'Astrophysique de Paris, F-75014, Paris, France

⁷ Department of Physics & Astronomy, University College London, Gower Street, London, WC1E 6BT, UK

⁸ Department of Astronomy and Astrophysics, University of Chicago, Chicago, IL 60637, USA

⁹ Kavli Institute for Particle Astrophysics & Cosmology, P.O. Box 2450, Stanford University, Stanford, CA 94305, USA

¹⁰ SLAC National Accelerator Laboratory, Menlo Park, CA 94025, USA

¹¹ Center for Astrophysical Surveys, National Center for Supercomputing Applications, 1205 West Clark Street, Urbana, IL 61801, USA

¹² Department of Astronomy, University of Illinois at Urbana-Champaign, 1002 West Green Street, Urbana, IL 61801, USA

¹³ Institut de Física d'Altes Energies (IFAE), The Barcelona Institute of Science and Technology, Campus UAB, E-08193 Bellaterra (Barcelona) Spain

¹⁴ Jodrell Bank Center for Astrophysics, School of Physics and Astronomy, University of Manchester, Oxford Road, Manchester, M13 9PL, UK

¹⁵ University of Nottingham, School of Physics and Astronomy, Nottingham NG7 2RD, UK

¹⁶ Astronomy Unit, Department of Physics, University of Trieste, via Tiepolo 11, I-34131 Trieste, Italy

¹⁷ INAF-Osservatorio Astronomico di Trieste, via G.B. Tiepolo 11, I-34143 Trieste, Italy

¹⁸ Institute for Fundamental Physics of the Universe, Via Beirut 2, I-34014 Trieste, Italy

¹⁹ Observatório Nacional, Rua Gal. José Cristino 77, Rio de Janeiro, RJ—20921-400, Brazil

²⁰ Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas (CIEMAT), Madrid, Spain

²¹ Department of Physics, IIT Hyderabad, Kandi, Telangana 502285, India

²² Faculty of Physics, Ludwig-Maximilians-Universität, Scheinerstr. 1, D-81679 Munich, Germany

²³ Santa Cruz Institute for Particle Physics, Santa Cruz, CA 95064, USA

²⁴ Institute of Theoretical Astrophysics, University of Oslo, P.O. Box 1029, Blindern, NO-0315 Oslo, Norway

²⁵ Institut d'Estudis Espacials de Catalunya (IEEC), E-08034 Barcelona, Spain

²⁶ Institute of Space Sciences (ICE, CSIC), Campus UAB, Carrer de Can Magrans, s/n, E-08193 Barcelona, Spain

²⁷ Kavli Institute for Cosmological Physics, University of Chicago, Chicago, IL 60637, USA

²⁸ Instituto de Física Teórica UAM/CSIC, Universidad Autónoma de Madrid, E-28049 Madrid, Spain

²⁹ Department of Astronomy, University of Michigan, Ann Arbor, MI 48109, USA

³⁰ Department of Physics, University of Michigan, Ann Arbor, MI 48109, USA

³¹ School of Mathematics and Physics, University of Queensland, Brisbane, QLD 4072, Australia

³² Center for Cosmology and Astro-Particle Physics, The Ohio State University, Columbus, OH 43210, USA

³³ Department of Physics, The Ohio State University, Columbus, OH 43210, USA

³⁴ Center for Astrophysics | Harvard & Smithsonian, 60 Garden Street, Cambridge, MA 02138, USA

³⁵ Australian Astronomical Optics, Macquarie University, North Ryde, NSW 2113, Australia

³⁶ Lowell Observatory, 1400 Mars Hill Road, Flagstaff, AZ 86001, USA

³⁷ Institució Catalana de Recerca i Estudis Avançats, E-08010 Barcelona, Spain

³⁸ Physics Department, 2320 Chamberlin Hall, University of Wisconsin-Madison, 1150 University Avenue, Madison, WI 53706-1390, USA

³⁹ Department of Physics, Stanford University, 382 Via Pueblo Mall, Stanford, CA 94305, USA

⁴⁰ Institute of Astronomy, University of Cambridge, Madingley Road, Cambridge CB3 0HA, UK

⁴¹ Department of Astrophysical Sciences, Princeton University, Peyton Hall, Princeton, NJ 08544, USA

⁴² Department of Physics and Astronomy, Pevensey Building, University of Sussex, Brighton, BN1 9QH, UK

⁴³ School of Physics and Astronomy, University of Southampton, Southampton, SO17 1BJ, UK

⁴⁴ Computer Science and Mathematics Division, Oak Ridge National Laboratory, Oak Ridge, TN 37831, USA

⁴⁵ Max Planck Institute for Extraterrestrial Physics, Giessenbachstrasse, D-85748 Garching, Germany

⁴⁶ Universitäts-Sternwarte, Fakultät für Physik, Ludwig-Maximilians Universität München, Scheinerstr. 1, D-81679 München, Germany

⁴⁷ Cerro Tololo Inter-American Observatory, NSFs NOIRLab, Casilla 603, La Serena, Chile

Received 2021 September 10; revised 2021 November 2; accepted 2021 November 10; published 2022 February 9

<figure><img src="image_2.png" /></figure>
Abstract

We present a search for outer solar system objects in the 6 yr of data from the Dark Energy Survey (DES). The DES covered a contiguous 5000 deg² of the southern sky with ≈80,000 3 deg² exposures in the *grizY* filters between 2013 and 2019. This search yielded 812 trans-Neptunian objects (TNOs), one Centaur and one Oort cloud comet, 458 reported here for the first time. We present methodology that builds upon our previous search on the first 4 yr of data. All images were reprocessed with an optimized detection pipeline that leads to an average completeness gain of 0.47 mag per exposure, as well as improved transient catalog production and algorithms for linkage of detections into orbits. All objects were verified by visual inspection and by the “sub-threshold significance,” the signal-to-noise ratio in the stack of images in which its presence is indicated by the orbit, but no detection was reported. This yields a pure catalog complete to r ≈ 23.8 mag and distances 29 < *d* < 2500 au. The TNOs have minimum (median) of 7 (12) nights’ detections and arcs of 1.1 (4.2) yr, and will have *grizY* magnitudes available in a further publication. We present software for simulating our observational biases for comparisons of models to our detections. Initial inferences demonstrating the catalog’s statistical power are: the data are inconsistent with the CFEPS-L7 model for the classical Kuiper Belt; the 16 “extreme” TNOs (a > 150 au, q > 30 au) are consistent with the null hypothesis of azimuthal isotropy; and nonresonant TNOs with q > 38 au, a > 50 au show a significant tendency to be sunward of major mean-motion resonances.

Unified Astronomy Thesaurus concepts: Small Solar System bodies (1469); Trans-Neptunian objects (1705); Comets (280); Astronomy data analysis (1858); Kuiper belt (893)

Supporting material: FITS file

**1. Introduction**

The population of small bodies orbiting beyond Neptune is a
remnant of events early in the formation of the solar system.
The current orbital distribution of these trans-Neptunian objects
(TNOs) is the result of the migration of the giant planets
(Fernández & Ip 1984; Tsiganis et al. 2005; Levison et al.
2008) and, since the discovery of the second Kuiper Belt object
by Jewitt & Luu (1993), numerous subsequent surveys of the
trans-Neptunian region have identified thousands of objects
(e.g., Jewitt & Luu 1995; Gladman et al. 1998; Allen et al.
2001, 2002; Trujillo et al. 2001; Millis et al. 2002; Bernstein
et al. 2004; Elliot et al. 2005; Fuentes & Holman 2008; Fraser
et al. 2010; Schwamb et al. 2010; Petit et al. 2011, 2017;
Rabinowitz et al. 2012; Brown et al. 2015; Bannister et al.
2016, 2018; Sheppard et al. 2016, 2019; Weryk et al. 2016;
Chen et al. 2018; Whidden et al. 2019; Bernardinelli et al.
2020a).<fnref n="48" /> The observed variety of dynamical classes (Gladman
et al. 2008) and surface compositions (Brown 2012) has led to
different hypotheses about the formation of this region.
Neptune’s migration can trap planetesimals into mean-motion
resonances (Malhotra 1993, 1995), and gravitational interac-
tions between Neptune and these planetesimals can further
excite their orbits (Gomes 2003). More detailed models of the
formation of this region include instabilities in Neptune’s
orbit (Batygin et al. 2012; Dawson & Murray-Clay 2012),
variations in Neptune’s migration timescale and smoothness
(Nesvorný 2015; Pike et al. 2017), the effects of a potential
unobserved giant planet in the outer solar system (Batygin et al.
2019), close stellar encounters (Pfalzner et al. 2018), the birth
environment of the solar system (Adams 2010), and Galactic
tides (Duncan et al. 2008).

The Dark Energy Survey (DES; The Dark Energy Survey Collaboration 2005, 2016) received an allocation of 575 observing nights on the 4 m Blanco Telescope in Cerro Tololo between 2013 and 2019, with the primary objective of measuring the nature of the accelerated expansion of the universe and the

spatial distribution of dark matter (The Dark Energy Survey Collaboration 2018a, 2018b, 2019a, 2019b, 2019c). The 3 deg², 520 Mpix Dark Energy Camera (DECam; Flaugher et al. 2015) was built for the survey, enabling the “wide” component of the survey to image a contiguous 5000 deg² area of the southern sky ≈10 times in each of the *grizY* bands over 6 yr. The interlaced “supernova” or “deep” survey imaged a smaller area (≈30 deg²) in *griz* at a ≈weekly cadence, aimed at detecting and characterizing Type Ia supernovae (Bernstein et al. 2012). DES has reported discoveries of individual TNOs of interest (Gerdes et al. 2016, 2017; Becker et al. 2018; Khain et al. 2018; Lin et al. 2019), as well as dynamical classifications for many detected objects (The Dark Energy Survey Collaboration 2016; Khain et al. 2020; Bernardinelli et al. 2020a), statistical analyses of the population of large-*a*, large-*q* “extreme” TNOs (Bernardinelli et al. 2020b; Napier et al. 2021), forecasts of future occultation events (Banda-Huarca et al. 2019), and a survey of machine-learning techniques for TNO searches (Henghes et al. 2021).

Bernardinelli et al. (2020b, hereafter Paper I) presented the
results of a uniform TNO search in the first four years of DES
wide-survey data (“Y4”) that yielded 316 objects. The Y4
search required the development of several algorithms for
moving object identification, orbit linking and recovery,
object confirmation, and simulations of object discovery and
observational biases (see also Bernardinelli et al. 2020b). In
this paper, we describe a search of the full six years of DES
wide-survey data (“Y6”) for TNOs. Methodological improve-
ments over the Y4 search are described in Section 2.
Section 3 characterizes the selection function for the Y6 DES
TNO search, and describes tools that we make available for
simulating DES discoveries given hypothetical TNO popula-
tions. These tools are required for statistical comparisons of
theoretical TNO populations to the Y6 catalog of DES TNOs.
We demonstrate that the detection efficiency for bound
objects at 30–2000 au distance is almost entirely dependent
on the mean *r*-band apparent magnitude of the TNO, and
independent of light-curve amplitude, color, or orbital
elements, aside from the considerations of having to reside
within the DES footprint for a sufficient fraction of the survey
epoch.

<footnote n="48">⁴⁸ Bannister (2020) presents a comprehensive review of these surveys and summarizes their discoveries.</footnote>
<table><caption>**Table 1**
Stages of the Y6 Processing</caption><thead><tr><th>Catalog/Processing Step</th><th>No. of Real Elements</th><th>No. of Injected Elements</th><th>Changes in Y6</th></tr></thead><tbody><tr><td>Single-epoch detections (Section 2.1)</td><td>1.60 × 10<sup>10</sup></td><td>1.10 × 10<sup>5</sup></td><td>Lower threshold</td></tr><tr><td>→Blinded fake injection (Section 2.3) →</td><td></td><td></td><td>Changed order</td></tr><tr><td>→Transient identification (Section 2.4) →</td><td></td><td></td><td>Changed order</td></tr><tr><td>→ Coadd avoidance (Section 2.5) →</td><td></td><td></td><td>New</td></tr><tr><td>→Pixel-level masking (Section 2.6) →</td><td></td><td></td><td>New</td></tr><tr><td>Transients (Section 2.7)</td><td>1.08 × 10<sup>8</sup></td><td>1.05 × 10<sup>5</sup></td><td></td></tr><tr><td>→Pair finding (Section 2.8) →</td><td></td><td></td><td></td></tr><tr><td>Pairs</td><td>≥ 10<sup>13</sup></td><td>...</td><td></td></tr><tr><td>→Triplet finding (Appendix) →</td><td></td><td></td><td>Modified</td></tr><tr><td>Triplets</td><td>≈10<sup>12</sup></td><td>...</td><td></td></tr><tr><td>→Orbit growing (Section 2.8) →</td><td></td><td></td><td></td></tr><tr><td>→Fake unblinding →</td><td></td><td></td><td></td></tr><tr><td>Sevenlets</td><td>31064</td><td></td><td></td></tr><tr><td>→Reliability cuts (Section 2.8) →</td><td></td><td></td><td></td></tr><tr><td>Candidates</td><td>9081</td><td>3937</td><td></td></tr><tr><td>→Sub-threshold significance test (Section 2.9) →</td><td></td><td></td><td>Lower threshold</td></tr><tr><td>Visual inspection</td><td>872</td><td>...</td><td></td></tr><tr><td>Confirmed objects</td><td>814</td><td>...</td><td></td></tr></tbody></table>

**Note.** The entries in *italics* indicate each intermediate step, going from the SE catalogs to the final catalog of TNOs. The columns also indicate the number of elements and the number of injected elements (when appropriate) at each stage of the processing. The final column describes changes, if present, between this step and the Y4 processing. The search required between 15 and 20 million CPU hours to be completed.

Section 4 presents the catalog of TNOs detected in the full DES wide-survey data. This catalog has 814 confirmed objects (458 first discovered in this work) with *grizY* photometry, high-quality multiyear orbital solutions yielding precise dynamical classifications, and well-characterized observational biases. This is the second largest TNO catalog from a single survey to date, as well as the largest catalog with multiband photometry. Section 5 discusses statistical comparisons between the DES objects and models of the trans-Neptunian region. Section 6 summarizes this release. The Appendix describes a new triplet search algorithm developed for this search.

## 2. The Year 6 Search

The search for TNOs in the DES data has three major stages:

1. The identification of single-night transients, that is, a source that appears in the sky in a given location for only a single night;

2. Orbit linking, where we “connect” all sets of transients that could potentially come from the same outer-solar-system object;

3. Verification of the linkage through the “sub-threshold significance” (STS) statistic, whereby we stack all images containing the putative TNO that did *not* yield a detection, thus yielding an independent check of whether the TNO actually exists.

These steps are summarized below, with more detailed descriptions of any aspects that have changed since the Y4 search. Table 1 outlines the search, and presents the number of sources relevant to each step of the processing.

### 2.1. Data Acquisition and Image Processing

The search presented here uses the DES Y6A2 internal release. Of 83,706 exposures processed from the wide survey, 76,217 pass quality cuts, and their detections are entered into

the “single-epoch” (SE) catalog. “Coadd” images are created as the average (per band) of all exposures in a given region of sky, and the detections and measurements derived from these images comprise the coadd catalog. The DES image-processing pipeline is described in detail in Morganson et al. (2018), and the coadd catalogs as well as images used in this work correspond to those in the DES Data Release 2 (The Dark Energy Survey Collaboration 2021).

The nominal survey strategy (Diehl et al. 2016, 2018) was such that the 5000 deg² footprint was tiled with 10 × 90 s exposures in each of the *griz* bands (and 6 × 45 + 2 × 90 s exposures in the *Y* band). Each point within the footprint has been imaged by working detector pixels in 7–10 exposures per band, typically eight (The Dark Energy Survey Collaboration 2021). The observation scheduler (Neilsen & Annis 2013; Neilsen et al. 2019) specifically avoids repeated exposures in the same region and in the same night, except for successive exposures with the same pointing in different bands, so the motion of a TNO is not readily detectable in our data. We have elected to remove the *Y*-band catalogs from the search, since most TNOs have very low signal-to-noise ratios (S/Ns) in this filter. We do, however, use these images for *Y*-band photometric measurements of any TNO discovered in the *griz* data, and for additional astrometric measurements of bright objects that have *Y*-band detections. This means that the catalog in which the search is conducted includes all 63,853 wide-field *griz* exposures, and we note that the linking process (Section 2.8) does not use the information of which filter was used for each exposure. Multiband *grizY* magnitudes are available for all of the objects reported herein, and these will be presented in a future publication.

Each DES exposure has its astrometric solution mapped to Gaia DR2 (Gaia Collaboration et al. 2018). The DECam astrometric model is fully described in Bernstein et al. (2017a), and all fixed optical distortions are known to ≈1 mas rms. This astrometric model includes two color-dependent effects— differential chromatic refraction in the atmosphere, and lateral
<figure><img src="image_3.png" /><figcaption>**Figure 1.** Left: distribution of single-exposure detection limit $m_{50}$, as defined in Equation (1), for all of the Y6 exposures in each band. Right: histogram of the difference in $m_{50}$ for between the year 4 and year 6 processing for each exposure. The mean gain in depth is ~0.46 for the r band. This gain comes from primarily from the optimized detection threshold, and to a lesser extent from the deeper DES Y6 coadd catalogs used to veto static sources.</figcaption></figure>

<figure><img src="image_4.png" /></figure>

color in the DECam prime-focus corrector. For the TNO search, we use astrometric positions that assume a nominal mean TNO color for all detections. Once the TNOs are linked and validated, we can determine a g-i color for each. We then recalculate the positions using the measured color, and refit the orbit. The dominant astrometric error for the brightest TNO (Eris) is from turbulent atmospheric distortions, with a typical rms value of 7 mas (Fortino et al. 2021). All other TNOs have astrometric accuracy limited by shot noise in the centroiding. The estimations of these measurement uncertainties are shown to be accurate to <10% by the $\chi$ extsuperscript2 values per degree of freedom of the orbital fits for Eris (62.9/58) and for the entire final catalog (28,832/25,236); that is, there is no overfitting in the data.

The DES detections are calibrated to a photometric system that is highly uniform across the focal plane (Bernstein et al. 2017b) and across all exposures (Burke et al. 2017), as demonstrated by an rms difference from the Gaia DR2 catalog of only 2.5 mmag (The Dark Energy Survey Collaboration 2021). The photometric measurements reported for our TNOs are obtained using a dedicated photometric pipeline, whose details will be presented in a future publication. We report the average point-spread function (PSF)-fitting flux over all images of a TNO in a given band, with the flux of each individual detection scaled to a nominal distance for that object.

## 2.2. Detection Threshold Optimization and Characterization

We have improved the image-processing pipeline presented in Morganson et al. (2018) to allow for fainter sources to be detected in each exposure. This is done by lowering the SOURCE EXTRACTOR (Bertin & Arnouts 1996) detection threshold to `DETECT_THRESHOLD = 0.8` instead of the previous 1.5 used in the DES First Cut (Morganson et al. 2018) processing, and using a detection kernel that better matches the typical PSF, changing from a boxcar (Morganson et al. 2018) to a Gaussian kernel. To optimize the threshold, we compared the catalogs of sources detected in DES exposures to the significantly deeper Hyper Suprime-Cam Subaru Strategic Program (HSC; Aihara et al. 2018) XMM-LSS catalog from their first data release. By using HSC, not only do we have an

independent set of measurements (as the exposures are in different times, single-night transients from any DES exposure are not present on the HSC images), but this also guarantees completeness of the reference catalog to fainter magnitudes than any DES exposure. The DES sources that do not match HSC detections correspond to astrophysical transients and/or to spurious detections coming from noise fluctuations or defects in the image. Our choice of `DETECT_THRESHOLD` led to an increase of ≈3× in the transient density, and incrementally smaller increases in transient density are found at higher detection thresholds. Any lower choice of `DETECT_THRESHOLD` leads to a substantial increase in this value (>10×), clearly due to detection of noise fluctuations, making it infeasible to conduct the search despite possible small gains in completeness for TNOs.

Once all DES images have been cataloged using the new Y6 detection parameters, we evaluate the detection threshold for unresolved sources as follows. For a given DES single exposure, we find all of the unresolved sources in the coadd catalog that overlap the exposure, and see whether each appears in the SE catalog. The list of (un)matched sources is used to constrain a completeness function for each DES exposure. The completeness model is defined by the magnitude of 50% completeness $m_{50}$, a scaling factor c, and transition sharpness k such that the probability for a source with magnitude m being detected in this exposure is given by

$$p(m) = \frac{c}{1 + \exp(k(m - m_{50}))} \quad (1)$$

The left-hand side of Figure 1 plots histograms of the $m_{50}$ values derived for all of the accepted Y6 exposures. The median depths are 24.2, 24.0, 23.5, 22.7, and 21.5 mag in the grizY bands, respectively. These are a significant improvement over the depths attained in the Y4 processing, even for the same images. Figure 1 shows the histogram of changes in $m_{50}$ between the Y4 and Y6 processing for all exposures contained in both. The median gain in completeness is 0.48, 0.46, 0.43, 0.39, and 0.41 mag in the grizY bands, respectively.
<figure><img src="image_5.png" /><figcaption>Figure 2. Trend lines for the g − r, g − i, and g − z colors for the Y4 sample (black points), as well as a linear regression (blue lines) for this data. Objects that were clipped are shown in red. We find that { α$_{b}$, β$_{b}$ } is { 1.49 , − 0 . 12 } for the i band and { 1.65 , − 0 . 133 } for the z band.</figcaption></figure>

2.3. Injection of Synthetic TNOs

We inject a synthetic population of TNOs (“fakes”) into the SE catalog in order to test for transient recovery and linking efficiencies. Unlike in Paper I, this synthetic population is injected before the production of the transient catalog, allowing the simulation to measure the (in)efficiency of the transient-selection process as well as those of the linking process and selection cuts.

Similar to Paper I, we construct an ensemble of orbits isotropically distributed in the Cartesian coordinates { x$_{0}$, y$_{0}$, z$_{0}$, ˙ x$_{0}$, ˙ y$_{0}$, ˙ z$_{0}$ } for heliocentric distances 25 au < d < 2500 au and velocities 0 ≤ v ≤ v$_{esc}$ ( d ) so as to cover a wide range of sky locations, semimajor axes, and eccentricities. The linking process, as will be seen later, aims to find sources at d ≳ 29 au, and the injection of sources at closer distances allows us to verify our completeness range. We construct a sample of ≈ 5000 fakes intersecting the DES observations, from three different populations:

1. 80% of the fakes are uniformly distributed in barycentric distance between 25 and 60 au, with the angular positions and velocities isotropically distributed using the method of Paper I;

2. 10% of the population is logarithmically distributed between 60 and 2500 au;

3. The remaining 10% are generated by constructing a larger isotropic population between 25 and 60 au, but keeping only those with inclination i < 20 $^{◦}$.

We have elected to include a larger number of objects at lower distances or inclinations than in Paper I, as these are more representative of the real trans-Neptunian population, and larger numbers of synthetic TNOs are needed in these regimes to keep the shot noise of the simulations subdominant to the shot noise of the real detections.

We also generate for each fake a mean r -band magnitude uniformly distributed in 20 ≤ m$_{r}$ ≤ 24 . 5, and a color 0 . 4 < g − r < 1 . 5. From this color index, we generated the other colors ( g − i and g − z ) by a linear parameterization for band b

$$
g - b = \alpha_{b}(g - r) + \beta_{b}, \quad (2)
$$

with the values for α$_{b}$ and β$_{b}$ found by fitting to the observed colors of the 316 TNOs from the Y4 search (see Figure 2). The objective of this fit is not to measure correlations of TNO colors, rather just to get representative color variation into the fakes. So we do not include measurement errors while fitting, and we apply a σ -clipping algorithm to remove outliers that are more than 3 σ away from the fits in any of the three dimensions. The fitted color trends are shown in Figure 2.

Differences in the surface reflectance as well as in the shape of TNOs can introduce apparent-magnitude variability on timescales shorter than the DES cadence (see Alexandersen et al. 2019, and references therein). To simulate this, each object is also assigned one of three possible light-curve scenarios:

1. Constant light curve (fixed magnitude, 50% of the sample);

2. Light curve with variation δm ≡ A sin ϕ and amplitude A = 0 . 2 (25%);

3. Light curve with the same variation as before and A = 0 . 5 (25%).

These amplitudes were chosen to represent a wide range of TNO variability (Alexandersen et al. 2019), with peak-to-peak excursions of up to 1 mag across the 30 to 50 potential detections of each fake. Since the time between DES observations of any given object is usually many times longer than a typical TNO rotation period, we do not assign a period to each object, but rather, for each observation, we draw a random
phase $\varphi \in [0, 2\pi)$. The percentages are chosen to provide accurate statistics for each of the three light-curve scenarios, and we do not attempt to represent the underlying light-curve amplitude distribution of TNOs.

For each observation of each fake TNO, we draw a random uniform deviate $u \in [0, 1]$, and use Equation (1) to find a probability of detection. If $p(m + \delta m) \ge u$ (where $\delta m$ is simulated with a random phase for each exposure), the synthetic source is considered as detected and is added to the SE catalog. Finally, each detected source receives an astrometric shot-noise error $\sigma_a$. For a given flux $f$, a background-limited source has $\sigma_a \propto f^{-1}$, while for brighter sources, $\sigma_a \propto f^{-1/2}$. Since most sources of interest are faint, we fixed $\sigma_a \propto f^{-1}$, and used the Y4 transient detections for magnitudes between $20 < m < 22$ for each band to find the expected amplitude of this shot noise. This corresponds to $\sigma_a(m_b = 20) = \{5.8, 4.9, 6.5, 9.6\}$ mas for $b = \{g, r, i, z\}$, respectively. These shot-noise errors then get added in quadrature to each exposure's atmospheric turbulence covariance matrix (see Section 2.3 of Paper I), and each object has its nominal positions shifted by drawing from the two-dimensional Gaussian defined by this sum of covariance matrices.

Section 3 presents a discussion of the survey completeness versus the color, light-curve, and orbital characteristics of the fake TNOs.

## 2.4. Transient Identification

The transient identification algorithm is the same as the one presented in Paper I. We begin by matching all detected SE and coadd sources by using a friends-of-friends (FoF) algorithm that links detections within 0.5" of each other. For each output group, we ask whether or not this source is matched to a coadd detection (COADD = 1, 0). Next we evaluate how far apart in time the first and the last detections of this match group are ($\Delta t \equiv t_{\text{last}} - t_{\text{first}}$). Then, for the groups that have both a coadd and an SE source, we ask what is the least negative deviation in magnitude between an SE detection and the coadd in the same band, $\Delta m \equiv m_{\text{SE}} - m_{\text{coadd}}$.

The coadd image at the location of an SE detection of a solar system object has this object's flux in only one exposure out of K in the same band and location; that is, the coadding process reduces the flux of a solar system object by a factor $1/K$, corresponding to a $\Delta m \lesssim -1.2$ for $K \ge 3$. We include in the transient catalog detections that satisfy

1. COADD = 1, $\Delta t < 2$ days, $\Delta m \le -1.2$; OR

2. COADD = 0, $\Delta t < 2$ days.

This process yields 160 million transient candidates, with the injections of fake TNOs detected with efficiency of 98.6%.

After these transients are identified, the first cleaning step is to remove large clusters of spurious detections. We apply an FoF algorithm with a linking radius of 30" to the catalog of all transients, and discard groups with over 20 detections, or groups with between 10 and 20 members whose detections have a tendency to lie along a line or a curve (see Section 2.4 of Paper I), masking 17 million transients (11%) and only 0.2% of the fakes.

## 2.5. Coadd Avoidance Radius

Static sources with complex shapes and/or low S/Ns can be split into multiple detections in ways that differ from exposure

<figure><img src="image_6.png" /><figcaption>**Figure 3.** Cumulative histogram of the distance between each identified transient and the closest coadd source, measured for all putative transients (blue), and for the fake detections injected in the catalog (red). The blue curve shows a steep growth for distances $\lesssim 1''$, while the fakes, which are uniformly distributed over the footprint, show much lower fractions with close coadd—24.3% (1.8%) of the transients (fakes) are closer than 1" to a coadd source. This indicates that the close pairs are primarily associated with static sources, and not true solar system transients, so we remove them from the transient catalog.</figcaption></figure>

to exposure. This can generate SE detections that spread in radius more than the 0.5" used for the friends-of-friends matching. Such detections can then be mistakenly identified as transients.

Figure 3 shows in blue the cumulative histogram of the distance between each transient and its closest coadd source (for those with COADD = 0). The red line shows the same quantity for the injected fakes, i.e., a truly randomly placed population. The significant excess of the former at low radii indicates that, for distances less than $\approx 1''$, the putative transients are not independent of the static coadd sources, and contamination such as that hypothesized above is present.

We therefore exclude from the transient catalog any SE detection that is within a “coadd avoidance radius” of 1" of a coadd catalog source, but had COADD = 0. This process removes 32 million (24.3%) of all identified transients, and only 1.8% of the fakes. This percentage has a weak dependence on galactic latitude $b$, increasing to 2.1% within $|b| \sim 20^\circ$, and dropping to 1.5% for sources farther than $|b| > 70^\circ$.

## 2.6. Pixel-level Masks

A final masking process has been applied to the transient catalog by combining the data from each CCD in an observing “epoch” (roughly two epochs per year; see Table 2 of Morganson et al. 2018) to search for additional unmasked spurious detections, such as unmasked bad columns. True astrophysical transients should be randomly located on each CCD, while signals due to CCD defects will be clustered in pixel coordinates. We count all detections in bins of size 8 × 16 pixels in each 2048 × 4096 pixels CCD image, and define a Poisson distribution with the mode of transient counts in each CCD bin. We mask all bins whose counts are above the 99.99% percentile of this Poisson distribution and at least twice the median bin count. Figure 4 plots an example CCD, showing that this process finds an excess of transients in low-x pixel coordinates in the CCDs, as well as bad columns that were not
<figure><img src="image_7.png" /><figcaption>**Figure 4.** Example of the transient density per 8 × 16 pixel bin in a CCD, with the threshold for masking set to 10 transients per bin. The bottom panel shows the location of the masks in yellow, while the blue region is the unmasked portion of the CCD. A few features can be seen here, with an excess of transients (and, therefore, of masked regions) in the low-x pixel counts of the CCD, as well as two bad columns in the CCD.</figcaption></figure>

<table><caption>**Table 2**

Completeness Fits for a Selection of Subsets of the Simulated Fake Population</caption><thead><tr><th>Subset of the Simulation</th><th>r<sub>50</sub></th><th>c</th><th>k</th><th>Subset of the Simulation</th><th>r<sub>50</sub></th><th>c</th><th>k</th></tr></thead><tbody><tr><td><strong>All fakes</strong></td><td><strong>23.77</strong></td><td><strong>0.943</strong></td><td><strong>6.37</strong></td><td>i &lt; 10°</td><td>23.67</td><td>0.944</td><td>6.54</td></tr><tr><td>Constant light curve</td><td>23.76</td><td>0.944</td><td>6.52</td><td>i &lt; 20°</td><td>23.68</td><td>0.938</td><td>6.24</td></tr><tr><td>A = 0.2 mag</td><td>23.77</td><td>0.943</td><td>6.38</td><td>i &lt; 50°</td><td>23.72</td><td>0.946</td><td>6.28</td></tr><tr><td>A = 0.5 mag</td><td>23.78</td><td>0.944</td><td>6.24</td><td>i &lt; 90°</td><td>23.77</td><td>0.951</td><td>6.43</td></tr><tr><td>0.4 &lt; g - r &lt; 0.7</td><td>23.75</td><td>0.950</td><td>6.73</td><td>d &lt; 45 au</td><td>23.76</td><td>0.930</td><td>6.21</td></tr><tr><td>0.7 &lt; g - r &lt; 1.0</td><td>23.75</td><td>0.952</td><td>6.80</td><td>d &lt; 60 au</td><td>23.77</td><td>0.943</td><td>6.37</td></tr><tr><td>1.0 &lt; g - r &lt; 1.3</td><td>23.81</td><td>0.951</td><td>6.70</td><td>d &lt; 90 au</td><td>23.77</td><td>0.943</td><td>6.37</td></tr><tr><td>1.3 &lt; g - r &lt; 1.5</td><td>23.88</td><td>0.948</td><td>6.92</td><td>d &lt; 500 au</td><td>23.77</td><td>0.943</td><td>6.37</td></tr></tbody></table>

**Note.** These subsets cover a different range of light-curve amplitudes A, g - r colors, inclinations i, and distances d. Other choices of parameter bins were tested (such as e), leading to no significant change among bins. The r$_{50}$ entry is the light-curve-averaged r-band magnitude where the probability of discovery of an implanted TNO is 50%, c is a scaling factor for the completeness curve, and k is its transition sharpness (see Equation (1)). The results for the complete set of fakes are shown in bold.

previously masked. We have not investigated the nature of the low-x transient excess, but both the directional nature of SExtractor's operation on an image, and CCD defects leading to charge smearing at the ends of columns, are potential explanations for this phenomenon.

This leads to 2.7% of the remaining transients being masked, and 1.2% of the fakes injected in the catalog. While this masking removes only twice as many transients as fakes, the masked transients had an outsized contribution to spurious linkages because of their concentration into small regions of

each image. Since the population of fakes is uniformly distributed in the images, we can conclude that 1.2% of the total image area is being masked by this process.

## 2.7. The Transient Catalog

The final transient catalog has 108 million sources (compare to 22 million in the Y4 catalog). Of these, 105,317 come from the fake detections injected into the SE catalog (of 110,246 originally injected), leading to an overall transient efficiency
<figure><img src="image_8.png" /><figcaption>**Figure 5.** Sinusoidal projection of the DES footprint in ecliptic coordinates; each dot corresponds to one *r*-band exposure, with its color representing transient counts per exposure and per square degree. The left panel counts transients with *r* < 23 mag. In this regime, most exposures are complete, and the transient density increases strongly toward the ecliptic plane (black line), reflecting density variations in the number of astrophysical transients (asteroids). The right panel corresponds to transients with *r* > 23. The faint end has three to five times more detections than the bright end. The lack of spatial dependence (except for a small decrease in transient counts for longitudes less than −30°) suggests that these are due primarily to noise detections, image artifacts, and/or other non-solar system transients.</figcaption></figure>

for moving object detection of 95.55%. Figure 5 shows the distribution of transients in the survey’s footprint.

The density of bright transients shows a strong concentration toward the ecliptic plane (to ≈200 per square degree for *r* < 23) at a level consistent with the expected density of asteroids (e.g., Gladman et al. 2009). The transient density increases significantly for fainter detections (*r* > 23), which dominate the catalog at all ecliptic latitudes. The mean transient density is a factor of 3–5 times the transient density of the Y4 catalog, primarily a consequence of the lower detection thresholds; we have lowered the purity of true solar system transients in an effort to increase completeness.

## 2.8. Orbit Linking

The process of linking of detections into orbits is very similar to the one described in Paper I. Pairs of detections are found by searching in bins of inverse distance $\gamma \equiv 1/d$ (Bernstein & Khushalani 2000), with $29 < d < 2500$ au. The detections are mapped to a frame that subtracts Earth’s parallax at that distance, and so the dominant motion becomes linear in time. The pairs are found by first finding all pairs of exposures within some time $\Delta t$ of each other that could contain a common bound TNO within the distance bin. Then kD trees of both exposures’ detections are searched for pairs of detections with separation that is consistent with the motion of a bound orbit. While the pair-finding algorithm remains unchanged, we have changed its implementation to use VAEX table software (Breddels & Veljanoski 2018). The VAEX out-of-core functionality works well for very large data sets, e.g., some of our distance bins generate billions of pairs.

A pair of exposures strongly constrains four out of the six orbital degrees of freedom, essentially leaving distance and line-of-sight motion weakly constrained. The triplet stage proceeds by determining the two-dimensional region of a future exposure spanned for plausible variations of these two parameters, and locates all detections within this region. The algorithm of Paper I used a kD-tree implementation for this stage, but for the Y6 search, we devise a linear-algebra-based “parallelogram” search, which is faster. The details are presented in the Appendix. We consider only triplets whose dates of observation $t_1, t_2$, and $t_3$

satisfy $|t_2 - t_1| < 60$ days and $|t_3 - t_2| < 60$ days. When searching for TNOs at distances $d > 50$ au, we increase these windows to 90 days, since the number of triplets decreases steeply with distance, and we can search larger time intervals without being overwhelmed by spurious triplets.

Once a triplet of detections is found, the search then proceeds by fitting these orbits using the procedures outlined in Bernstein & Khushalani (2000) and Paper I. An “n-let” of linked detections is fit to an orbit with the six orbital elements left free, but with a tight Gaussian prior on the inverse distance. All exposures are then searched for transient detections lying within the $4\sigma$ predicted error ellipse of its position in that exposure. If a new transient is found to be consistent with the orbit, an $(n+1)$-let is created, and the orbit is refit. This is iterated until no new transients are consistent with the orbit. Such “terminal” n-lets are retained as TNO linkage candidates.

For each linkage, we compute the following quantities:

1. The $\chi^2$ of the orbit fit, following the routines of Bernstein & Khushalani (2000). Here, the number of degrees of freedom $\nu = 2 \times \text{NDETECT} - 6$, with NDETECT being the number of detections in this orbit. We reject all orbits with $\chi^2/\nu > 4$.

2. The number of unique nights NUNIQUE $\le NDETECT$ on which detections were made. This is more indicative than NDETECT of the chances of accidental linkage of asteroid apparitions into a TNO orbit, because the DES observing algorithm (Neilsen & Annis 2013) occasionally chooses to take successive exposures (at 2 minute intervals) with the same pointing in a single night. The short interval between these repeated pointings means that these intra-night exposures have a highly correlated chance of containing an asteroid or image defect. We keep only orbits with NUNIQUE $\ge 7$ (this choice is discussed in Sections 2.9 and 3.2).

3. The false-positive rate (FPR) of the linkage. With *j* being an index over all exposures, $A_{j,search}$ being the area in that exposure consistent with the orbit fit at $4\sigma$, and transient density $n_j$ in exposure *j*, the FPR for a spurious linkage is

$$ \text{FPR} = \sum_j A_{j,search} n_j. \qquad (3) $$
<figure><img src="image_9.png" /><figcaption>**Figure 6.** Results from the STS evaluation of 3918 sources with NUNIQUE $\ge 7$. The left panel shows a histogram of the STS value for all sources with NUNIQUE = 7 (>7) in blue (red). The right panel plots the number of images used in each STS stack vs. the STS values, with color and shape encoding the NUNIQUE value and eyeball classification, respectively. The STS = 1.2$\sqrt{N_{images}}$ curve separates the bulk of low-STS sources from the ones considered as real in the eyeball test. The long tail of the blue histogram in the left panel and the separation of the two groupings of the right panel indicate that real sources do reliably have signals in all griz images, even when they are not detected in individual exposures. The 814 objects to the right of the curve or with NUNIQUE $\ge 9$ comprise our final detection catalog.</figcaption></figure>

<figure><img src="image_10.png" /></figure>

We calculate the FPR for the detection last linked to the
orbit, and retain only those with FPR < 0.02. A larger
minimum value of FPR would increase substantially the
number of linkages to be tested, with no significant gain
in recovery rate of the synthetic objects.

4. The ARC, corresponding to the time between the first and last detections of the orbit, and ARCCUT < ARC, the shortest arc that remains after any single night of detections is removed from this orbit. We keep only orbits with ARCCUT > 6 months—that is, orbits with $\ge 2$ detections outside the season of the discovery triplet.

After making these cuts, all duplicate linkages are merged,
and every orbit is refit with no priors on distance or binding
energy. This leaves 760 candidates with NUNIQUE $\ge 8$ and
another 8321 with NUNIQUE = 7.

## 2.9. Sub-threshold Significance

To assess the reliability of a linkage, we implement the sub-threshold significance (STS) test from Paper I: we stack all images in the griz bands where the object is not detected, as a real object will have some signal lurking below the detection limit of the individual images, such that stacking will lead to a significant signal. The STS is the significance of the flux peak in a stack of SE images centered on the positions predicted from the orbit fit, measured inside a 1" FWHM aperture. The key is that all exposures that are already linked into the orbit are excluded from the stack, as are images taken on the same night as linked detections. This leaves only images that are statistically independent of the original linkage in terms of both shot noise and presence of asteroids or defects. With Nimages griz exposures kept in the STS stack, we examine by eye all objects with STS > 1.0√Nimages to eliminate any images where the STS aperture is contaminated by static sources or unmasked artifacts. Three of the authors (P.B., G.B., and M.S.) examined by eye postage stamps and stacked images of all candidates, and attributed to each a score of R(real), M(maybe), and F(false); if two of the scores agree, this is this object's final classification. As seen in the right-hand plot of

Figure 6, the curve STS = 1.2√Nimages perfectly separates
those graded R from those graded F. We therefore consider as a
confirmed object:

1. any candidate with NUNIQUE $\ge 9$ (independent of its STS),

2. and any candidate with STS > 1.2√Nimages (after removing contaminated images from the STS stack).

Of the 20 candidates graded M, 13 are above this STS
threshold.

Our final catalog consists of those 814 objects that pass at least one of these criteria. Figure 6 shows the results of the STS of all 9081 sources with NUNIQUE $\ge 7$ found in the search.

We also computed the STS statistic for ≈100,000 orbits with NUNIQUE = 6 found at $d > 50$ au. Visual inspection of all sources with STS > 1.5√Nimages yielded only one source whose STS value was not spuriously high due to contamination by unmasked artifacts or static sources. This is contrary to the result in Paper I, where most NUNIQUE = 6 orbits came from real linkages. Because examination of all NUNIQUE = 6 candidates would be a large undertaking with very little gain in secure TNOs, we elect to set NUNIQUE = 7 as our minimum threshold for the Y6 TNO search. The modest effects of this choice on TNO completeness are discussed in Section 3.2.

## 2.10. False Negatives

The combination of the NUNIQUE and STS criteria gives us
confidence that the final catalog does not contain false-positive
detections. During the peer review process, however, two of
our orbits with NUNIQUE = 7 and one with NUNIQUE = 8
were noted by the Minor Planet Center<fnref n="49" /> to contain detections
that could be linked to previously known TNOs, on orbits
inconsistent with our linkages. In each case, a string of real
detections in a single DES season were falsely linked to two
spurious detections in other DES seasons. Closer investigation
of whether this problem could be more prevalent in our sample

<footnote n="49">⁴⁹ https://www.minorplanetcenter.net/</footnote>
shows that these three cases were the lowest-quality ones from the Y6 search and therefore unlikely to represent a larger set of erroneous linkages. In all three cases, there was disagreement among the visual inspection grades, but a passing value STS. Two of these cases only had images for the STS in the year of the discovery triplet, before the erroneous linkage “drifted off” the true orbit. The final case had only one image used for the STS verification, also in the same year as the discovery triplet. These were the only three cases in all of our catalog where these conditions occurred.

Now we ask whether we have false negatives, i.e., objects we should have detected but missed. Our proposition is that our final catalog is almost 100% complete for TNOs whose appearances in the transient catalog satisfy all of the following criteria for being “discoverable”:

1. At least one triplet that fits the timing criteria in Section 2.8;

2. NUNIQUE ≥ 7;

3. ARCCUT > 6 months.

4. Heliocentric distance in 2016 > 29 au.

One caveat (discussed further in Section 3) that prevents us from achieving 100% completeness is that regions around relatively bright stars that are bright enough to be removed as part of the DES processing are not taken into account in our estimates, leading to a potential overestimation of our completeness near the detection limits of the survey.

### 2.10.1. Synthetic TNOs

The efficiency of our linking process can be judged, first, by seeing whether the linker successfully found all of the TNOs whose synthetic detections injected into the SE catalogs (as per Section 2.3) created a set of transients meeting the discoverability criterion. We find that 3710 out of 3749 (≈99%) such objects were in fact successfully linked and “discovered.” All 39 of the missed linkages were in regions of the footprint with high transient density, suggesting that these did not satisfy the FPR criteria during the linking process.

There is one caveat here, which is that we do not have postage stamp images for the fakes and thus cannot calculate an STS value. The fairly clear separation between high- and low-STS populations in Figure 6 strongly suggests that >99% of real objects do lie above the threshold.

### 2.10.2. Known Objects

We can also address the false-negative question by asking whether all previously known TNOs that crossed the DES footprint made it into our catalog. We start by comparing the Y6 catalog to the Y4 catalog. Six objects out of the 316 that were found in the Y4 search are missed in the Y6 search. The objects, and the reasons they are missing in this search, are:

1. 2003 QT90 and 2013 VQ46: these objects had NUNIQUE = 6 in this search (as in the Y4 search), and are on the edges of the DES footprint, not allowing for further observations in observing years 5 and 6;

2. 2010 SB41 and (534315) 2014 SK349, which were discovered at distances *d* = 28.2 au and *d* = 27.8 au, respectively, are too close to be found in the smallest characterized distance bin for this search;

3. 2013 RL124: the only possible discovery triplet for this object in the Y6 search included a pair of detections 73

days apart, and the maximum time difference for pairs of detections in this search for *d* < 50 au is 60 days (see Section 2.8; note that this restriction was not present in the Y4 search), meaning that triplets of this object were never searched for further detections in Y6;

4. 2013 SL106: this is one of our objects that was mis-linked in this search. The true orbit for this TNO yields NUNIQUE = 7, but ARCCUT < 6 months. One of the detections from the Y4 search that led to a successful ARCCUT was masked due to its proximity to a coadd source.

We thus see that none of these five met the discoverability criterion for the Y6 search.

We also compare our list of objects to an updated list of known objects in the Minor Planet Center using the methodology of Banda-Huarca et al. (2019). This list leads to a total of 591 objects that cross the DES footprint, including objects that do not satisfy our detection criteria. Of these, 356 have detection sets that are “discoverable” by our criteria, and 100% of these have indeed been retrieved in this processing.

In summary, the evidence from both simulated and real TNO observations is that our linking process is highly complete (>99%) for TNOs whose brightness and geometry produce a discoverable set of entries in the transient catalog.

## 3. Survey Simulation

### 3.1. Simulator Methods

In order to quantify our observational biases, we have developed software that allows observations of population models to be simulated, similar to Jones et al. (2006) and Lawler et al. (2018). Each simulated object requires a set of orbital elements, magnitudes (either apparent or absolute) in griz and potential light curve, allowing for a wide range of parameters to be varied in these tests. As an example, the populations of fakes described in Paper I and Section 2.3 as well as the extreme TNO ensemble of Bernardinelli et al. (2020b) were realized within this framework.

For each set of orbital elements, we find all potential observations of such a TNO; that is, we find every DES pointing in which such an object's position would lie inside a functional DECam CCD. We have recorded completeness estimates {$m_{50}$, *c*, *k*} for all DES exposures used in the search, and so we can iterate over the exposures and calculate the probability of the TNO being detected using Equation (1) given its magnitude and light curve. An additional factor of 95.5% gives the probability that a detection will be correctly identified as a transient (as determined in Section 2.7). A random-number draw determines whether this TNO yields a transient detection on each of its observable exposures.

From the list of transient detections for a posited TNO, we determine whether the TNO is discoverable by the criteria listed in Section 2.10. Given the results of that discussion, we can safely assume that any posited TNO that is discoverable would, in fact, have been discovered by the survey.

There is one caveat to our simulation algorithm, which is that the stellar objects used to determine the detection probabilities (Section 2.2) are, by definition, only drawn from parts of the sky that yield usable images. In particular, regions of the sky near very bright stars or galaxies, in which neither TNOs nor static stellar sources would have been found, are not properly accounted for by this method. We are thus slightly
<figure><img src="image_11.png" /><figcaption>**Figure 7.** Measured completeness as a function of r-band magnitude for the same subsets of the simulated fake population presented in Table 2. The left panel presents subsets of photometric properties (light curves and colors), and the right panel shows subsets of orbital properties (inclination and distance).</figcaption></figure>

<figure><img src="image_12.png" /></figure>

overestimating the parameter *c* in Equation (1), perhaps by a few percent. The brightest stars, galaxies, and globular-cluster regions are totally removed from DES processing and are not considered as “observed”—these are properly excluded. But the “deserts” around less-bright stars are not. We expect this issue to cause, at worst, a few percent misestimation of our total efficiency for TNOs near the completeness limits.

This survey simulation software, as well a tutorial for its use, is available on GitHub.<fnref n="50" /> The software works with the user stipulating a set of orbital elements or phase space coordinates, an absolute or apparent magnitude and colors for each object. Routines to include the color ranges of Equation (2) and light curves are also included. The software then evaluates if this object would be detected or not, and returns this information.

## 3.2. Completeness Testing

We show here the results of a simulation of a large population of synthetic TNOs, and the effects of selecting distinct subsets from this sample. The population of fakes is simulated in a similar manner as in Section 2.3, but this population is a factor of 100 larger (yielding ≈450,000 fakes) than those injected into the SE catalogs, and covers the range $20 < m_r < 25$. This time, however, we use the detectability simulations described above rather than actually linking the detections.

We aim to summarize the detection probability as a function of mean apparent *r*-band magnitude of each simulated source (averaging over any light curves), subdividing the simulated population into subsets of color, light-curve behavior, distance, and inclination, as well as the “overall” efficiency that combines all of the simulated objects. Figure 7 shows the completeness fraction of the simulation as a function of mean $m_r$ within each subset. We fit the logit function (Equation (1)), see also Section 2.6 of Paper I) to the detection efficiency curves and report the fitted parameters in Table 2. This choice of function slightly overestimates the completeness for $m_r$ slightly brighter than $r_{50}$, but the logit is an acceptable approximation for our purposes.

It is clear from this exercise that the detection efficiency of the DES Y6 TNO survey is essentially a function only of mean

$m_r$ for objects that do not move out of the DES footprint during the survey. For the population as a whole, the survey is 50% complete at $m_r = r_{50} = 23.77$ mag, a significant improvement on the value of $r_{50} = 23.3$ mag in the Y4 search. This procedure also allows us to determine a magnitude of 50% completeness in the *g*, *i*, and *z* bands, yielding $g_{50} = 24.71$, $i_{50} = 23.42$, and $z_{50} = 23.28$ mag. We note, however, that the underlying assumption of a uniform *g*–*r* color distribution and *r*–*i*, *r*–*z* colors is derived as in Equation (2).

The total footprint area is 5090 deg², and the value of *c* = 0.943 leads to an effective search area of 4800 deg². The amplitude of a TNO light curve changes $r_{50}$ by only ±0.01 mag for $0 < A < 0.5$ mag, and the distance to the TNO also makes no appreciable difference once the TNO is in the search region of $d > 29$ au. There is a mild dependence of $r_{50}$ on color, changing by 0.13 mag from the bluest to reddest simulated objects. The “pivot point” of the four-band DES survey, i.e., the wavelength at which detection efficiency is nearly independent of TNO color, is slightly redward of the *r* band.

The completeness drops slightly (≈0.1 mag brighter) for $i < 10, 20^{\circ}$ when compared to all of the other subsets. This is a direct effect of the shape of the DES footprint (Figure 5), which is very narrow for ecliptic latitudes $|\beta| \lesssim 15^{\circ}$: low-latitude objects (especially at closer distances) have a greater chance of moving in or out of the footprint during the DES duration. This reduces the number of opportunities for the TNO to meet the detectability criteria, causing loss of fainter sources.

The inclusion of orbits with NUNIQUE = 6 would lead to $r_{50} = 23.81$ and increase the effective search area to 4812 deg², a minor gain compared to the increased burden of verification both by STS and visual inspection of potentially hundreds of thousands of orbits.

A very good approximation to the DES Y6 TNO search selection function, therefore, is that TNOs follow the logit function with $r_{50} = 23.77$ mag, for any TNO that spent most of 2013–2019 inside the DES footprint at a distance beyond 29 au.

# 4. Catalog of DES TNOs

## 4.1. Dynamical Classification

We divide the catalog of confirmed Y6 TNOs into dynamical classes following a methodology similar to Gladman et al.

<footnote n="50">⁵⁰ https://github.com/bernardinelli/DESTNOSIM</footnote>
(2008) and Khain et al. (2020). We sample the barycentric Cartesian phase space position and covariance matrix of the orbital solution at its solution epoch ($t_0$) to generate 30 clones of each orbit, and integrate these for 10 Myr in time steps of 30 days using the WHFAST implementation of the symplectic Wisdom-Holman mapping (Wisdom & Holman 1991; Rein & Tamayo 2015), part of the REBOUND N-body simulator (Rein & Liu 2012). The simulation is such that Jupiter, Saturn, Uranus, and Neptune are treated as active particles, and the mass of the terrestrial planets is added to the solar mass, as in the orbital solutions (Bernstein & Khushalani 2000).

Similarly to Elliot et al. (2005), for each clone, we check all resonant arguments of the form

$$ \sigma = p\lambda - q\lambda_N + m\varpi - n\varpi_N + r\Omega - s\Omega_N, \quad (4) $$

where $\varpi \equiv \Omega + \omega$ is the longitude of perihelion, and $\lambda \equiv \varpi + M$ is the mean longitude. The subscript “N” refers to Neptune’s orbital elements, and we check all $p, q, |m|, |n|$, and $|r| \in [1, 29]$ such that $p - q + m - n + r = 0$. For convenience, we define $s = 0$, as the number of possible resonant arguments becomes significantly larger if this parameter is also left free.

We use an automated resonance identification algorithm similar to the one outlined in Khain et al. (2020): for each resonance whose nominal orbital period $P_{res}$ is $|P_{res} - P(t_0)| < 0.15 P_N$, that is, whose orbital period is up to 15% of Neptune’s orbital period away from the resonance, we construct a two-dimensional histogram of $(t, \sigma)$, and, for each time slice, we check whether the values of $\sigma$ are bounded; that is, we check if there are multiple empty $\sigma$ bins in each $t$ slice. To achieve this, the bins require 1000 points in each window of 200 kyr, and we manage to identify resonances independently of the libration center. If the largest contiguous interval of libration spans more than 90% of the full integration, we consider this clone to be resonant.

We classify the nonresonant objects as follows:

1. Objects with $a(t_0) > 2000$ au are classified as Oort cloud objects;

2. Objects with $a(t_0) < a_N$ are considered Centaurs;

3. As in Khain et al. (2020), an object is considered to be scattering if it experiences excursions in semimajor axis as a function of time $a(t)$ such that

$$ \frac{\max|a(t) - a(t_0)|}{a(t_0)} > 0.0375; \quad (5) $$

4. Following Gladman et al. (2008) and Khain et al. (2020), detached objects are nonscattering, nonresonant orbits that have $e(t_0) > 0.24$; and

5. Objects that do not fall into these categories are classified as Classical.

Each object gets assigned the dynamical class indicated by the behavior of the majority of its clones. As in Khain et al. (2020), we require at least 50% of the clones to present scattering behavior to assign this classification, and resonant objects are called “candidates” if over 50% of their clones exhibit resonant behavior, while those with over 80% of the clones showing libration of the same resonant argument are identified as securely resonant. A minority of objects do not satisfy either dynamical requirement for the majority of their clones, and receive an “insecure” classification.

## 4.2. Sample of TNOs

The final catalog contains 814 objects. Of these, 503 have not been found in the Y4 processing, including the 458 newly identified objects. One of the discoveries is a Centaur, and one is an Oort cloud comet discovered interior to the $d > 29$ au region in which we did a complete search for bound sources. Therefore, there are 812 TNOs in the fully characterized sample. This is the second largest catalog of TNOs to date (OSSOS; Bannister et al. 2018, has 840 objects, 818 in the characterized TNO phase space), and the largest with multi-band photometry accompanied by orbital arcs of multiple years. Figure 8 shows the semimajor axes $a$, eccentricities $e$, and inclinations $i$ for all bound objects, as well as their dynamical classifications. This figure excludes the one object with $a > 10, 000$ au, discussed below. Table 3 describes the parameters included in the released catalog for each object, and the full table is provided in machine-readable format. A summary of the dynamical classification is presented in Table 4.

The following properties of the sample, as well as particularly interesting objects, are highlighted:

1. Object C/2014 UN$_{271}$ (Bernardinelli-Bernstein), inbound on a near-parabolic cometary orbit, with $e = 0.999419 \pm 0.000011$, $q = 10.95$ au, $i = 95.56^{\circ}$, and nominal discovery distance $d = 27.53$ au, incoming from the Oort cloud. This is also the only one of the objects discovered at $d < 29$ au, and is therefore not part of our well-characterized search space. This object is discussed in detail in Bernardinelli et al. (2021).

2. Two new discoveries at $d > 70$ au: the first is 2014 US$_{277}$, our largest-$q$ object with $q = 54.31$ au, discovered at a distance of 79 au (the third most distant object from this search). While most clones of this orbit are nonresonant and this object is classified as detached, a portion of its clones were identified to be in the 6:1 resonance. The second is 2014 YC$_{92}$, a scattering object with low-$q$ discovered far from its perihelion, at $d = 73$ au.

3. Several new objects that securely occupy distant, high-order mean-motion resonances: for example, we report one object in each of the 8:1 (2014 RV$_{86}$) and 25:3 (2014 RL$_{86}$) mean-motion resonances (see Volk et al. 2018 for a discussion of a sample of 9:1 resonators). An object is also classified as a resonant candidate for the 10:1 resonance (that is, over half of its clones are in the 10:1 resonance, but not more than 80% as defined above).

4. Nine new “extreme” TNOs (eTNOs; $a > 150, q > 30$ au), increasing the total measured by DES to 16. Of these, nine (four newly reported here) also have $a > 230$ au and are of interest to the Planet 9 hypothesis (Batygin et al. 2019; see also Shankman et al. 2017; Bernardinelli et al. 2020b; Napier et al. 2021). These objects are discussed in Section 5.2.

5. Four new Neptune Trojans, increasing our sample size to ten, as well as the unrecovered 2013 RL$_{124}$.

6. Several new objects with high inclinations, including a detached object with $a = 53.28$ au and $i = 62.84^{\circ}$, the most inclined object with $q > 30$ au to date, and two new objects in classical orbits with inclination $i > 54^{\circ}$.

7. Two objects in low-$e$ (classical) orbits with $a > 55$ au, on orbits with high inclinations ($i > 29^{\circ}$) and perihelia
<figure><img src="image_13.png" /><figcaption>**Figure 8.** Semimajor axes, eccentricities, and inclinations of the 813 TNOs (excluding the object found in a cometary orbit) found in this search, color-coded by dynamical class. The dotted lines represent constant perihelion $q = a(1 - e)$, and the vertical dashed lines show the approximate locations of some *p:*q mean-motion resonances. Since many of these have similar *a*, we do not indicate all occupied resonances. The solid circles represent objects with secure dynamical classifications, and the triangles represent objects whose classification is insecure (including resonant candidates). Table 4 presents the number of objects found in each dynamical class.</figcaption></figure>

($q > 48$ au), being similar to the detached population, despite their assigned dynamical classes.

8. A total of 154 objects securely classified as detached, the largest sample of this population found by a single survey to date.

The discovery distance and magnitude of these objects are shown in Figure 9, as well as the 50% completeness limits of this search and the Y4 search. Most objects new to this search are found at $m_r > r_{50,Y4} \approx 23.3$, as expected. However, 30 new “bright” ($m_r < 23$) objects were found in the Y6 search, since the two additional years of data allow for new objects to enter the survey’s footprint for long enough to meet our selection criteria (see Figure 9).

**5. Initial Implications**

**5.1. Comparison to CFEPS Model**

The most accurate model to date for the classical TNOs is that derived by the CFEPS project (Kavelaars et al. 2009; Petit et al. 2011; Gladman et al. 2012). In their “L7” model, the classical population is described as a combination of a low-*e*, low-*i* “kernel” in a narrow *a* range, a “stirred” component of stable,

low-*i* orbits and a high-*i* “hot” component. The semimajor axis distribution of the CFEPS-L7 model is found to be consistent with the classical TNOs detected in the early data from the OSSOS survey (Bannister et al. 2016), but the inclination component has been found to be inconsistent with the CFEPS high latitude extension data (HiLat; Petit et al. 2017). Here we ask whether the main belt classical TNOs ($40 < a < 47$ au, roughly the region between the 3:2 and 2:1 resonances) in the DES Y6 catalog are consistent with the CFEPS-L7 model, or if these improved data require the model to be revised.

The CFEPS-L7 model consists of a simulated population of 26,031 members of the main belt with $H_g < 8.5$. We simulate the DES observation of this proposed population as in Section 3. The orbital elements and $H_g$ provided for each body by CFEPS-L7 need to be assigned colors for the DES simulation. For each member of the kernel and stirred components, we draw a synthetic $g-r$ color from the observed distribution in the DES classicals with $i < 5°$, and for the hot component, we draw colors from the DES classicals at $i > 5°$. The $g-i$ and $g-z$ colors are obtained following Equation (2). The simulation predicts 695 DES detections, about 3% of the total CFEPS-L7 population in this region.
<table><caption>**Table 3**
DES Trans-Neptunian Objects for the Y6 Release</caption><thead><tr><th>Column Name</th><th>Unit</th><th>Description</th></tr></thead><tbody><tr><td>MPC</td><td></td><td>Minor Planet Center object designation</td></tr><tr><td>a (a)</td><td>au</td><td>Semimajor axis of the best-fit orbit</td></tr><tr><td>σ<sub>a</sub> (sigma_a)</td><td>au</td><td>Uncertainty in a</td></tr><tr><td>e (e)</td><td></td><td>Eccentricity</td></tr><tr><td>σ<sub>e</sub> (sigma_e)</td><td></td><td>Uncertainty in e</td></tr><tr><td>i (i)</td><td>deg</td><td>Inclination</td></tr><tr><td>σ<sub>i</sub> (sigma_i)</td><td>deg</td><td>Uncertainty in i</td></tr><tr><td>ω (aop)</td><td>deg</td><td>Argument of perihelion</td></tr><tr><td>σ<sub>ω</sub> (sigma_aop)</td><td>deg</td><td>Uncertainty in ω</td></tr><tr><td>Ω (lan)</td><td>deg</td><td>Longitude of ascending node</td></tr><tr><td>σ<sub>Ω</sub> (sigma_lan)</td><td>deg</td><td>Uncertainty in Ω</td></tr><tr><td>T<sub>p</sub> (T_p)</td><td>UTC Modified Julian date</td><td>Time of perihelion passage</td></tr><tr><td>σ<sub>T</sub> (sigma_T)</td><td>days</td><td>Uncertainty in T<sub>p</sub></td></tr><tr><td>q (q)</td><td>au</td><td>Perihelion distance</td></tr><tr><td>σ<sub>q</sub> (sigma_q)</td><td>au</td><td>Uncertainty in q</td></tr><tr><td>d (d)</td><td>au</td><td>Discovery distance (geocentric)</td></tr><tr><td>σ<sub>d</sub> (sigma_d)</td><td>au</td><td>Uncertainty in d</td></tr><tr><td>Δ (delta)</td><td>au</td><td>Discovery distance (heliocentric)</td></tr><tr><td>σ<sub>Δ</sub> (sigma_delta)</td><td>au</td><td>Uncertainty in Δ</td></tr><tr><td>m<sub>r</sub> (m_r)</td><td>mag</td><td>Mean r-band magnitude</td></tr><tr><td>σ<sub>m</sub> (sigma_m)</td><td>mag</td><td>Uncertainty in m<sub>r</sub></td></tr><tr><td>H<sub>r</sub> (H_r)</td><td>mag</td><td>Absolute magnitude in band r</td></tr><tr><td>σ<sub>H</sub> (sigma_H)</td><td>mag</td><td>Uncertainty in H<sub>r</sub></td></tr><tr><td>NUNIQUE</td><td></td><td>Number of unique nights of detections</td></tr><tr><td>NDETECT</td><td></td><td>Number of detections</td></tr><tr><td>χ extsuperscript2 (CHI2)</td><td></td><td>χ extsuperscript2 of the orbit fit, where ν = 2 × NDETECT − 6</td></tr><tr><td>x, y, z (x, y, z; 3 columns)</td><td>au</td><td>ICRS-oriented positions</td></tr><tr><td>v<sub>x</sub>, v<sub>y</sub>, v<sub>z</sub> (v_x, v_y, v_z; 3 columns)</td><td>au yr<sup>−1</sup></td><td>ICRS-oriented velocities</td></tr><tr><td>Σ<sub>μ,ν</sub> (Sigma_mu_nu; 21 columns)</td><td>(au,au/yr)<sup>2</sup></td><td>μ,ν element of the state vector covariance matrix.</td></tr><tr><td>Class</td><td></td><td>Dynamical classification</td></tr><tr><td>Notes</td><td></td><td>Notes on the object<sup>a</sup></td></tr></tbody></table>

Notes. The description of each column is given here. All of the elements reported are barycentric and refer to epoch 2016.0, and their uncertainties correspond to the 1σ uncertainty marginalized over other orbital parameters.

<sup>a</sup> Insecure dynamical classifications include ins in this column, noncharacterized objects have nc, and resonant objects also have their {p, q, m, n, r, s} resonant arguments listed.

(This table is available in its entirety in FITS format.)

We compare this simulated DES classical TNO population to the 251 real DES Y6 objects in the same a range with $H_g < 8.5$. The tests presented here only use the normalized distributions in each parameter, so the difference in the total number of objects is not considered in the tests. Figure 10 shows the cumulative distributions in a, e, i, and $H_g$ expected from the L7 model in comparison to the observed DES detections. A Kolmogorov-Smirnov test (K-S test; e.g., Press et al. 2007) leads to p-values = {0.0257, 0.0382, 0.0046, $\le 10^{-5}$} for the {a, e, i, $H_g$} distributions, respectively. The observed i and $H_g$ distributions are incompatible with the CFEPS-L7 model, while the model a and e distributions are

<table><caption>**Table 4**
Number of Objects Per Dynamical Classification for the 814 Objects</caption><thead><tr><th>Dynamical Class</th><th>Number of Objects</th><th>Dynamical Class</th><th>Number of Objects</th></tr></thead><tbody><tr><td>Classical belt</td><td>381 (+1 insecure)</td><td>Scattering</td><td>50 (+1 insecure)</td></tr><tr><td>Detached</td><td>154</td><td>Centaur</td><td>1</td></tr><tr><td>Oort cloud</td><td>1</td><td></td><td></td></tr><tr><th>Mean-motion resonance</th><th>Number of objects</th><th>Mean-motion resonance</th><th>Number of objects</th></tr><tr><td>1:1 (30.1 au)</td><td>10</td><td>13:6 (50.9 au)</td><td>1</td></tr><tr><td>5:4 (34.9 au)</td><td>6</td><td>11:5 (50.4 au)</td><td>1</td></tr><tr><td>4:3 (36.3 au)</td><td>14</td><td>9:4 (51.7 au)</td><td>1</td></tr><tr><td>7:5 (37.7 au)</td><td>1</td><td>7:3 (52.9 au)</td><td>7</td></tr><tr><td>3:2 (39.4 au)</td><td>69</td><td>5:2 (55.4 au)</td><td>15</td></tr><tr><td>26:17 (39.9 au)</td><td>1</td><td>3:1 (62.6 au)</td><td>4</td></tr><tr><td>23:15 (40.0 au)</td><td>1</td><td>19:6 (64.9 au)</td><td>1 candidate</td></tr><tr><td>11:7 (40.7 au)</td><td>1</td><td>16:5 (65.4 au)</td><td>2</td></tr><tr><td>19:12 (40.9 au)</td><td>1 candidate</td><td>10:3 (67.1 au)</td><td>1</td></tr><tr><td>5:3 (42.3 au)</td><td>14+1 candidate</td><td>7:2 (69.4 au)</td><td>8</td></tr><tr><td>22:13 (42.7 au)</td><td>1 candidate</td><td>4:1 (75.8 au)</td><td>2</td></tr><tr><td>12:7 (43.1 au)</td><td>1</td><td>21:5 (78.4 au)</td><td>1</td></tr><tr><td>7:4 (43.7 au)</td><td>25+2 candidates</td><td>13:3 (80.0 au)</td><td>1</td></tr><tr><td>29:16 (44.8 au)</td><td>1</td><td>22:5 (80.8 au)</td><td>1 candidate</td></tr><tr><td>11:6 (45.1 au)</td><td>1+1 candidate</td><td>11:2 (93.8 au)</td><td>1</td></tr><tr><td>13:7 (45.5 au)</td><td>1 candidate</td><td>8:1 (120.4 au)</td><td>1</td></tr><tr><td>2:1 (47.7 au)</td><td>21</td><td>25:3 (123.7 au)</td><td>1</td></tr><tr><td>23:11 (49.2 au)</td><td>1</td><td>10:1 (139.7 au)</td><td>1 candidate</td></tr><tr><td>21:10 (49.3 au)</td><td>1 candidate</td><td></td><td></td></tr><tr><td>Total</td><td>814</td><td></td><td></td></tr></tbody></table>

Note. The resonant objects are presented in order of increasing semimajor axis, with the approximate value presented in parentheses.

only marginally rejectable by the data, which is not surprising considering the large increase in information on high-inclination TNOs in the DES data compared to the CFEPS data to which the L7 model was initially fit. Since all four distributions show a low p-value from this test, this indicates a significant deviation between the model and the data.

We investigate whether a simple reweighting of the L7 subpopulations can improve agreement with the DES observations. We reparameterize the L7 classical model by changing the fraction $f_μ$ of the stirred, hot, and kernel components μ as

$$ \frac{dn}{di} = f_s \frac{dn_s}{di} + f_h \frac{dn_h}{di} + f_k \frac{dn_k}{di} \quad (6) $$

with the constraint that $f_s+f_h+f_k=1$. A least-squares fit of $dn/di$ to the measured inclination distribution leads to $f_h = 0.653$, $f_s = 0.222$ and $f_k = 0.125$, compared to $f_h = 0.51$, $f_s = 0.38$ and $f_k = 0.11$ in the original CFEPS model.<fnref n="51" /> Repeating the K-S test with these new fractions leads to K-S test p-values = {0.0018, 0.7792, 0.0134, $\le 10^{-5}$}. The model reproduces the measured e distribution satisfactorily, but rejects the new a, i, and $H_g$ distributions. No simple subpopulation

<footnote n="51">⁵¹ Note that the stirred and kernel components share the same inclination distribution in the model, and so any difference in $dn/di$ comes from the selection functions for each component.</footnote>
<figure><img src="image_14.png" /></figure>

<figure><img src="image_15.png" /><figcaption>**Figure 9.** Left: absolute magnitude $H_r$ versus barycentric discovery distance for the 814 objects reported here. The blue dots correspond to objects new to this search, while the red dots correspond to those that were present in the Y4 catalog. The dotted curve shows the $r_{50} = 23.3$ of the Y4 search, and the solid black curve shows the $r_{50} = 23.77$ of this search. The absolute magnitudes are found by taking the mean flux in each exposure corrected to a nominal geocentric distance $d$ and heliocentric distance $\Delta$. Right: histogram of apparent magnitudes $m_r$ for the objects found in this search. The solid purple histogram shows the full sample, while the lined histograms show the sample found in the Y4 search and the sample new to this Y6 search. The black lines (dotted and solid) are the same as in the left panel. While the majority of the new objects were discovered at $m_r > r_{50, Y4}$, 56 were found at magnitudes brighter than $r_{50, Y4}$.</figcaption></figure>

<figure><img src="image_16.png" /></figure>

<figure><img src="image_17.png" /><figcaption>**Figure 10.** Comparison between the CFEPS-L7 classical objects with $40 \text{ au} < a < 47 \text{ au}$ simulated into the DES Y6 and the DES-detected classical TNOs in this range. Each panel presents the cumulative distribution in {$a, e, i, H_g$}, with the blue lines representing the baseline model and the red lines representing the reweighted models. The black line corresponds to the DES observations.</figcaption></figure>

reweighting can bring the L7 model into agreement with the observations. This result is in agreement with the CFEPS-HiLat result from Petit et al. (2017), and so the DES Y6 catalog, then, will enable considerable refinement of TNO population models, particularly in combination with the OSSOS data.

Repeating these tests with shallower $H_g < 8.0$ samples (351 simulated DES detections, 130 real ones)—this being the limit of the CFEPS data to which the model was fit—yields consistency between the data and the model. The p-values are {$0.079, 0.338, 0.308, 0.067$} in $a, e, i$, and $H_g$.

## 5.2. Isotropy of Extreme TNOs Revisited

We repeat the test of eTNOs' isotropy presented in Bernardinelli et al. (2020b), which asked whether the population of TNOs with $a > 150 \text{ au}$ and $q > 30 \text{ au}$ shows deviation from an underlying uniform distribution in $\omega, \Omega$, and $\varpi \equiv \Omega + \omega$. Deviations from such isotropy serve as the primary motivation for the Planet 9 hypothesis (Batygin & Brown 2016; Batygin et al. 2019). The underlying population model chosen

for these objects is such that we clone the detected {$a_j, e_j, i_j, H_j$} for each object $j$, and randomize their $\omega, \Omega$, and mean anomaly $M$. DES Y6 observations of these objects are then simulated as in Section 3, and we compute the probabilities $p(\theta|s)$ of detecting an object with angle $\theta \in \{\omega, \Omega, M\}$ conditioned on a successful detection $s$. We compare $p(\theta|s)$ to the empirical distributions coming from the true detected eTNOs using Kuiper's variant of the K-S test (Kuiper 1960), as described in Bernardinelli et al. (2020b).

We present here the results for the population with $a > 150 \text{ au}$ and $q > 30 \text{ au}$ (16 objects), as in Shankman et al. (2017) and Case 1 of Bernardinelli et al. (2020b), as well as a more restrictive $a > 230 \text{ au}$ selection (nine objects), as in Brown & Batygin (2019) and Napier et al. (2021), where the clustering signal is supposed to be stronger.<fnref n="52" /> The selection functions and the angles for the detected objects for this more restrictive selection are shown in Figure 11. We obtain, for

<footnote n="52">⁵² We note that Brown & Batygin (2021) present a similar test on objects with $a > 150 \text{ au}$ and $q > 42 \text{ au}$.</footnote>
<figure><img src="image_18.png" /><figcaption>Figure 11. Histograms of the relative detection probabilities for the a > 230 au eTNOs (in blue) for ω, Ω, and ϖ constructed as in Equation (1) of Bernardinelli et al. (2020b). The black vertical dashed lines represent the detected objects.</figcaption></figure>

{ω, Ω, ϖ}, p-values = {0.9224, 0.9514, 0.8788} for the a > 150 au sample, and p-values = {0.1355, 0.0144, 0.5851} for the a > 230 au sample. The a > 150 au sample is consistent with an underlying isotropic distribution, that is, the null hypothesis that the observed values come from the parent isotropic distribution cannot be rejected. The a > 230 au case shows poorer agreement with isotropy in Ω—but note that the probability of a truly isotropic population yielding a p-value = 0.0144 among six uncorrelated tests (see discussion in Bernardinelli et al. 2020b) is 1 − (1 − 0.0144)⁶ = 8% (or 4% considering only the 3 a > 230 au tests). If the tests are correlated, the odds of a chance occurrence of one such low p-value are reduced. Furthermore, the isotropic hypothesis is fully acceptable in ϖ, the variable in which the clustering signal is supposed to be the strongest. We repeat the f-test of Bernardinelli et al. (2020b), where the total likelihood of detection for all objects is compared to the likelihood of an ensemble of clones drawn from each object's selection function. We have that f(Ω) = 0.287; that is, 28.7% of the sets of clones yield a lower likelihood than one measured for the detected objects. This test indicates that the Kuiper statistics are not being driven by individual highly unlikely points. This result, then, remains in agreement with Shankman et al. (2017), Bernardinelli et al. (2020b), and Napier et al. (2021): the apparent clustering in orbital element space of these objects is consistent with the selection functions of the surveys and does not conclusively demand a massive perturber or other mechanism to break the eTNO isotropy in the DES survey. There is, however, a suggestive tendency for avoidance of − 90° < Ω < 90° in the DES sample, at 4%–8% significance, which motivates continuation of this exercise with the Vera C. Rubin Observatory Large Synoptic Survey Telescope or other future data.

5.3. Resonant Dropouts

We test the tendency of high-q objects to be found preferentially “sunward” of distant resonances (a lower than the nominal resonance center) rather than “outward” (a higher than the nominal resonance center; Kaib & Sheppard 2016; Nesvorny et al. 2016; Pike & Lawler 2017). Such asymmetry is expected in models where “resonance sweeping” is carrying TNOs outward as Neptune migrates, but objects drop out of the resonance en route. This process is sensitive to the smoothness and speed of Neptune’s migration (see detailed discussion in

Kaib & Sheppard 2016). We conduct a similar test to the one presented in Lawler et al. (2019), who found marginal inconsistency (p-value ≈ 0.029) of the OSSOS TNOs near the 5:2 and 4:1 resonances with q > 40 au arising from uniform distribution in period P ≡ a³/².

Following the notation introduced in Bernardinelli et al. (2020b), we produce a model for the high-q TNOs that is uniform in ω, Ω, and M and predicts a population p(e, i, H|s) conditioned on a successful detection s that reproduces the observed distribution in {e, i, H}:

$$
\begin{equation}
\begin{split}
p(a, e, i, \Omega, \omega, \mathcal{M}, H) \\
& \propto \sum_j \frac{\delta(e - e_j)\delta(i - i_j)\delta(H - H_j)}{p(s|e_j, i_j, H_j)} \\
& \times u(\Omega)u(\omega)u(\mathcal{M})p(a),
\end{split}
\tag{7}
\end{equation}
$$

where δ is the Dirac delta function, u(θ) is a uniform distribution in θ ∈ [0, 2π), p(a) is a chosen distribution in a, and the sum is carried over the detected objects j. We simulate these detections following the procedure defined in Section 3, and obtain an ensemble of detected TNOs, where each object is weighted by the fraction of successful detections p(s|e, i, H).

We select the 16 nonresonant, nonscattering objects with q > 38 au (so a larger sample can be obtained) within ± 2 au of the 5:2 (two objects, all sunward), 3:1 (seven objects, five sunward), 7:2 (three objects, all sunward), and 4:1 (four objects, all sunward) resonances, and test whether the sample is consistent with a distribution uniform in P (as in Lawler et al. 2019), as well as a distribution p(a) ∝ a − 5 / 2 (as in the hot and stirred components of the CFEPS model). We define

$$
\Delta P \equiv \frac{P - P_{\text{res}}}{P_{\text{N}}}, \qquad (8)
$$

that is, the fractional difference between an object's period and the nominal resonance location in units of Neptune's period. Objects with ΔP < 0 (ΔP > 0) are sunward (outward) of the resonance. We apply a K-S test between the observed distribution in ΔP and the distribution coming from the simulated ensemble with smooth p(a), as shown in Figure 12. The p-value of the K-S statistic is determined by sampling 10⁶ different realizations of the smooth-p(a) reference distribution. As an even simpler statistic, we calculate the expected fraction f of the detected ± 2 au population expected to be sunward under
<figure><img src="image_19.png" /><figcaption>**Figure 12.** Selection functions in fractional period difference $\Delta P$ for the uniform (left) and $a^{-5/2}$ (right) models for the high (blue lines) and low (red lines) perihelion sample of objects near the *n*:1 and *n*:2 resonances. The solid lines show the cumulative fraction of detections in each population, and the dashed lines represent the detections expected from the model distribution. The solid black line divides the “sunward” ($\Delta P < 0$) from the “outward” ($\Delta P > 0$) detections. The distinction between these two groups is visually apparent, with the high-$q$ curves deviating significantly from the expected distribution under the assumed model.</figcaption></figure>

<figure><img src="image_20.png" /></figure>

the smooth-$p(a)$ reference model, and use binomial statistics to calculate the probability of matching or exceeding the 14 out of 16 real objects that are sunward.

Both the uniform-period and smooth power-law $p(a)$ distributions can be rejected at very high significance for the first time, with $p = 0.0003$, 0.0007 or $p = 0.0008$, 0.0010 from the K-S (indicating deviation from uniformity) and binomial (indicating a preference for sunward detections) tests, respectively, for this set of objects near the *n*:1 and *n*:2 resonances. We repeat these tests for the $q < 38$ au nonresonant, non-scattering objects near the 5:2 (17 objects, nine sunward), 3:1 (eight objects, all sunward), 7:2 (six objects, two sunward), and 4:1 (two objects, zero sunward) resonances. Dynamical expectations are that perturbations from Neptune should weaken or erase any signature of resonant dropouts in this population (Kaib & Sheppard 2016). The 33 low-$q$ detections indeed show no evidence of preference for being sunward of resonance, having $p$-values of 0.92 (uniform and smooth power law) and 0.58 (uniform), 0.69 (smooth power law) in the K-S and binomial statistics, respectively.

These statistics confirm the clear visual impression from the upper panel of Figure 8 that in fact nearly all nonresonant TNOs beyond the 2:1 resonance with $q \gtrsim 40$ are found just sunward of *n*:1 or *n*:2 resonances.

## 6. Summary

We describe a catalog of 814 objects found in a complete search of the 6 yr, 5000 deg² coverage DES, yielding a sample 50% complete at $m_r \approx 23.8$. All detected objects have astrometry over multiyear arcs tied to Gaia DR2, and photometric measurements in the grizY bands, with all measurements at the shot-noise limits. We demonstrate using both synthetic detections as well as previously discovered objects that the survey is highly complete given our thresholds of NUNIQUE $\ge 7$, ARCCUT $> 6$ months, and heliocentric distance $> 29$ au. The STS technique discards virtually all accidental linkages, yielding a powerful methodology for discovery confirmation of faint sources. The Y6 sample is complete to $r_{50} = 23.77$, a gain of ~0.5 mag over the Y4 sample, leading to the recovery of 506 objects that were not found in Y4, all but

45 of which were new discoveries. We believe we have reached the limits of the DES data for detection of solar system objects at $d > 29$ au. A search for closer objects is currently beyond our computational means, as the burden grows as a high inverse power of $d$, and the search reported here required an estimated 15–20 million CPU hours. It could be made feasible by reducing the transient catalog density significantly, e.g., by imposing a minimum S/N. This release includes only the *r*-band photometry—a future publication will describe our dedicated TNO photometric processing and report detailed photometry for all bands and exposures.

Searching for TNOs in the DES images presents challenges since it was not optimized for this purpose, such as its temporally sparse data, multiband survey coverage, and extensive coverage at high ecliptic latitudes with lower TNO densities. Despite this, the DES sample is comparable in a number of objects to the largest predecessor TNO-targeted surveys. Its well-characterized observational biases and publicly available survey simulator lead to a sample suited for statistical tests of solar system population models. Indeed the fact the DES footprint is *not* concentrated at the ecliptic will make the combination of the DES and OSSOS particularly powerful for constraining TNO population models. The DES catalog represents a significant increase in the number of known dynamically detached objects, including new extreme and high-$q$ TNOs, as well as objects in distant resonances or high inclinations, which are relevant to more detailed hypotheses of the formation of the outer solar system.

The data set presented here also contains a number of new notable objects, in particular the Oort cloud comet C/2014 UN271 (Bernardinelli-Bernstein), as well as the detached object 2014 YX$_{91}$ with $i = 62^\circ 8$ and the $q = 54.3$ au object 2014 US$_{277}$.

Initial applications of the statistical power of the DES Y6 TNO data are presented in three results for distinct populations of the trans-Neptunian region:

1. We show that our population of classical TNOs is formally incompatible with the CFEPS-L7 model for the $H_g < 8.5$ limit even if we change the relative sizes of the hot, stirred, and kernel components. Such results indicate
suggesting further substructure and/or a wider inclination distribution for the classical population, in agreement with recent results (e.g., Petit et al. 2017).

2. The sample of extreme TNOs measured by DES is consistent with an underlying uniform distribution in $\omega$ and $\varpi$, agreeing with previous results from OSSOS (Shankman et al. 2017), DES, (Bernardinelli et al. 2020b), and their combination with the Sheppard-Trujillo surveys (Napier et al. 2021). There is a low significance disagreement between the observed $\Omega$ distribution and uniformity, and future surveys with large sky coverage, such as the upcoming Rubin Observatory's Legacy Survey of Space and Time (Ivezić et al. 2019), will enable these two scenarios to be tested even further.

3. We measure with high significance the preference of high-$q$ objects to be found sunward of $n:1$ and $n:2$ resonances beyond 50 au, as expected from solar system formation models with grainy or slow migration.

The DES data, then, hold a unique set of objects comprising ≈20% of all currently known TNOs, with well-characterized observational biases and covering ≈1/8 of the sky. These will be valuable for further detailed statistical tests of formation models for the trans-Neptunian region.

University of Pennsylvania authors have been supported in this work by grants AST-1515804 and AST-2009210 from the National Science Foundation, and grant DE-SC0007901 from the Department of Energy.

Funding for the DES Projects has been provided by the U.S. Department of Energy, the U.S. National Science Foundation, the Ministry of Science and Education of Spain, the Science and Technology Facilities Council of the United Kingdom, the Higher Education Funding Council for England, the National Center for Supercomputing Applications at the University of Illinois at Urbana-Champaign, the Kavli Institute of Cosmological Physics at the University of Chicago, the Center for Cosmology and Astro-Particle Physics at the Ohio State University, the Mitchell Institute for Fundamental Physics and Astronomy at Texas A&M University, Financiadora de Estudos e Projetos, Fundação Carlos Chagas Filho de Amparo à Pesquisa do Estado do Rio de Janeiro, Conselho Nacional de Desenvolvimento Científico e Tecnológico and the Ministério da Ciência, Tecnologia e Inovação, the Deutsche Forschungsgemeinschaft, and the Collaborating Institutions in the Dark Energy Survey.

The Collaborating Institutions are Argonne National Laboratory, the University of California at Santa Cruz, the University of Cambridge, Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas-Madrid, the University of Chicago, University College London, the DES-Brazil Consortium, the University of Edinburgh, the Eidgenössische Technische Hochschule (ETH) Zürich, Fermi National Accelerator Laboratory, the University of Illinois at Urbana-Champaign, the Institut de Ciències de l'Espai (IEEC/CSIC), the Institut de Física d'Altes Energies, Lawrence Berkeley National Laboratory, the Ludwig-Maximilians Universität München and the associated Excellence Cluster Universe, the University of Michigan, NSF's NOIRLab, the University of Nottingham, The Ohio State University, the University of Pennsylvania, the University of Portsmouth, SLAC National Accelerator Laboratory, Stanford University, the University of

Sussex, Texas A&M University, and the OzDES Membership Consortium.

Based in part on observations at Cerro Tololo Inter-American Observatory at NSF's NOIRLab (NOIRLab Prop. ID 2012B-0001; PI: J. Frieman), which is managed by the Association of Universities for Research in Astronomy (AURA) under a cooperative agreement with the National Science Foundation.

The DES data management system is supported by the National Science Foundation under grant Nos. AST-1138766 and AST-1536171. The DES participants from Spanish institutions are partially supported by MICINN under grants ESP2017-89838, PGC2018-094773, PGC2018-102021, SEV-2016-0588, SEV-2016-0597, and MDM-2015-0509, some of which include ERDF funds from the European Union. IFAE is partially funded by the CERCA program of the Generalitat de Catalunya. Research leading to these results has received funding from the European Research Council under the European Union's Seventh Framework Program (FP7/2007-2013) including ERC grant agreements 240672, 291329, and 306478. We acknowledge support from the Brazilian Instituto Nacional de Ciência e Tecnologia (INCT) do e-Universo (CNPq grant 465376/2014-2).

This manuscript has been authored by Fermi Research Alliance, LLC under Contract No. DE-AC02-07CH11359 with the U.S. Department of Energy, Office of Science, Office of High Energy Physics.

**Software:** The software developed in this work will be made public shortly after the publication. This work made use of the following public codes: NUMPY (Oliphant 2015), SCIPY (Jones et al. 2001), ASTROPY (Astropy Collaboration 2013, 2018), MATPLOTLIB (Hunter 2007), IPYTHON (Pérez & Granger 2007), EASYACCESS (Carrasco Kind et al. 2019), WCSFIT and PIXMAPPY (Bernstein et al. 2017a), SExTRACTOR (Bertin & Arnouts 1996), CFITSIO (Pence 1999), EIGEN (Guennebaud & Jacob 2010), CSPICE (Acton 1996; Acton et al. 2018), VAEX (Breddels & Veljanoski 2018).

## Appendix

### Changes to the Linking Algorithm

We have updated the triplet finding algorithm presented in Paper I. Here, we will present only the details required to understand the changes; we refer the readers to Bernstein & Khushalani (2000) and Paper I for the full details.

Following the formalism of Bernstein & Khushalani (2000), we assume a coordinate system whose z-axis points to the center of a field that contains the object, with the origin being the location of the observatory at the midpoint of the observing season ($t_0$). We use the six orbital parameters {$\alpha = x/z, \beta = y/z, \gamma = 1/z, \dot{\alpha} = \dot{x}/z, \dot{\beta} = \dot{y}/z, \dot{\gamma} = \dot{z}/z$}. Under the approximation of inertial motion for the TNO (which is adequate for this purpose), the observed angular coordinates of the TNO at time $t$ are

$$ \theta(t) = \frac{1}{1 + \dot{\gamma}(t - t_0) - \gamma z_{\text{obs}}(t)} \times (\alpha + \dot{\alpha}(t - t_0) - \gamma x_{\text{obs}}(t), \beta + \dot{\beta}(t - t_0) - \gamma y_{\text{obs}}(t)). \quad (A1) $$

Nominal orbital parameters for a given pair are determined by the four observed coordinates of the pair, plus a choice of nominal (inverse) distance $\gamma_0$ and an assignment $\dot{\gamma}_0 = 0$ to the line-of-
sight velocity. For a given exposure being searched for a potential third detection, this defines a nominal search position $θ_0$. We define two (nonorthogonal) vectors on the sky that track deviations of the predicted position as we vary the two degrees of freedom $γ$ and $\dot{γ}$ while keeping $α, β, \dot{α},$ and $\dot{β}$ constant:

$$ \mathbf{u} = \delta\gamma \frac{\partial\theta}{\partial\gamma}, \qquad (A2) $$

$$ \mathbf{v} = \dot{\gamma}_{\text{bind}} \frac{\partial\theta}{\partial\dot{\gamma}}. \qquad (A3) $$

These vectors span the maximum deviation of the pair's third detection from $θ_0$ if the TNO is assumed to lie in an (inverse) distance bin bounded by $γ_0 \pm \delta\gamma$ and have line-of-sight velocity bounded by $|\dot{\gamma}| \le \dot{\gamma}_{\text{bind}}$. As the name suggests, $\dot{\gamma}_{\text{bind}}$ is chosen to be the maximum allowed for a bound orbit given $\dot{α}, \dot{β}$, and $γ_0$.

Now the position of any potential third detection of the pair in this exposure can be written as

$$ \theta = \theta_0 + \mu\mathbf{u} + \nu\mathbf{v}, \qquad (A4) $$

with $|\mu| < 1, |\nu| < 1$. Conversely, each detection's position can be mapped to a pair $(\mu, \nu)$ by a linear transformation. To test whether a detection is inside the parallelogram defined by a bound orbit within the distance bin, it suffices to ask whether $|\mu|, |\nu| \le 1$. This simple linear test leads to a speed gain of $100 \times$, on average, over the kD-tree triplet search of Paper I. Another benefit of the parallelogram search is a much simpler, parallelized implementation using NUMBA (Lam et al. 2015).

## ORCID iDs

Pedro H. Bernardinelli https://orcid.org/0000-0003-0743-9422

Gary M. Bernstein https://orcid.org/0000-0002-8613-8259

Masao Sako https://orcid.org/0000-0003-2764-7093

Brian Yanny https://orcid.org/0000-0002-9541-2678

M. Aguena https://orcid.org/0000-0001-5679-6747

E. Bertin https://orcid.org/0000-0002-3602-3664

J. Carretero https://orcid.org/0000-0002-3130-0204

C. Conselice https://orcid.org/0000-0003-1949-7638

M. Costanzi https://orcid.org/0000-0001-8158-1449

H. T. Diehl https://orcid.org/0000-0002-8357-7467

K. Eckert https://orcid.org/0000-0002-1407-4700

I. Ferrero https://orcid.org/0000-0002-1295-1132

P. Fosalba https://orcid.org/0000-0002-1510-5214

J. García-Bellido https://orcid.org/0000-0002-9370-8360

D. W. Gerdes https://orcid.org/0000-0001-6942-2736

D. Gruen https://orcid.org/0000-0003-3270-7644

R. A. Gruendl https://orcid.org/0000-0002-4588-6517

J. Gschwend https://orcid.org/0000-0003-3023-8362

S. R. Hinton https://orcid.org/0000-0003-2071-9349

D. L. Hollowood https://orcid.org/0000-0002-9369-4157

S. Kent https://orcid.org/0000-0003-4207-7420

K. Kuehn https://orcid.org/0000-0003-0120-0808

N. Kuropatkin https://orcid.org/0000-0003-2511-0946

F. Menanteau https://orcid.org/0000-0002-1372-2534

R. Miquel https://orcid.org/0000-0002-6610-4836

J. Myles https://orcid.org/0000-0001-6145-5859

R. L. C. Ogando https://orcid.org/0000-0003-2120-1154

A. Palmese https://orcid.org/0000-0002-6011-0530

F. Paz-Chinchón https://orcid.org/0000-0003-1339-2683

A. Pieres https://orcid.org/0000-0001-9186-6042

A. K. Romer https://orcid.org/0000-0002-9328-879X

A. Roodman https://orcid.org/0000-0001-5326-3486

S. Serrano https://orcid.org/0000-0002-0211-2861

I. Sevilla-Noarbe https://orcid.org/0000-0002-1831-1953

M. Soares-Santos https://orcid.org/0000-0001-6082-8529

G. Tarle https://orcid.org/0000-0003-1704-0781

A. R. Walker https://orcid.org/0000-0002-7123-8943

## References

Acton, C., Bachman, N., Semenov, B., & Wright, E. 2018, P&SS, 150, 9

Acton, C. H. 1996, P&SS, 44, 65

Adams, F. C. 2010, ARA&A, 48, 47

Aihara, H., Armstrong, R., Bickerton, S., et al. 2018, PASJ, 70, 8

Alexandersen, M., Benecchi, S. D., Chen, Y.-T., et al. 2019, ApJS, 244, 19

Allen, R. L., Bernstein, G. M., & Malhotra, R. 2001, ApJ, 549, L241

Allen, R. L., Bernstein, G. M., & Malhotra, R. 2002, AJ, 124, 2949

Astropy Collaboration 2013, A&A, 558, A33

Astropy Collaboration 2018, AJ, 156, 123

Banda-Huarca, M. V., Camargo, J. I. B., Desmars, J., et al. 2019, AJ, 157, 120

Bannister, M. T. 2020, in The Trans-Neptunian Solar System, ed. D. Prialnik, M. A. Barucci, & L. Young (Amsterdam: Elsevier), 439

Bannister, M. T., Gladman, B. J., Kavelaars, J. J., et al. 2018, ApJS, 236, 18

Bannister, M. T., Kavelaars, J. J., Petit, J.-M., et al. 2016, AJ, 152, 70

Batygin, K., Adams, F. C., Brown, M. E., & Becker, J. C. 2019, PhR, 805, 1

Batygin, K., & Brown, M. E. 2016, AJ, 151, 22

Batygin, K., Brown, M. E., & Betts, H. 2012, ApJL, 744, L3

Becker, J. C., Khain, T., Hamilton, S. J., et al. 2018, AJ, 156, 81

Bernardinelli, P. H., Bernstein, G. M., Montet, B. T., et al. 2021, ApJL, 921, L37

Bernardinelli, P. H., Bernstein, G. M., Sako, M., et al. 2020a, ApJS, 247, 32

Bernardinelli, P. H., Bernstein, G. M., Sako, M., et al. 2020b, PSJ, 1, 28

Bernstein, G., & Khushalani, B. 2000, AJ, 120, 3323

Bernstein, G. M., Abbott, T. M. C., Desai, S., et al. 2017b, PASP, 129, 114502

Bernstein, G. M., Armstrong, R., Plazas, A. A., et al. 2017a, PASP, 129, 074503

Bernstein, G. M., Trilling, D. E., Allen, R. L., et al. 2004, AJ, 128, 1364

Bernstein, J. P., Kessler, R., Kuhlmann, S., et al. 2012, ApJ, 753, 152

Bertin, E., & Arnouts, S. 1996, A&AS, 117, 393

Breddels, M. A., & Veljanoski, J. 2018, A&A, 618, A13

Brown, M. E. 2012, AREPS, 40, 467

Brown, M. E., & Batygin, K. 2019, AJ, 157, 62

Brown, M. E., & Batygin, K. 2021, AJ, 162, 219

Brown, M. E., Bannister, M. T., Schmidt, B. P., et al. 2015, AJ, 149, 69

Burke, D. L., Rykoff, E. S., Allam, S., et al. 2017, AJ, 155, 41

Carrasco Kind, M., Drlica-Wagner, A., Koziol, A., & Petravick, D. 2019, JOSS, 4, 1022

Chen, Y.-T., Lin, H.-W., Alexandersen, M., et al. 2018, PASJ, 70, S38

Dawson, R. I., & Murray-Clay, R. 2012, ApJ, 750, 43

Diehl, H. T., Neilsen, E., Gruendl, R., et al. 2016, Proc. SPIE, 9910, 99101D

Diehl, H. T., Neilsen, E., Gruendl, R. A., et al. 2018, Proc. SPIE, 10704, 10704

Duncan, M. J., Brasser, R., Dones, L., & Levison, H. F. 2008, in The Solar System Beyond Neptune, ed. M. A. Barucci et al. (Tucson, AZ: Univ. of Arizona Press), 315

Elliot, J. L., Kern, S. D., Clancy, K. B., et al. 2005, AJ, 129, 1117

Fernández, J., & Ip, W.-H. 1984, Icar, 58, 109

Flaugher, B., Diehl, H. T., Honscheid, K., et al. 2015, AJ, 150, 150

Fortino, W. F., Bernstein, G. M., Bernardinelli, P. H., et al. 2021, AJ, 162, 106

Fraser, W. C., Brown, M. E., & Schwamb, M. E. 2010, Icar, 210, 944

Fuentes, C. I., & Holman, M. J. 2008, AJ, 136, 83

Gaia Collaboration, Brown, A. G. A., Vallenari, A., et al. 2018, A&A, 616, A1

Gerdes, D. W., Jennings, R. J., Bernstein, G. M., et al. 2016, AJ, 151, 39

Gerdes, D. W., Sako, M., Hamilton, S., et al. 2017, ApJ, 839, L15

Gladman, B., Kavelaars, J. J., Nicholson, P. D., Loredo, T. J., & Burns, J. A. 1998, AJ, 116, 2042

Gladman, B., Marsden, B. G., & Vanlaerhoven, C. 2008, in The Solar System Beyond Neptune, ed. M. A. Barucci et al. (Tucson, AZ: Univ. of Arizona Press), 43

Gladman, B., Lawler, S. M., Petit, J. M., et al. 2012, AJ, 144, 23

Gladman, B. J., Davis, D. R., Neese, C., et al. 2009, Icar, 202, 104

Gomes, R. S. 2003, Icar, 161, 404

Guennebaud, G., & Jacob, B. 2010, Eigen, v3, https://eigen.tuxfamily.org

Henghes, B., Lahav, O., Gerdes, D. W., et al. 2021, PASP, 133, 014501
Hunter, J. D. 2007, CSE, 9, 90

Ivezić, Ž, Kahn, S. M., Tyson, J. A., et al. 2019, ApJ, 873, 111

Jewitt, D., & Luu, J. 1993, Natur, 362, 730

Jewitt, D. C., & Luu, J. X. 1995, AJ, 109, 1867

Jones, E., Oliphant, T., Peterson, P., et al. 2001, SciPy: Open source scientific tools for Python, http://www.scipy.org/

Jones, R. L., Gladman, B., Petit, J.-M., et al. 2006, Icar, 185, 508

Kaib, N. A., & Sheppard, S. S. 2016, AJ, 152, 133

Kavelaars, J. J., Jones, R. L., Gladman, B. J., et al. 2009, AJ, 137, 4917

Khain, T., Becker, J. C., Adams, F. C., et al. 2018, AJ, 156, 273

Khain, T., Becker, J. C., Wen Lin, H., et al. 2020, AJ, 159, 133

Kuiper, N. H. 1960, Nederl. Akad. Wetensch. Proc. Ser. A, 63, 38

Lam, S. K., Pitrou, A., & Seibert, S. 2015, Proc. of the Second Workshop on the LLVM Compiler Infrastructure in HPC, LLVM '15 (New York: Association for Computing Machinery)

Lawler, S. M., Kavelaars, J. J., Alexandersen, M., et al. 2018, FrASS, 5, 14

Lawler, S. M., Pike, R. E., Kaib, N., et al. 2019, AJ, 157, 253

Levison, H. F., Morbidelli, A., VanLaerhoven, C., Gomes, R., & Tsiganis, K. 2008, Icar, 196, 258

Lin, H. W. W., Gerdes, D. J., Hamilton, S., et al. 2019, Icar, 321, 426

Malhotra, R. 1993, Natur, 365, 819

Malhotra, R. 1995, AJ, 110, 420

Millis, R. L., Buie, M. W., Wasserman, L. H., et al. 2002, AJ, 123, 2083

Morganson, E., Gruendl, R. A., Menanteau, F., et al. 2018, PASP, 130, 074501

Napier, K. J., Gerdes, D. W., Lin, H. W., et al. 2021, PSJ, 2, 59

Neilsen, E., & Annis, J. 2013, in ASP Conf. Ser. 485, Astronomical Data Analysis Software and Systems XXIII, ed. N. Manset & P. Forshay (San Francisco, CA: ASP), 77

Neilsen, E. H., Annis, J. T., Diehl, H. T., et al. 2019, Dark Energy Survey's Observation Strategy, Tactics, and Exposure Scheduler, Fermilab Report, FERMILAB-TM-2714-AE-CD-PPD

Nesvorný, D. 2015, AJ, 150, 73

Nesvorný, D., Vokrouhlický, D., & Roig, F. 2016, ApJL, 827, L35

Oliphant, T. E. 2015, Guide to NumPy (2nd ed.; USA: CreateSpace Independent Publishing Platform)

Pence, W. 1999, in ASP Conf. Ser. 172, Astronomical Data Analysis Software and Systems VIII, ed. D. M. Mehringer, R. L. Plante, & D. A. Roberts (San Francisco, CA: ASP), 487

Petit, J. M., Kavelaars, J. J., Gladman, B. J., et al. 2011, AJ, 142, 131

Petit, J.-M., Kavelaars, J. J., Gladman, B. J., et al. 2017, AJ, 153, 236

Pfalzner, S., Bhandare, A., Vincke, K., & Lacerda, P. 2018, ApJ, 863, 45

Pike, R. E., & Lawler, S. M. 2017, AJ, 154, 171

Pike, R. E., Fraser, W. C., Schwamb, M. E., et al. 2017, AJ, 154, 101

Press, W., Teukolsky, S., Vetterling, W., & Flannery, B. 2007, Numerical Recipes: The Art of Scientific Computing (3rd ed.; Cambridge: Cambridge Univ. Press)

Pérez, F., & Granger, B. E. 2007, CSE, 9, 21

Rabinowitz, D., Schwamb, M. E., Hadjiyska, E., & Tourtellotte, S. 2012, AJ, 144, 140

Rein, H., & Liu, S. F. 2012, A&A, 537, A128

Rein, H., & Tamayo, D. 2015, MNRAS, 452, 376

Schwamb, M. E., Brown, M. E., Rabinowitz, D. L., & Ragozzine, D. 2010, ApJ, 720, 1691

Shankman, C., Kavelaars, J. J., Bannister, M. T., et al. 2017, AJ, 154, 50

Sheppard, S. S., Trujillo, C., & Tholen, D. J. 2016, ApJ, 825, L13

Sheppard, S. S., Trujillo, C. A., Tholen, D. J., & Kaib, N. 2019, AJ, 157, 139

The Dark Energy Survey Collaboration 2005, arXiv:astro-ph/0510346

The Dark Energy Survey Collaboration 2016, MNRAS, 460, 1270

The Dark Energy Survey Collaboration 2018a, PhRvD, 98, 043526

The Dark Energy Survey Collaboration 2018b, MNRAS, 480, 3879

The Dark Energy Survey Collaboration 2019a, MNRAS, 483, 4866

The Dark Energy Survey Collaboration 2019b, ApJ, 872, L30

The Dark Energy Survey Collaboration 2019c, PhRvL, 122, 171301

The Dark Energy Survey Collaboration 2021, ApJS, 255, 20

Trujillo, C. A., Jewitt, D. C., & Luu, J. X. 2001, AJ, 122, 457

Tsiganis, K., Gomes, R., Morbidelli, A., & Levison, H. F. 2005, Natur, 435, 459

Volk, K., Murray-Clay, R. A., Gladman, B. J., et al. 2018, AJ, 155, 260

Weryk, R. J., Lilly, E., Chastel, S., et al. 2016, arXiv:1607.04895

Whidden, P. J., Kalmbach, J. B., Connolly, A. J., et al. 2019, AJ, 157, 119

Wisdom, J., & Holman, M. 1991, AJ, 102, 1528