# The PAndAS View of the Andromeda Satellite System.
## III. Dwarf Galaxy Detection Limits

Amandine Doliva-Dolinsky, Nicolas F. Martin, Guillaume F. Thomas, Annette M. N. Ferguson, Rodrigo A. Ibata, Geraint F. Lewis, Dougal Mackey, Alan W. Mcconnachie, Zhen Yuan

► To cite this version:

Amandine Doliva-Dolinsky, Nicolas F. Martin, Guillaume F. Thomas, Annette M. N. Ferguson, Rodrigo A. Ibata, et al.. The PAndAS View of the Andromeda Satellite System. III. Dwarf Galaxy Detection Limits. The Astrophysical Journal, 2022, 933, 10.3847/1538-4357/ac6fd5. insu-03764635

HAL Id: insu-03764635

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

Submitted on 1 Sep 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
OPEN ACCESS

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

CrossMark

# The PAndAS View of the Andromeda Satellite System. III. Dwarf Galaxy Detection Limits

Amandine Doliva-Dolinsky¹, Nicolas F. Martin¹,², Guillaume F. Thomas³,⁴, Annette M. N. Ferguson⁵, Rodrigo A. Ibata¹, Geraint F. Lewis⁶, Dougal Mackey⁷, Alan W. McConnachie⁸, and Zhen Yuan¹

¹Université de Strasbourg, CNRS, Observatoire astronomique de Strasbourg, UMR 7550, F-67000, France; amandine.doliva-dolinsky@astro.unistra.fr

²Max-Planck-Institut für Astronomie, Königstuhl 17, D-69117, Heidelberg, Germany

³Instituto de Astrofísica de Canarias, Calle Vía Láctea s/n, E-38206 La Laguna, Tenerife, Spain

⁴Universidad de La Laguna, Avda. Astrofísico Fco. Sánchez, E-38205 La Laguna, Tenerife, Spain

⁵Institute for Astronomy, University of Edinburgh Royal Observatory, Blackford Hill, Edinburgh EH9 3HJ, UK

⁶Sydney Institute for Astronomy, School of Physics, A28, The University of Sydney, NSW 2006, Australia

⁷Research School of Astronomy and Astrophysics, Australian National University, Canberra, ACT 2611, Australia

⁸NRC Herzberg Astronomy and Astrophysics, 5071 West Saanich Road, Victoria, BC V9E 2E7, Canada

Received 2021 October 1; revised 2022 April 29; accepted 2022 May 3; published 2022 July 8

## Abstract

We determine the detection limits of the search for dwarf galaxies in the Pan-Andromeda Archaeological Survey (PAndAS) using the algorithm developed by the PAndAS team. The recovery fractions of artificial dwarf galaxies are, as expected, a strong function of physical size and luminosity and, to a lesser extent, distance. We show that these recovery fractions vary strongly with location in the surveyed area because of varying levels of contamination from both the Milky Way foreground stars and the stellar halo of Andromeda. We therefore provide recovery fractions that are a function of size, luminosity, and location within the survey on a scale of ~1 × 1 deg² (or ~14 × 14 kpc²). Overall, the effective surface brightness for a 50% detection rate ranges between 28 and 30 mag arcsec⁻². This is in line with expectations for a search that relies on photometric data that are as deep as the PAndAS survey. The derived detection limits are an essential ingredient on the path to constraining the global properties of Andromeda’s system of satellite dwarf galaxies and, more broadly, to providing constraints on dwarf galaxy formation and evolution in a cosmological context.

Unified Astronomy Thesaurus concepts: Local Group (929); Dwarf galaxies (416); Andromeda Galaxy (39)

Supporting material: machine-readable table

## 1. Introduction

In cosmological models that include dark matter, hundreds of subhalos of dark matter are orbiting around a central halo (Klypin et al. 1999; Moore et al. 1999). These subhalos are the birthplace of dwarf galaxies, but not all will “light up.” The expected number of dwarf galaxies is highly sensitive to the cosmology (Spergel & Steinhardt 2000; Bode et al. 2001) but also to diverse physical phenomena that can stop star formation, such as stellar feedback or reionization (Bullock et al. 2000; Somerville 2002; Mashchenko et al. 2008; Wheeler et al. 2015). Dark matter simulations adding baryon physics are compared to dwarf galaxy observations in order to constrain those parameters (Koposov 2009; Kim et al. 2018; Nadler et al. 2019).

The search for dwarf galaxies was revolutionized by large, homogeneous photometric surveys such as the Sloan Digital Sky Survey (SDSS; Abazajian et al. 2003), the Panoramic System Telescope and Rapid Response System 1 3π survey (Pan-STARRS1; Chambers et al. 2016), and the Dark Energy Survey (DES; Abbott et al. 2018). Those new data sets are searched for dwarf galaxies thanks to detection algorithms that look for overdensities of stars compatible with an old stellar population (Belokurov et al. 2007; Martin et al. 2013; Laevens et al. 2015). Thanks to these improvements, numerous new

dwarf galaxies were detected. Indeed, before those surveys, there were a dozen known dwarf galaxies around the Milky Way (MW) with an absolute magnitude $M_V < -8.8 \pm 0.2$ (Draco; Mateo 1998). By now, the number of known MW satellites is 59, with $M_V < -0.8 \pm 0.9$ (Virgo I; Homma et al. 2016). Even if the nature of those objects is sometimes debated (Conn et al. 2018a; Jerjen et al. 2018; Mutlu-Pakdil et al. 2018), the increase in the number of known dwarf galaxies allows us to determine the dwarf galaxy detection limits of each survey and, hence, to obtain more precise statistical estimations of the real number of dwarf galaxies around the MW. This was already done for SDSS, Pan-STARRS1, and DES searches (Koposov et al. 2008; Drlica-Wagner et al. 2020) and used to derived an estimate of the luminosity function for the MW.

The Local Group hosts another large galaxy that provides us with a different test sample for cosmological and galaxy formation models. Large photometric surveys also changed our knowledge on M31 and its satellites. Searches focused on the SDSS led to the discovery of three new dwarf galaxies around M31, as faint as $M_V = -8.1 \pm 0.5$ with And X (Zucker et al. 2004, 2007; Bell et al. 2011). The INT/WFC imaging of the surroundings of M31 led to the discovery of And XVII (Irwin et al. 2008). But the main survey leading to the discovery of numerous new M31 dwarf galaxies is the Pan-Andromeda Archaeological Survey (PAndAS), a dedicated survey of the halo of the Andromeda galaxy. Searches within this survey led to the discovery of 19 satellite dwarf galaxies of M31 with $M_V < -5.9 \pm 0.7$ (Martin et al. 2006; Ibata et al. 2007; McConnachie et al. 2008; Richardson et al. 2011). For this

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

Original content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.
survey, an algorithm was created by the PAndAS team to
search for stellar populations of dwarf galaxies in the full
survey by looking for overdensities of stars both spatially and
in the color-magnitude diagram (CMD; Martin et al. 2013,
hereafter M13). Thanks to the algorithm, we can now
determine the survey completeness limit.

In this paper, we determine the detection limits of Andromeda’s dwarf galaxies using the PAndAS survey. In parallel to similar studies already conducted on the MW, this will provide a second independent test for cosmological or galaxy formation simulations, using a different satellite system. Recovery fractions of dwarf galaxies are obtained by adding artificial galaxies to the PAndAS photometric catalog and by determining whether they are detected by blindly running the search algorithm developed by the PAndAS team (M13).

This paper is structured as follows: Section 2 briefly describes PAndAS and its data. Section 3 details the generation of artificial dwarf galaxies and how we determine their recovery fractions. Section 4 presents the modeled dwarf galaxy completeness for each field of the PAndAS survey. The impact of the galaxy’s distance and metallicity is also investigated. Finally, we discuss our results in Section 5.

## 2. Preliminaries

PAndAS was a Large Program conducted from 2008 to 2011 at the Canada-France-Hawaii Telescope (CFHT), using the wide-field imager MegaCam to map the surroundings of the Andromeda galaxy (McConnachie et al. 2009, 2018). Combined with previous observations obtained through PI time, it includes over 400 fields of ~1 deg² each and probes a region within ~150 kpc of M31 and ~50 kpc of its companion galaxy M33. Each field was observed in both the *g* and *i* bands for photometry that includes the brightest three magnitudes of the red giant branch (RGB) at the distance of M31. The survey is presented in detail in McConnachie et al. (2018), along with the data reduction and catalog creation steps, and we refer the reader to this publication for the full description of the catalogs that we use here. In a nutshell, these reach median, 5σ *g* and *i* depths of 26.0 and 24.8, respectively (Ibata et al. 2014).

The algorithm to generate dwarf galaxies, described and used below, requires two field-specific ingredients: a model of the photometric uncertainties and a model of the completeness. Thomas et al. (2021) present models for these two components that are determined carefully for a small set of fields and propagated to the full survey by scaling these initial models with an anchor magnitude that corresponds to the magnitude at which the uncertainties are equal to 0.1 mag.

For the model of the photometric uncertainties, we model them with an exponential model as a function of magnitude (see Ibata et al. 2007 for more details). The photometric uncertainties in the *i* and *g* bands are shown in Figure 1 for three different fields of the survey with the corresponding models. The exponential fits describe the photometric uncertainties reasonably well, despite small differences that arise between fields.

For the completeness, we use the work of Thomas et al. (2021), who determined, independently, the *i*- and *g*-band completeness of the PAndAS data by comparison with regions of deep photometry observed with the Hubble Space Telescope (HST) data. This process is based on 14 HST fields that are ~2.5–3 mag deeper than PAndAS and propagated to the full

survey via changes in the anchoring magnitude, as described
above.

## 3. Methods

### 3.1. Generating Artificial Dwarf Galaxies

To determine the dwarf galaxy detection limits of the PAndAS survey, we need to generate artificial dwarf galaxies that can be added to the PAndAS photometric catalog to measure their probability of being detected. For any artificial dwarf galaxy, we aim to generate a catalog of fake positions and *g* and *i* magnitudes that, combined, create a system of known properties (total magnitude, metallicity, radius, distance, and/or age).

The first step to generate an artificial dwarf galaxy is to determine the photometry of its individual stars. For simplicity, we assume that the artificial dwarf galaxy can be parameterized by a single age and metallicity. While this is clearly not accurate, the different stellar populations of the faint dwarf galaxies that will be simulated here are not clearly separated in the PAndAS photometry (Martin et al. 2016). For every artificial dwarf galaxy, we go through the following steps:

1. For a given choice of metallicity, age, distance, and total magnitude in the *V* band ($M_V$), we randomly draw stars from the corresponding luminosity function taken from the PARSEC library (Marigo et al. 2017), using the pre-2014 MegaCam photometric system and a Kroupa initial mass function (Kroupa 2001). The corresponding isochrone from the same library provides the color that the randomly drawn stars must follow. Artificial stars are drawn with an absolute magnitude $-5 < M_i < 20$ until the total flux of the system reaches the target total magnitude $M_V$.<fnref n="9" /> We then keep only bright stars, as the selection box of the search algorithm will only use upper RGB stars with $i_0 < 23.5$. The exact shape of the horizontal branch is not modeled in detail, as its stars are fainter than the magnitude limit used by the algorithm, but their flux is properly taken into account to determine the total luminosity of the artificial dwarf galaxy.

2. To best represent the PAndAS observations, we randomize the perfect magnitudes drawn from the isochrone and luminosity function by adding noise to the magnitudes, following the models of photometric uncertainties of the field in which the artificial galaxy will be placed. The added observed photometric uncertainties broaden the locus of stars in the CMD.

3. Finally, we use the (in)completeness model determined by Thomas et al. (2021) to test for the observability of any star drawn in the artificial galaxy. Every star generated in the previous step is tested against the completeness model. For every star drawn from the isochrone and luminosity function, we draw random deviates between 0 and 1 and independently test those against the completeness at the star's *g*- and *i*-band magnitudes.

The CMD of an artificial dwarf galaxy with $M_V = -8.5$, $[Fe/H] = -1.7$, and an age of 10 Gyr (and $r_h = 265$ pc) is shown in Figure 2 as it is pushed through the different steps of

<footnote n="9">⁹ MegaCam magnitudes are transformed into V-band magnitudes using the color equations presented in Ibata et al. (2014).</footnote>
<figure><img src="image_4.png" /><figcaption>**Figure 1.** Uncertainties in the *i* (left) and *g* (right) magnitudes for three typical fields of the survey. The exponential fits (green lines) reasonably follow the photometric uncertainty distributions of each field.</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>**Figure 2.** CMD of the stars of an artificial dwarf galaxy with $M_V = -8.5$, [Fe/H] = -1.7, an age of 10 Gyr, and $r_h = 265$ pc, as it is pushed through the different steps of the generation process. The first panel corresponds to the CMD of stars drawn from the chosen luminosity function and corresponding isochrone. As the search algorithm is only looking at RGB stars, we only represent the bright stars that are of interest here. The second panel shows the CMD after taking the PAndAS photometric uncertainties into account. The CMD of the third panel folds in the incompleteness of the PAndAS data. Finally, the fourth panel shows the CMD of this artificial system within $2r_h$ overlaid on top of the typical field contamination from PAndAS (gray points). A direct comparison with the CMD within $2r_h$ of And XIV (fifth panel) that has similar parameters shows the similarities between the two CMDs. The sixth panel displays the perfect agreement between the radial surface density of stars in the artificial dwarf galaxy (blue points) and the chosen radial density model (green line).</figcaption></figure>

the process. We can see that the resulting CMD in the fourth panel, which is the one of an artificial dwarf galaxy added on a random location of the PAndAS survey without any known satellite or stellar stream, is very similar to the CMD of the real dwarf galaxy And XIV, which has a similar total luminosity ($M_V = -8.5$; Martin et al. 2016) and is shown in the fifth panel.

Now that the number of stars in the galaxy and their photometry are known, the next step is to determine their position with respect to the center of the galaxy. The radial stellar density is assumed to follow an exponential law with scale radius $r_e = r_h/1.68$, where $r_h$ is the half-light radius of the system in pc. The probability density function of stars, $P(r)$, at
radius $r$ given in pc, is

$$P(r) = \frac{1.68^2 r}{2\pi r_h^2} \exp\left(-1.68 \frac{r}{r_h}\right). \quad (1)$$

For every star in the artificial dwarf galaxy that was assigned a color and a magnitude in the previous step and passed the completeness test, we draw its radius from this probability distribution function (pdf) and assign it a random angle between 0 and $2\pi$. These polar coordinates are then transformed into sky coordinates using the center of the artificial dwarf galaxy. The right panel of Figure 2 shows the perfect agreement between the radial density profile of stars in this test artificial dwarf galaxy and the chosen model density profile.

## 3.2. Determining the Recovery Rate of Artificial Galaxies

There are clear differences in the properties of observed stars throughout the PAndAS survey: for instance, the contamination from MW foreground stars increases significantly toward the north (M13), and the density of M31 stellar halo stars changes with the density of stellar stream stars (Ibata et al. 2014). We therefore aim to determine the recovery fractions of artificial dwarf galaxies (and therefore the detection limits of dwarf galaxies in the survey) as a function of the location in the survey. A degree, at the distance of M31, corresponds to ~14 kpc. This means that the $1^\circ \times 1^\circ$ MegaCam field is a natural areal subdivision to consider since it corresponds to a sizable area of the M31 halo without being too large. The fact that there are subtle field-to-field differences in the photometry depth and completeness, despite the survey being very homogeneous, further comforts us in our choice to determine the artificial dwarf galaxy recovery fractions on a field-to-field basis.

### 3.2.1. Choice of Parameters for the Artificial Dwarf Galaxies

We first need to carefully decide which parameters have the most significant impact on the detection limits so as to limit as much as possible the significant amount of computation required to determine the dwarf galaxy recovery fractions.

We know that the parameters with the most impact on the detection of a dwarf galaxy (or absence thereof) are the size and the total luminosity of a system because these are directly related to the surface brightness of a dwarf galaxy (e.g., Koposov 2009; Walsh et al. 2009; Drlica-Wagner et al. 2020). The distance to the system is very important for searches of MW dwarf galaxies but is not so important for distant systems. At the distance of M31, a change in distance of $\pm 300$ kpc only leads to a change in the location of stars in the CMD by $-1.1/+0.5$ mag, respectively. While this is not subtle, the fact that PAndAS only observes RGBs that are sparsely populated, and so very noisy, means that a shift in distance is similar to a shift in metallicity, and so the distance does not have a very significant impact on the recovery rate of a dwarf galaxy (M13).

To limit the required calculations, we therefore consider the distance to an artificial system as a secondary parameter and, for the moment, assume that all artificial dwarf galaxies are located at the distance of M31, with a distance modulus of 24.47 (McConnachie et al. 2005).

The range of interest for a dwarf galaxy's total magnitude is $-4.5 \lesssim M_V \lesssim -8.5$ (Martin et al. 2016) and straddles the total magnitude of the faintest dwarf galaxy detections. In this regime, all known dwarf galaxies are invariably metal-poor with a metallicity contained between $-1.5$ and $-2.3$ (Tollerud et al. 2012; Collins et al. 2013). In addition, the PAndAS photometry is not very sensitive to metallicity variations in this range, as the tracks followed by RGB stars in this metallicity range almost overlap in the $(g-i, i)$ CMD. To save on unnecessary computing time, we therefore assume a fixed metallicity for all dwarf galaxies, with $[Fe/H] = -1.7$.

Finally, we know that the ages of the stellar populations in a dwarf galaxy have little impact on the color and magnitude of its RGB stars as long as they are at least moderately old ($\gtrsim 2$ Gyr). Consequently, we assume an age of 10 Gyr for all artificial dwarf galaxies (Weisz et al. 2019).

We check the impact of these three assumptions on the distance, metallicity, and age in Section 4 and confirm that it is indeed minimal, or can easily be modeled in the case of the distance.

In summary, for each field in the PAndAS footprint, we generate artificial galaxies with a metallicity $[Fe/H] = -1.7$, an age of 10 Gyr, at a distance modulus of M31 ($m-M=24.47$). To test the impact of the total magnitude ($M_V$) and the size (parameterized by the half-light radius $r_h$ given in pc) of the system, we bin the $M_V-\log_{10}(r_h)$ plane over the generous range $-8.5 < M_V < -4.5$ and $1.8 < \log_{10}(r_h) < 3$, with bin sizes of 0.25 and 0.1, respectively. From known M31 dwarf galaxies, we know that the transition from detected to undetected systems is within this range (Martin et al. 2016). For each of these bins, we generate an artificial galaxy with a random size and magnitude within the limits of this bin and ingest the photometry and position of its stars in the PAndAS catalog at a random location within the considered field. Stars that fall in neighboring fields are kept, but we remove any star of the artificial galaxy that falls in a CCD gap, in a hole between two fields, or outside of the PAndAS footprint, and we then test for the detection significance of this system (see below). This step is usually repeated 5 times so we can determine a recovery fraction for the studied magnitude-size bin and field. Ideally, we would like to ingest more than five artificial galaxies, but, given the $\sim 400$ fields and the number of bins in the magnitude-size plane, this already corresponds to about half a million ingested artificial galaxies and running the search algorithm 1000 times since ingesting all the artificial dwarf galaxies at once would lead to overlapping systems, change the properties of the survey, and bias the recovery rates. However, we will confirm below from a small subset of representative fields that ingesting only five artificial dwarf galaxies per $M_V-\log_{10}(r_h)$ bin and per field still yields reliable results.

Finally, we remove regions near known dwarf galaxies where the search algorithm would invariably return a detection, irrespective of the ingested system, and the regions close to M31 and M33, for which the search algorithm does not work well because of the complex mix of stellar populations (M13). We therefore mask regions within $4r_h$ of all known dwarf galaxies, as well as within $\sim 2^\circ$ and $\sim 1^\circ$ of M31 and M33, respectively. In the south of the survey, a group of background nearby elliptical galaxies with globular clusters that masquerade as M31 RGB stars also needs to be masked out. We choose not to mask M31's globular clusters, as they are rarely
<figure><img src="image_6.png" /><figcaption>**Figure 3.** Local density of PAndAS stars in the CMD selection box (left panel). The white ellipses correspond to the regions we mask because of known stellar systems, and the three fields chosen to be representative of the various regions of the survey are highlighted in red. The additional four fields used to test the effect of the distance on the recovery fractions are represented in purple. The right panel shows the distribution of PAndAS fields, following the numbering scheme of Ibata et al. (2014) and McConnachie et al. (2018). Masked regions are here shown in cyan.</figcaption></figure>

recovered by the search algorithm and only one of these is detected in M13 with a significance above 6. The resulting mask is shown in Figure 3. The detection limits of fields that are entirely within the masked regions are not determined, and for partially masked fields we make sure to insert artificial dwarf galaxies in the nonmasked part of the field.

### 3.2.2. Setup of the Search Algorithm

To test for the recovery of an artificial galaxy, we use the algorithm developed by M13 to search for dwarf galaxies in the PAndAS survey. The algorithm is computationally costly to run, as it determines the likelihood of there being a dwarf galaxy at all locations of PAndAS on a 0.5 grid (115 pc at the distance of M31), given the distribution of local RGB stars on the sky and in the color–magnitude space. The assumption is that this distribution can be well reproduced by three different components: one that simulates a metal-poor, compact dwarf galaxy, one that simulates the locally constant contamination of M31’s halo stars, and a third part that represents the contamination from MW stars, following an empirically built model. We refer the reader to M13 for the full description of the search algorithm and its parameterization, but we mention here our choices for the model parameter grid that we let the algorithm explore. Ideally, it would be best to explore a broad range of choices for all parameters of the model, but, as mentioned in M13, this can be quite costly, and we restrict the exploration to likely model parameters.

The five model parameters we focus on are as follows:

1. the metallicity of the dwarf galaxy, $[Fe/H]_{dw}$;

2. the half-light radius of the dwarf galaxy, $r_h$;

3. the number of stars in the dwarf galaxy component of the model, $N$;

4. the metallicity of the M31 halo contamination, $[Fe/H]_{halo}$; and

5. the fractional contribution of MW stars to the total local contamination, $\eta$.

The five-dimensional grid of parameters for which the algorithm determines a likelihood value is listed in Table 1. Similarly to M13, the favored model is the one that maximizes

<table><caption>**Table 1**
Ranges of Parameters Used to Run the Search Algorithm</caption><thead><tr><td>Parameter</td><td>Minimal Value</td><td>Maximal Value</td><td>Step</td></tr></thead><tbody><tr><td>$[Fe/H]_{halo}$</td><td>-1.3</td><td>-0.6</td><td>0.1</td></tr><tr><td>$[Fe/H]_{dw}$</td><td>-2.3</td><td>-1.1</td><td>0.3</td></tr><tr><td>$r_h$</td><td>0.5</td><td>3.5</td><td>1'</td></tr><tr><td>$\log_{10}(N)$</td><td>-0.5</td><td>4</td><td>0.5</td></tr><tr><td>$\eta$</td><td>0</td><td>1</td><td>0.1</td></tr></tbody></table>

the likelihood, and we use the comparison of this maximal likelihood with the likelihood of the best model with $N \simeq 0$ to determine the significance $S$ of a detection.

Since we know the location of the artificial dwarf galaxies ingested in the survey catalog, we can significantly save on computing time by only running the search algorithm near these locations. Therefore, we only test for the significance of a dwarf galaxy detection within $\pm 1.5$ of the known centers of ingested dwarf galaxies in both the $\xi$ and $\eta$ directions (or $\pm 3$ steps on the spatial grid over which the algorithm is run).$^{10}$ The corresponding locations are highlighted in Figure 4 for a run with artificial galaxies that have $r_h = 1$ kpc.

### 3.2.3. Choice of the Detection Threshold

The search algorithm outputs the significance $S$ of there being a dwarf galaxy at all tested locations, but we still need to determine what constitutes a detection. If the threshold we consider for a detection is too low, it will lead to the detection of systems that would be too faint to be reliably classified as dwarf galaxies in PAndAS, and it will not be possible to use the resulting detection limits to understand the M31 satellite system as seen by PAndAS. The threshold is derived from the distribution of significances determined by M13 for an area with low M31 and MW contamination near And XI–XIII (see their Figure 9). To determine expectations of detections stemming from noise in the data, we linearly fit the low-

<footnote>$^{10}$ We tested on field 81 that the recovery fractions are similar when the significance calculations are made within a distance of $\pm 2.5$ and $\pm 1.5$ of the artificial galaxy center.</footnote>
<figure><img src="image_7.png" /><figcaption>**Figure 4.** Example of the distribution of artificial dwarf galaxy stars in the survey (red points). The large polygons correspond to the limits of the PAndAS fields in this region. The dark squares each correspond to the 7 × 7 test locations for the search algorithm and are chosen to be around the known location of an artificial galaxy's center.</figcaption></figure>

significance tail of the significance distribution and choose as the significance threshold $S_{\text{th}}$ the value of $S$ that corresponds to less than one detection in this area of ~9 deg$^2$. This yields $S_{\text{th}} = 6$. It is coherent with the significance from the dwarf galaxy with the lowest value of $S$ that is unambiguously detected, and XXVI, for which the algorithm determines $S_{\text{th}} = 5.9$ (M13).

The dwarf galaxy recovery fractions for three representative fields are shown in Figure 5 for this threshold. These three fields are representative of an outer M31 halo region (field 81), a region contaminated by dense M31 stellar streams (field 227), and a region heavily impacted by the MW foreground contamination (field 375), as highlighted in Figure 3. We will discuss the recovery fractions in more detail in the next subsection, but at this stage we wish to focus on the fact that, overall, they behave as expected, with faint/compact dwarf galaxies recovered at a high fraction and larger and/or fainter systems showing lower detection rates. However, we note that for field 227 the recovery fractions are nonzero for systems that the search algorithm has no hope of detecting in PAndAS, as they are too faint and too extended (region highlighted in red in the figure). This behavior is typical of fields that are heavily contaminated by M31 stellar structures. In these fields, the M31 streams have a color-magnitude distribution that is similar to what the algorithm searches for in a dwarf galaxy, which biases upward the values of $S$. This was already shown by M13, and it forces us to increase the detection threshold for those fields to avoid false detections and unrealistically high recovery fractions. In practice, we increase $S_{\text{th}}$ by steps of 0.5 until there are only at most two galaxies detected in the region of faint and extended dwarf galaxies highlighted in red in the upper left corner of the recovery fraction panels of Figure 5. Our aim is to ensure a detection threshold that leads to only

a very small number of false positives, but to not have a threshold that is so conservative that only bright dwarf galaxies are recovered. The resulting recovery rate for field 227 is shown in the second column of panels of the figure, along with the map of final values of $S_{\text{th}}$ throughout the survey in the right panel. As expected, the nominal detection threshold value, $S_{\text{th}} = 6$, is used for most fields, and only fields close to M31 and significantly affected by M31 stellar halo structures require higher detection threshold values.

### 3.3. The Recovery Fraction Model

The large amount of computing time required to determine the recovery fractions of artificial dwarf galaxies limits our calculations to only five simulated dwarf galaxies per field and per magnitude-size pixel, which leads to somewhat noisy recovery fractions. In this subsection, we aim to build an analytical recovery fraction model for each field to bypass these limitations.

The leftmost column of panels in Figure 6 shows the outcome of the algorithm for the same three fields 81, 227, and 375 already discussed above. In this case, we run the algorithm with 20 galaxies per ($M_V$, $\log_{10}(r_h$)) pixel and per field, to validate that the model we build on the regular simulations with only five dwarf galaxies per pixel is representative. In these panels, we also show the properties of known M31 dwarf galaxies located in the PAndAS survey (Martin et al. 2016). As expected, these known galaxies are located in the area of high recovery fractions but also cover the area of transition, in agreement with the results of the search algorithm (M13) and our choice for the detection threshold values. For comparison, the second column of panels displays the recovery fractions of dwarf galaxies with our regular setup of five dwarf galaxies per pixel. It is reassuring that, even though the recovery fractions are here more noisy, as expected with fewer artificial galaxies, they do look similar to the more accurate results with 20 galaxies per pixel. We will provide a quantitative comparison below.

The recovery fractions behave as expected: high surface brightness systems, i.e., compact and/or bright dwarf galaxies in the lower right portion of the plots, are recovered at high efficiency, while more extended and/or fainter systems are progressively missed by the search algorithm. The field contamination (mainly from the M31 stellar halo in field 227 or the MW in field 375) has an impact on the recovery fractions and, in general, degrades them. We note that, in the case of a high level of contamination from M31 stellar halo structures, the algorithm recovers a small number of large and faint dwarf galaxies below the transition region. As mentioned above, these are false detections and an artifact of the search algorithm that is sometimes struggling to discriminate contaminating M31 stream stars from the dwarf galaxy stars at the same location of the CMD because they all follow similar RGB tracks. However, this effect remains small since, as mentioned above, we adapt the threshold limit to minimize the presence of these false detections. This effect is not present for fields heavily contaminated by foreground MW stars (field 375 in Figure 6) since these contaminating stars have a color-magnitude distribution that is quite different from that of the M31 (artificial) dwarf galaxies (M13). In this case, the added contamination leads to recovery fractions that are shifted toward brighter and/or more compact systems.
<figure><img src="image_8.png" /><figcaption>**Figure 5.** Left: recovery fractions for three fields of the PAndAS survey, both with the common detection threshold value, $S_{th} = 6$, and with the updated threshold value determined to minimize the number of false positives in the region of large and faint dwarf galaxies that is highlighted in red. Right: detection thresholds that are used for all considered PAndAS fields. Only fields with a high level of M31 stellar halo contamination have a higher $S_{th}$, as the search algorithm struggles to always distinguish the structures of M31 and dwarf galaxies.</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>**Figure 6.** Recovery fractions for fields 81, 227, and 375. The first column of panels presents the recovery fractions when simulating 20 galaxies per field and per ($M_V$, log$_{10}$(r$_h$)) pixel. The red circles correspond to known dwarf galaxies around M31. These straddle the transition region between no recovery and full recovery, as expected. The second column of panels presents the recovery fractions for only five galaxies per field and per ($M_V$, log$_{10}$(r$_h$)) pixel. The third column shows the modeled recovery fractions. In these panels, the red line represents the 50% detection threshold. The fourth column displays the residuals that remain small, confirming that our model reproduces the recovery fractions well.</figcaption></figure>

We fit an analytical model to the recovery fractions to overcome the noise and the false detections. Given the shape of the recovery fraction distributions, our analytical model is built around a quadratic transition region in the $M_V$-log$_{10}$(r$_h$) plane.

We allow for the parameters of this quadratic transition to vary from field to field, and we parameterize the width of the smooth transition between regions of full to no recovery.
<figure><img src="image_10.png" /><figcaption>**Figure 7.** The top panels show the pdf's for each parameter of the model for the case with 20 and five artificial dwarf galaxies per field and per ($M_v$, log$_{10}$(r$_h$)) pixel for field 81. The bottom panels show the pdf's for the case of five galaxies and for different fields. While using only five galaxies per pixel produces wider pdf's, as expected, these show no systematic bias compared to the pdf's resulting from the 20-galaxy case.</figcaption></figure>

We define the model for the recovery fraction of artificial dwarf galaxies as

$$f(M_v, \log_{10}(r_h)) = F\left(\frac{M_v - M_{v,\lim}}{\sigma}\right), \\ \text{with } M_{v,\lim} = \mathbf{a} \log_{10}(r_h)^2 + \mathbf{b}. \qquad (2)$$

Here a and b do not carry any significant physical meaning and are just used to parameterize the space, σ is the width of the transition, and F is the complementary error function, defined as

$$F(x) = \frac{1}{2} \text{erfc}\left(\frac{x}{\sqrt{2}}\right). \qquad (3)$$

For each field, we fit the model to the recovery fraction distribution in the $M_v-\log_{10}(r_h)$ plane by determining the likelihood of the data given the model thanks to a Metropolis-Hastings method. The variables a, b, and σ are only used to parameterize the model; the important result is the threshold in surface brightness (which we will use for further analyses). Then, the best set of parameters is obtained by calculating the binomial likelihood of the number of recovered dwarf galaxies compared to the model. For the three representative fields and the case with five artificial galaxies per pixel, the resulting best-fit models are displayed in the third column of panels in Figure 6, with the residuals shown in the fourth column. The models provide a good representation of the recovery fractions, and the residuals remain small.

To confirm that we do not introduce any bias by simulating only five galaxies per ($M_v$, log$_{10}$(r$_h$)) pixel, in Figure 7 we compare the marginalized pdf's obtained for the three parameters in field 81 when using 20 (dashed lines) and subsamples of five (solid lines) simulated galaxies per pixel. Even though they are, as expected, wider, we find no systematic bias in the pdf's obtained when using only five dwarf galaxies per pixel. Therefore, we go ahead and determine the recovery fractions with five artificial galaxies per pixel and per field, as it already represents more than 90,000 CPU hours on the Strasbourg University High Performance Computing center for the simulations over the ~400 PAndAS fields.

# 4. Results

## 4.1. Dwarf Galaxy Recovery Fractions over the PAndAS Survey

Figure 8 summarizes the recovery fractions of all studied PAndAS fields after ingesting more than 350,000 artificial dwarf galaxies in the stellar catalog. The top three panels of the figure show, for each field, the values of the three model parameters over the survey. There are some clear changes to the values of parameters a and b that can be tracked to changes in the properties of the survey: the north-south gradient is linked to the increased MW foreground contamination toward the north, and the presence of the Giant Stream south of M31 (0 ≲ ξ ≲ 5° and -5 ≲ η ≲ 0°) clearly impacts the values of a. On the other hand, the speed of the transition, parameterized by σ, remains fairly constant over the survey.

The slope of the transition region in the $\log_{10}(r_h)-M_v$ plane does not map lines of constant surface brightness. While this is common (e.g., Koposov et al. 2008; Walsh et al. 2009; Drlica-Wagner et al. 2020) and is likely related to the complex nature of the data and of the search algorithm, it can make it difficult to interpret our results in terms of surface brightness limits. We therefore show, in the bottom panel of Figure 8, the corresponding surface brightness, $\mu_{250}$, within the half-light radius of a system for a fixed $r_h = 250$ pc. Most fields show $29 \text{ mag arcsec}^{-2} < \mu_{250} < 30 \text{ mag arcsec}^{-2}$, with a clear impact from the MW contamination that is significantly more important on the northern side of the survey that reaches Galactic latitude $b - 12°$. One can also note the impact of the M31 stellar halo contamination in the regions nearest M31. The necessity to increase the detection threshold $S_{th}$ to avoid false positives is partly responsible (Figure 6), but the lower surface brightness limit is also linked to a lower contrast because of a higher density of contaminating M31 stars. The model parameters for all fields are listed in Table 2.

## 4.2. Impact of Fixed Artificial Dwarf Galaxy Parameters

As mentioned in Section 3, we fixed some of the supposedly less impactful parameters of the artificial dwarf galaxies to save
<figure><img src="image_11.png" /><figcaption>**Figure 8.** Values of the favored model parameters, *a*, *b*, and *σ*, for each field in the PAndAS survey (top three panels). The bottom panel translates those into the surface brightness within the half-light radius 50% detection limit for a dwarf galaxy with $r_h = 250$ pc (see text for more details). The fields highlighted in red are those chosen to be representative of the different regions of the survey, while those highlighted in purple are the additional fields used to test the impact of the distance on the recovery fractions. The surface brightness threshold varies depending on the position: in regions with low contamination fainter galaxies can be detected compared to regions with significant M31 or MW contamination.</figcaption></figure>

on the computing time, which is already quite large. Here we test the impact of these parameters on the recovery fraction model parameters for seven representative fields of the survey. Those are chosen to be dispersed throughout the survey to represent different types of contamination.

To determine the impact of the distance to the artificial dwarf galaxies, we conduct a set of simulations with five galaxies per pixel and distances that vary between ±300 kpc of the distance to M31 that is assumed for the main run. The resulting variations in the recovery model parameters are shown in Figure 9 for steps of 150 kpc. One can expect that for a closer/further distance the recovery fractions shown in Figure 6 would simply shift to the upper left/lower right. This is indeed the main effect that we seen in the left panel of Figure 9, with *b* typically shifting by the distance modulus change that corresponds to the distance change. We note, however, that changes in the distance to the artificial dwarf galaxies also have some impact on *a* of the recovery fraction model (middle

panel). While the reason for this correlation is not readily evident, it may be related to the algorithm needing a minimum number of stars in the M31 CMD box to yield a reliable detection, with a smaller impact from the size of the system in this regime. Such an effect would tend to make the slope of the transition region of the recovery fractions steeper as the number of stars in the CMD box becomes smaller. It could explain why more distant artificial galaxies, which have fewer stars in the magnitude-limited CMD box, yield higher values of *a*.

Irrespective of the reason for the changes of *a* and *b* with distance, these anticorrelations are regular and can easily be modeled. This has the benefit of preventing numerous additional simulations while still taking into account that dwarf galaxies significantly in front of or behind M31 have somewhat different recovery fractions. We therefore fit a linear model to the values of *a* and *b* shown in Figure 9, using the data from all fields together, but allowing for different intercepts for each field. We define as $a_0$ and $b_0$ the intercepts (their values are
<figure><img src="image_12.png" /><figcaption>Figure 9. Changes to the best parameters of the recovery fraction model as a function of the distance modulus to the dwarf galaxy. Each point is slightly shifted in distance, in order to make the graph more legible. The black dashed line represents the theoretical change of $M_V$ as a function of the distance modulus for the field 50. Parameters a and b smoothly change with distance, which allows us to easily model the influence of the distance on the recovery fractions.</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>Figure 10. Changes to the best parameters of the recovery fraction model as a function of a dwarf galaxy’s metallicity. Variations to the metallicity within this metal-poor range do not significantly change the parameters of the model.</figcaption></figure>

accessible in the Appendix) and $\mu$ the distance modulus. The resulting fits yield

$$ 
\begin{aligned}
a &= 0.06 \mu + a_0, \\
b &= -0.89 \mu + b_0.
\end{aligned}
\qquad (4) $$

The third parameter of the model, $\sigma$, which represents the width of the transition region of the recovery fractions in the $\log_{10}(r_h)-M_V$ space, shows no significant changes with distance. We therefore assume that it is not impacted by changes in the distance to the dwarf galaxies.

In Section 3, we also fixed the metallicity of the metal-poor artificial dwarf galaxies to [Fe/H] = -1.7, as the color of the RGB only shifts slightly in color in this metal-poor regime. The impact of this choice of metallicity is explored by determining the recovery fraction for artificial galaxies in the range -2.3 < [Fe/H] < -1.5, with steps of 0.1 dex, for fields 81, 227, and 375. The results are statistically similar over the metallicity range, as shown in Figure 10. We conclude that the impact of [Fe/H] on the recovery fraction is negligible and that our assumption of [Fe/H] = -1.7 does not affect our results.

### 4.3. Using the Recovery Fractions

To facilitate the use of the previous results and, for instance, easily apply the PAndAS recovery fraction to sets of simulated dwarf galaxy satellite systems, we make public a Python module that is accessible at https://github.com/dolivadolinsky/Recovery_dwarf_galaxy_M31. Given the position of a dwarf galaxy around its host, its half-light radius, and its V-band magnitude, this module returns the probability that such a dwarf galaxy would have been discovered in PAndAS. The module first determines the fields in which the tested dwarf galaxy would be located or if it is outside the PAndAS footprint. If there is an overlap of fields, the parameters values for the galaxy are the ones of the deepest photometry. Then, the relations described in Section 4.2 are used to determine the values of a and b at the galaxy’s distance in order to retrieve the recovery fractions.

## 5. Discussion and Conclusion

In this paper we determined the recovery fraction of Andromeda’s dwarf galaxy based on the algorithm of M13 applied to the photometric data from the PAndAS survey. We repeatedly simulated artificial dwarf galaxies, blindly added them to the PAndAS photometric catalog, and checked the significance of their recovery by running the search algorithm. We determined and modeled the recovery fractions of satellites on a field-to-field basis for the 390 $1^\circ \times 1^\circ$ fields of the survey and as a function of a dwarf galaxy’s magnitude, half-light radius, and heliocentric distance. We publish the resulting recovery fraction as a Python module that provides the probability of observing a dwarf galaxy of a given magnitude, half-light radius, and distance in the PAndAS survey.

As expected, we found that the recovery fractions are highly sensitive to the size and luminosity of a dwarf galaxy but also to its position owing to the significant changes in the level of MW and M31 stellar halo contamination throughout the PAndAS footprint. Constant recovery fractions do not map constant surface brightness levels, likely because of the complex nature of the data and of the search algorithm. For galaxies with $r_h = 250$ pc, the surface brightness threshold that corresponds to a 50% detection rate mainly varies between ~29.0 and 30.0 mag arcsec$^{-2}$, with a median value of 29.8. These results are in good agreement with the surface brightness of genuine dwarf galaxies discovered in PAndAS (Martin et al. 2016). Depending on the field of the survey, these limits can be as low as 28.0 mag arcsec$^{-2}$ or as high as 30.5 mag arcsec$^{-2}$. The completeness of the survey also depends on the distance to a dwarf galaxy, but this variation of the surface brightness threshold roughly maps expectations from changes to the distance modulus.

Our results bear similarities to studies based on other data sets but using the same methodology of inserting stars from artificial dwarf galaxies at the photometric catalog level. Even if the given surface brightness threshold values are not always
comparable, we can calculate the $\mu_{250}$ values based on the results from other studies. Koposov et al. (2008) and Walsh et al. (2009) derive $\mu_{250} \sim 28$ mag arcsec$^{-2}$ at the distance of M31, but for the significantly shallower SDSS data and without considering the specific region around M31 that becomes significantly contaminated by foreground MW stars. Recently, Drlica-Wagner et al. (2020) determined dwarf galaxy recovery fractions for a large part of the sky from both the PAndAS STARRS1 and DES photometric data. The surface brightness thresholds obtained for the two surveys in this study are 28 and 30 mag arcsec$^{-2}$, respectively. It may seem surprising that these results based on the DES are similar to those we obtain for the deeper PAndAS data. First, we note that the PAndAS footprint (latitude $|b| \sim 22^{\circ}$) is more contaminated by MW and M31 stars than a typical DES field (latitude $|b| \ge 30^{\circ}$). In addition, the surface brightness limits are sensitive to the distance and the depth of surveys, as well as the algorithm used to search for dwarf galaxies; those can be an explanation for the similarity between the DES thresholds and the ones obtained in this work. Finally, the detection threshold ($S_{th}$ in our work) used by Drlica-Wagner et al. (2020) to determine whether an artificial dwarf galaxy is recovered is particularly low ($S_{th} = 4.9$). In fact, it is lower than the value $S_{th} = 6.0$ used by Drlica-Wagner et al. (2015) to detect candidate satellites like Tucana V or Cetus II that were not confirmed as well-defined, localized stellar overdensities from deeper photometry (Conn et al. 2018a, 2018b). A particularly aggressive detection threshold would naturally lead to deeper detection limits, but these would be overestimated compared to the detection limits of genuine dwarf galaxies, as the low threshold would lead to false positives in the list of candidate detections.

It is also interesting to compare our results with those of Garling et al. (2021) for the more distant galaxies NGC 3077 and NGC 4214, located more than 3 Mpc away. While the data of the LBT-SONG survey used by these authors are comparable to or deeper than the PAndAS data, the larger distance forces them to ingest artificial dwarf galaxy stars directly in the images to properly account for crowding that becomes an issue because of the smaller angular size of the systems. Their analysis yields $\mu_{250} \sim 28$ mag arcsec$^{-2}$, toward the brighter end of the detection limits we determined from M31, and shows one of the limitations of ground-based surveys beyond the Local Group.

Finally, we compare our results with the ones of Huxor et al. (2014), who determined the globular cluster completeness of the PAndAS survey. They found a 50% detection surface brightness inside the half-light radius threshold of $\sim 26$ mag arcsec$^{-2}$. Compared to our results, the brighter limit can be explained by the fact that they looked for more compact objects and that the search for globular clusters was performed through visual inspection of the images.

The work presented here complements our buildup of a thorough understanding of the dwarf galaxy satellite system of M31 (Martin et al. 2013, 2016) by clearly expressing the limits

to our dwarf galaxy searches within the PAndAS survey. Our results can be seen as an intermediate step toward the dwarf galaxy searches that will be possible within $\sim 1$ Mpc with the next generation of ground-based photometric surveys like the Legacy Survey of Space and Time (LSST; Ivezic et al. 2019). More immediately, we aim to use these recovery fractions to infer the underlying properties of the M31 dwarf galaxy satellite system (number of dwarf galaxies, distribution, size–luminosity relation, etc.) to provide constraints on dwarf galaxy formation and evolution in a cosmological context (A. Doliva-Dolinsky et al. 2022, in preparation).

Based on observations obtained with MegaPrime/MegaCam, a joint project of CFHT and CEA/DAPNIA, at the Canada-France-Hawaii Telescope (CFHT), which is operated by the National Research Council (NRC) of Canada, the Institut National des Sciences de l'Université de la Centre Nationale de la Recherche Scientifique (CNRS) of France, and the University of Hawaii. The authors would like to acknowledge the High Performance Computing Center of the University of Strasbourg for supporting this work by providing scientific support and access to computing resources. Part of the computing resources were funded by the Equipex Equip@Meso project (Programme Investissements d'Avenir) and the CPER Alsacacul/Big Data. N.F.M., R.I., and Z.Y. gratefully acknowledge support from the French National Research Agency (ANR) funded project "Pristine" (ANR-18-CE31-0017) and from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation program (grant agreement No. 834148). G.T. acknowledges support from the Agencia Estatal de Investigación (AEI) of the Ministerio de Ciencia e Innovación (MCINN) under grant with reference (FJC2018-037323-I).

## Appendix
Model Parameters

<table><caption>Table 2
Model Parameters for Each Field of the Survey</caption><thead><tr><th>Field</th><th>a</th><th>b</th><th>σ</th><th>a<sub>0</sub></th><th>b<sub>0</sub></th></tr></thead><tbody><tr><td>1</td><td>-0.27</td><td>-4.69</td><td>0.48</td><td>-1.74</td><td>14.89</td></tr><tr><td>2</td><td>-0.31</td><td>-4.40</td><td>0.34</td><td>-1.78</td><td>15.18</td></tr><tr><td>3</td><td>-0.29</td><td>-4.41</td><td>0.41</td><td>-1.76</td><td>15.17</td></tr><tr><td>4</td><td>-0.28</td><td>-4.49</td><td>0.38</td><td>-1.75</td><td>15.09</td></tr><tr><td>5</td><td>-0.32</td><td>-4.42</td><td>0.38</td><td>-1.79</td><td>15.16</td></tr><tr><td>6</td><td>-0.28</td><td>-4.30</td><td>0.48</td><td>-1.75</td><td>15.28</td></tr><tr><td>7</td><td>-0.27</td><td>-4.44</td><td>0.51</td><td>-1.74</td><td>15.14</td></tr><tr><td>8</td><td>-0.31</td><td>-4.44</td><td>0.37</td><td>-1.78</td><td>15.14</td></tr><tr><td>9</td><td>-0.26</td><td>-4.44</td><td>0.46</td><td>-1.73</td><td>15.14</td></tr><tr><td>10</td><td>-0.26</td><td>-4.5</td><td>0.43</td><td>-1.73</td><td>15.08</td></tr></tbody></table>

(This table is available in its entirety in machine-readable form.)
ORCID iDs

Amandine Doliva-Dolinsky https://orcid.org/0000-0001-9775-9029

Nicolas F. Martin https://orcid.org/0000-0002-1349-202X

Guillaume F. Thomas https://orcid.org/0000-0002-2468-5521

Annette M. N. Ferguson https://orcid.org/0000-0001-7934-1278

Rodrigo A. Ibata https://orcid.org/0000-0002-3292-9709

Geraint F. Lewis https://orcid.org/0000-0003-3081-9319

Dougal Mackey https://orcid.org/0000-0002-6529-8093

Alan W. McConnachie https://orcid.org/0000-0003-4666-6564

References

Abazajian, K., Adelman-McCarthy, J. K., Agüeros, M. A., et al. 2003, AJ, 126, 2081

Abbott, T. M. C., Abdalla, F. B., Allam, S., et al. 2018, ApJS, 239, 18

Bell, E. F., Slater, C. T., & Martin, N. F. 2011, ApJL, 742, L15

Belokurov, V., Zucker, D. B., Evans, N. W., et al. 2007, ApJ, 654, 897

Bode, P., Ostriker, J. P., & Turok, N. 2001, ApJ, 556, 93

Bullock, J. S., Kravtsov, A. V., & Weinberg, D. H. 2000, ApJ, 539, 517

Chambers, K. C., Magnier, E. A., Metcalfe, N., et al. 2016, arXiv:1612.05560

Collins, M. L. M., Chapman, S. C., Rich, R. M., et al. 2013, ApJ, 768, 172

Conn, B. C., Jerjen, H., Kim, D., & Schirmer, M. 2018a, ApJ, 852, 68

Conn, B. C., Jerjen, H., Kim, D., & Schirmer, M. 2018b, ApJ, 857, 70

Drlica-Wagner, A., Bechtol, K., Luque, E., et al. 2015, ApJ, 813, 109

Drlica-Wagner, A., Bechtol, K., Mau, S., et al. 2020, ApJ, 893, 47

Garling, C. T., Peter, A. H. G., Kochanek, C. S., Sand, D. J., & Crnojević, D. 2021, MNRAS, 507, 4764

Homma, D., Chiba, M., Okamoto, S., et al. 2016, ApJ, 832, 21

Huxor, A. P., Mackey, A. D., Ferguson, A. M. N., et al. 2014, MNRAS, 442, 2165

Ibata, R. A., Lewis, G. F., McConnachie, A. W., et al. 2014, ApJ, 780, 128

Ibata, R., Martin, N. F., Irwin, M., et al. 2007, ApJ, 671, 1591

Irwin, M. J., Ferguson, A. M. N., Huxor, A. P., et al. 2008, ApJL, 676, L17

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

Jerjen, H., Conn, B., Kim, D., & Schirmer, M. 2018, arXiv:1809.02259

Kim, S. Y., Peter, A. H. G., & Hargis, J. R. 2018, PhRvL, 121, 211302

Klypin, A., Kravtsov, A. V., Valenzuela, O., & Prada, F. 1999, ApJ, 522, 82

Koposov, S. 2009, PhD thesis, Max-Planck Institute for Astronomy

Koposov, S., Belokurov, V., Evans, N. W., et al. 2008, ApJ, 686, 279

Kroupa, P. 2001, MNRAS, 322, 231

Laevens, B. P. M., Martin, N. F., Bernard, E. J., et al. 2015, ApJ, 813, 44

Marigo, P., Girardi, L., Bressan, A., et al. 2017, ApJ, 835, 77

Martin, N. F., Ibata, R. A., Irwin, M. J., et al. 2006, MNRAS, 371, 1983

Martin, N. F., Ibata, R. A., Lewis, G. F., et al. 2016, ApJ, 833, 167

Martin, N. F., Ibata, R. A., McConnachie, A. W., et al. 2013, ApJ, 776, 80

Mashchenko, S., Wadsley, J., & Couchman, H. M. P. 2008, Sci, 319, 174

Mateo, M. L. 1998, ARA&A, 36, 435

McConnachie, A. W., Irwin, M. J., Ferguson, A. M. N., et al. 2005, MNRAS, 356, 979

McConnachie, A. W., Huxor, A., Martin, N. F., et al. 2008, ApJ, 688, 1009

McConnachie, A. W., Irwin, M. J., Ibata, R. A., et al. 2009, Natur, 461, 66

McConnachie, A. W., Ibata, R., Martin, N., et al. 2018, ApJ, 868, 55

Moore, B., Ghigna, S., Governato, F., et al. 1999, ApJL, 524, L19

Mutlu-Pakdil, B., Sand, D. J., Carlin, J. L., et al. 2018, ApJ, 863, 25

Nadler, E. O., Gluscevic, V., Boddy, K. K., & Wechsler, R. H. 2019, ApJL, 878, L32

Richardson, J. C., Irwin, M. J., McConnachie, A. W., et al. 2011, ApJ, 732, 76

Somerville, R. S. 2002, ApJL, 572, L23

Spergel, D. N., & Steinhardt, P. J. 2000, PhRvL, 84, 3760

Thomas, G. F., Martin, N. F., Fattahi, A., et al. 2021, ApJ, 910, 92

Tollerud, E. J., Beaton, R. L., Geha, M. C., et al. 2012, ApJ, 752, 45

Walsh, S. M., Willman, B., & Jerjen, H. 2009, AJ, 137, 450

Weisz, D. R., Martin, N. F., Dolphin, A. E., et al. 2019, ApJL, 885, L8

Wheeler, C., Oñorbe, J., Bullock, J. S., et al. 2015, MNRAS, 453, 1305

Zucker, D. B., Kniazev, A. Y., Bell, E. F., et al. 2004, ApJL, 612, L121

Zucker, D. B., Kniazev, A. Y., Martínez-Delgado, D., et al. 2007, ApJL, 659, L21