Extensions to the halo occupation distribution model for more accurate clustering predictions

TL;DR

This study evaluates various implementations of the halo occupation distribution (HOD) model to improve galaxy clustering predictions, emphasizing the importance of satellite HOD scatter for star formation rate-selected samples.

Contribution

It introduces and tests three HOD schemes with increasing complexity, highlighting the role of satellite HOD scatter in modeling galaxy clustering more accurately.

Findings

All models accurately reproduce the two-halo clustering component.

Differences observed mainly in the one-halo term, especially for SFR-selected samples.

Using a negative binomial distribution for satellite counts yields the best results for SFR selections.

Abstract

We test different implementations of the halo occupation distribution (HOD) model to reconstruct the spatial distribution of galaxies as predicted by a publicly available semi-analytical model (SAM). We compare the measured two-point correlation functions of the HOD mock catalogues and the SAM samples to quantify the fidelity of the reconstruction. We use fixed number density galaxy samples selected according to stellar mass or star formation rate (SFR). We develop three different schemes to populate haloes with galaxies with increasing complexity, considering the scatter of the satellite HOD as an additional parameter in the modelling. We first modify the SAM output, removing assembly bias and using a standard Navarro-Frenk-White density profile for the satellite galaxies as the target to reproduce with our HOD mocks. We find that all models give similar reproductions of the two-halo contribution to the clustering signal, but there are differences in the one-halo term. In particular, the HOD mock reproductions work equally well using either the HOD of central and satellites separately or using a model that also accounts for whether or not the haloes contain a central galaxy. We find that the HOD scatter does not have an important impact on the clustering predictions for stellar mass selected samples. For SFR selections, we obtain the most accurate results assuming a negative binomial distribution for the number of satellites in a halo. The scatter in the satellites HOD is a key consideration for HOD mock catalogues that mimic ELG or SFR selected samples in future galaxy surveys.

Tables

Table 1. Table 1: The first column shows the abundance of galaxies in the three density samples used here. The second and third columns show the cuts applied to G13 galaxies in stellar mass and star formation rate, respectively, to achieve these abundances.

$n/h^3{\mathrm{Mpc}}^{-3}$	$M_{\min }^{\ast }/h^{-1}{M}_{\odot }$	${\mathrm{SFR}}_{\min }/{\mathrm{yr}}^{-1}{M}_{\odot }$
${10}^{-3.0}$	$5.95\times {10}^{10}$	$5.25$
${10}^{-2.5}$	$3.38\times {10}^{10}$	$2.53$
${10}^{-2.0}$	$1.25\times {10}^{10}$	$0.70$

Table 2. Table 2: Satellite fractions of the galaxy samples used. The first column indicates their number densities. Column 2,3,4 and 5 show the satellite fraction in the SAM samples and in the HOD mock built using the 1-HOD, 2-HOD and 4-HOD models, respectively.

		Stellar Mass
$n/h^3{\mathrm{Mpc}}^{-3}$	SAM	1-HOD	2-HOD	4-HOD
${10}^{-3.0}$	$0.280$	$0.317$	$0.279$	$0.280$
${10}^{-2.5}$	$0.322$	$0.405$	$0.324$	$0.324$
${10}^{-2.0}$	$0.381$	$0.517$	$0.383$	$0.383$
		SFR
$n/h^3{\mathrm{Mpc}}^{-3}$	SAM	1-HOD	2-HOD	4-HOD
${10}^{-3.0}$	$0.171$	$0.084$	$0.175$	$0.172$
${10}^{-2.5}$	$0.195$	$0.192$	$0.197$	$0.197$
${10}^{-2.0}$	$0.244$	$0.334$	$0.246$	$0.246$

Equations

$${\sigma }_{\mathrm{NB}}=\sigma +\beta \sigma .$$

$$P(N|r,p)=\frac{\Gamma (N+r)}{\Gamma (r)\Gamma (N+1)}{p}^r(1-p{)}^N.$$

$$p=\frac{⟨N⟩}{{\sigma }_{\mathrm{NB}}^2},r=\frac{{⟨N⟩}^2}{{\sigma }_{\mathrm{NB}}^2-⟨N⟩}.$$

$$⟨{\mathrm{N}}_{sat\_cen}({\mathrm{M}}_{\mathrm{h}})⟩$$

$$⟨{\mathrm{N}}_{sat\_nocen}({\mathrm{M}}_{\mathrm{h}})⟩$$

Full text

††pubyear: 2019††pagerange: Extensions to the halo occupation distribution model for more accurate clustering predictions–11

1 Introduction

In the current cosmological paradigm, the Universe is composed of a filamentary network of structures shaped by gravity. In this framework, dark matter haloes correspond to overdense regions that evolve by gravitational instability due to mergers and interactions with other haloes. Galaxy formation occurs inside haloes where baryons collapse in the gravitational potentials and the condensation of cold gas allows the formation of stars and the evolution of galaxies (White & Rees, 1978). A detailed description of the halo-galaxy connection enables using the galaxies to constrain the cosmological model.

The evolution of dark matter haloes can be followed, to high accuracy, using N-body simulations which use a set of cosmological parameters as inputs. In contrast, the evolution of galaxies in haloes involves many physical processes that are still poorly understood. The fate of baryons within dark matter haloes has been modelled using different approaches, such as hydrodynamical simulations that provide an insight into the formation and evolution of galaxies (e.g. Vogelsberger et al., 2014; Schaye et al., 2015). However, these models are computationally expensive and cannot be run over the large volumes needed for cosmological studies. Alternatively, the effect of baryons can be probed in such large volumes using semi-analytical models (SAMs) of galaxy formation. These start from haloes extracted from a large volume dark matter only simulation and use simplified physical models of the processes that shape the evolution of baryons (Cole et al., 2000; Baugh, 2006; Benson, 2010; Somerville & Davé, 2015). Hence, SAMs make predictions for the abundance and clustering of galaxies that can be compared and tested with large surveys.

Another way to describe the galaxy population is with the halo occupation distribution (HOD) framework (Benson et al., 2000; Peacock & Smith, 2000; Scoccimarro et al., 2001; Yang et al., 2003). This is an empirical approach that provides a relation between the mass of haloes and the number of galaxies hosted by them. This is expressed as the probability distribution $P(N|M_{h})$ that a halo of virial mass $M_{h}$ hosts $N$ galaxies which satisfy some selection criteria. This approach provides insight into the halo-galaxy connection and can be used to study galaxy clustering (Berlind & Weinberg, 2002; Zheng et al., 2005; Conroy et al., 2006; Zehavi et al., 2011; Wechsler & Tinker, 2018). Furthermore, the HOD parameters can be tuned in detail because they only aim to reproduce a limited set of observables such as the galaxy number density and clustering. Thus, HOD modelling is one of the most efficient ways to populate very large volumes or to produce many realizations required for, e.g., estimating covariance matrices using mock galaxy catalogues (e.g. Norberg et al., 2009; Manera et al., 2013). These mock catalogues can then be used to test and develop new algorithms that will be used for the next generation of surveys.

The study of star forming emission line galaxies (ELGs) has gained interest over the last decade as they will be targeted by surveys such as Euclid and the Dark Energy Spectroscopic Instrument (DESI) surveys (Laureijs et al., 2011; DESI Collaboration et al., 2016). The luminosity of an emission line depends on a number of factors, including the star formation rate (SFR), gas metallicity and the conditions in the HII regions (e.g. Orsi et al., 2014). Even though ELGs samples are related to star formation, they are not the same as SFR-selected samples. Still, a similar HOD approach can be used to study both galaxy populations (Geach et al., 2012; Cochrane et al., 2017; Cochrane & Best, 2018). In particular, the shape of the HOD in SFR selected samples is more complex than the case of the more widely studied stellar mass selected samples (e.g. Contreras et al., 2013; Gonzalez-Perez et al., 2018). For example, the occupation function of central galaxies in ELG samples does not follow the canonical step-like form. Accurate modelling of the HOD will provide the more realistic mock catalogues needed for the analysis of future observational samples.

Here, we use the HOD formalism to test three different ways to populate dark matter haloes with galaxies. The prescriptions of these models aim to replicate as accurately as possible the target galaxy populations of a SAM sample. The comparison between the galaxy population in the mock catalogues and SAM samples is done via the analysis of their two-point correlation function (2PCF), which is related to the power spectra of density fluctuations and is sensitive to cosmology (e.g. DeRose et al., 2019). We also include the scatter of the HOD of satellites in our modelling, and quantify the impact of using this additional parameter on the clustering.

The outline of this paper is as follows. The definition of galaxy samples used and the basic properties of the N-body simulation and the SAM are given in Section 2. The correlation functions and the HODs of the samples are presented in Section 3. In Section 4 we introduce the HOD models used to build the mock catalogues and the recipes employed to perform this procedure. The main results and analysis are discussed in Section 5 while in Section 6 we present our conclusions. Appendix A shows the predicted occupation functions for a particular HOD model.

2 Simulation data

In this section we give a brief overview of the galaxy formation model used (§ 2.1) and the N-body simulation in which it is implemented (§ 2.2).

2.1 Galaxy formation model

A galaxy formation model needs to take into account a variety of physical processes such as radiative cooling of gas; AGN, supernovae and photoionisation feedback; chemical evolution; star formation; disc instabilities; collapse and merging of dark matter haloes; and galaxy mergers. These affect the fate of baryons in haloes which lead to the formation and evolution of galaxies. Several physical processes such as star formation and gas cooling are not fully understood due to their complexity. As a consequence, a set of free parameters are used in the equations that model these processes. These free parameters are tuned in order to reproduce observations such as the luminosity functions, colours and the distribution of morphological types. In this context, different SAMs usually have their own implementations to model these physical processes, predicting different galaxy populations. Here we use the outputs at $z=0$ from the SAM of Guo et al. (2013) (hereafter G13) which is a version of the L-GALAXIES code from the Munich group (De Lucia et al., 2004; Croton et al., 2006; De Lucia & Blaizot, 2007; Guo et al., 2011; Henriques et al., 2013). The outputs are publicly available from the Millennium Archive111http://gavo.mpa-garching.mpg.de/Millennium/.

The samples used here are defined according to three different number densities where we rank the galaxies in the SAM by their stellar mass or SFR in a decreasing way (hereafter the SAM samples). These samples are useful in order to compare with observational catalogues with similar space densities. Table 1 shows the three number densities and the cuts in stellar mass and SFR used in each case.

2.2 The Millennium simulation

The distribution of dark matter haloes used in this work is drawn from the Millennium-WMAP7 simulation Guo et al. (2013) which is identical with the Millennium Simulation Springel et al. (2005), but with updated cosmological parameters that match the results from the WMAP7 observations. This version assumes a flat $\Lambda$CDM universe considering $\Omega_{m}=0.27$, $\Omega_{\Lambda}=0.73$, $h=H_{0}/(\rm 100\ km\ s^{-1}\ Mpc^{-1})=0.704$ and $\sigma_{8}=0.81$. The simulation was carried out in a box-size of $500\ h^{-1}\rm Mpc$ following $2160^{3}$ particles of mass $9.31\times 10^{8}h^{-1}\rm M_{\odot}$. The run produced 61 simulation snapshots from $z=50$ up to $z=0$. G13 use a friends-of-friends group finding algorithm (FOF) to identify dark matter haloes in each snapshot (Davis et al., 1985) and then run SUBFIND to identify the subhaloes (Springel et al., 2001). Halo merger trees are constructed for each output and track the evolution of haloes through cosmic time. These trees are the starting points for the SAM.

3 Characterization of the SAM galaxy samples

This section introduces the statistics used to characterize the distribution of galaxies, starting with the measurement of the correlation function (§ 3.1), the form of the HOD predicted by the SAM (§ 3.2) and the scatter in the HOD (§ 3.3).

3.1 Clustering measurement: two-point galaxy correlation function

The spatial two-point correlation function, $\xi(r)$, measures the excess probability of finding a pair of galaxies at a given separation with respect to a random distribution. We compute the 2PCF of the galaxy samples with the Corrfunc code (Sinha & Garrison, 2017).

Fig. 1 shows the 2PCF of the stellar mass (top) and SFR (bottom) selected samples for the different space densities. For the former, the amplitude of the clustering increases with decreasing number density, as we consider more massive galaxies. The impact of the inclusion of these massive galaxies is stronger at small scales and is weaker at large scales. In contrast, for the SFR selected galaxies the amplitude of the 2PCF for the different samples remains largely unchanged except for small scales where the satellite-satellite pairs make an important contribution to the clustering amplitude. For both selections, the satellite fraction increases with increasing number density.

In the 2PCF, we can distinguish between the contribution from galaxy pairs in the same halo and from different haloes. The former are the main contributors to the amplitude of the 2PCF on small scales, namely the one-halo term which dominates up to $\sim 1\ h^{-1}{\rm Mpc}$, while galaxy pairs between different haloes contribute mostly to the two-halo term which determines the clustering on large scales. In this regime the total number of galaxies in the halo, regardless of whether they are satellites or the central, drives the amplitude of the clustering, acting as a weighting for the bias of each halo in computing an overall “effective” bias for the sample (see e.g. Baugh et al. 1999). The one-halo term is sensitive not only to the number of satellites, but also depends on their spatial distribution.

3.2 The halo occupation function predicted by the SAMs

The galaxy populations in the SAMs depend on the choices adopted for the modelling of the baryonic processes. Hence, depending on the SAM employed, different galaxy catalogues with different luminosity functions, stellar mass functions or correlation functions can be obtained for the same dark matter simulation. For example, Contreras et al. (2013) studied the effects on the clustering predicted from different SAMs and found some differences particularly in galaxy samples selected by SFR and cold gas mass. Moreover, they show that the shapes of the HODs are model-dependent which reflects the differences in the implementation of physical processes in each SAM. For example, the specific modelling of dynamical friction affects the satellite population in SAMs. Here we are not interested in the detailed shape of the HOD predicted by a particular SAM, but on how best to use the occupation functions to populate dark matter haloes with galaxies to produce a similar spatial distribution to that resulting from a SAM.

The HOD is usually broken down into the contribution from central and satellite galaxies. Fig. 2 shows these two components for stellar mass and SFR selected samples with the same number density for the G13 SAM. Here, each HOD is computed in bins of width 0.08 dex in the logarithm of the halo mass where the position of each $\langle N\rangle$ value is plotted at the median value within each bin. The striking difference in the shape of the HOD of centrals between the two selections is due to the different galaxies that are included. Massive centrals tend to be red galaxies hosted by massive dark haloes. Such centrals are included in stellar mass selected samples but not when selecting by SFR. The galaxies in the SFR samples correspond mainly to blue star-forming galaxies excluding luminous red galaxies with high stellar mass but low SFR. It is noteworthy that the fraction of halos that contain a central passing the SFR selection never reaches unity for the sample plotted in Fig. 2. These features of the HOD of SFR galaxies have been noted in SAMs before (e.g. Contreras et al., 2013, 2019; Gonzalez-Perez et al., 2018) and inferred for blue galaxies in the SDSS (e.g. Zehavi et al., 2011). This shows that a significant number of haloes in the SFR selected samples do not host a central as their SFR is below the threshold. The same situation is found in the other number density samples. Note that, in observational samples, the ranking of galaxies in order of their emission line luminosity may not correspond to the ranking in SFR due to dust attenuation, which means that the highest SFR galaxies may not necessarily have the brightest emission lines.

We estimate the uncertainties of the HOD values using jackknife resampling (Norberg et al., 2009), dividing the simulation volume into 10 slices. We use the position of the centre of the potential of haloes to classify the galaxies within each halo. The resulting errors are shown as the shaded regions in Fig. 2, and they are negligible for all halo masses except at the high mass end and for the HOD of centrals selected by SFR.

Because of the simple relation between halo mass and occupation number, the HOD represents a useful approach for the construction of mock galaxy catalogues. Here we have described the first moment of the HOD, the main ingredient to building-mocks recipes. Nevertheless, it is important to also consider the second moment, i.e the dispersion in the HOD of satellites.

3.3 The predicted dispersion in the halo occupation number

When the simplest HOD approach is used to build mock catalogues, the mean of the distribution is the main parameter. Central galaxies are assumed to follow a nearest-integer distribution where the mean $\rm\langle N_{cen}\rangle$ is between zero and one. For satellites, a Poisson distribution with mean $\rm\langle N_{sat}\rangle$ is the most widely assumed distribution (e.g Kravtsov et al., 2004; Zheng et al., 2005).

In G13, satellites are classified as type-1 if they are hosted by a resolved subhalo, and type-2 or orphans if the subhalo has been destroyed by tidal effects and is no longer identified. Boylan-Kolchin et al. (2010) found that the number of low mass subhaloes in main haloes in the Millennium-II Simulation (Boylan-Kolchin et al., 2009) is well described by a negative binomial distribution which corresponds to a super-Poissonian statistic as its scatter is larger than a Poisson distribution. This suggests that the type-1 satellite population can also be described by this distribution. Based on the outputs from the SAM presented in Jiang & van den Bosch (2016), and using the Bolshoi (Klypin et al., 2011) and MultiDark (Prada et al., 2012) simulations, Jiang & van den Bosch (2017) showed that ignoring this non-Poissonity in the HOD of subhaloes results in systematic errors in the predicted clustering of galaxies. Here we extend the application of the negative binomial distribution by checking whether the HOD of G13 satellites, that is including type-1 and type-2, is well described by this statistic. We expect that the HOD scatter is model-dependent because of the different treatments of dynamical friction. Moreover, as some galaxy properties, such as SFR and stellar mass, have a model-dependent scatter it is reasonable to assume the same for HODs. For example, Guo et al. (2016) showed that different galaxy formation models do not have the same dispersion in the stellar mass-halo mass relation. Therefore our results are specific to the G13 model. It is likely that a different SAM would require an adjustment to to the value of $\beta$ to describe the scatter of the satellite HOD. Nevertheless we expect our general results to hold for any SAM, and that the satellite distribution displays more scatter than Poisson.

The Poisson and negative binomial distributions differ in their shapes so it is useful to parametrise the departure from the Poisson scatter. We use the parameter $\beta$ (defined below) to denote this departure. For a Poisson distribution the variance is given by the mean value of the random variable, namely $\rm\langle N_{sat}\rangle$, with the standard deviation given by $\sigma=\sqrt{\rm\langle N_{sat}\rangle}$. The negative binomial distribution has the same mean as the Poisson distribution, but a larger scatter which can be expressed as

[TABLE]

where $0<\beta<1$. Then, $\beta$ indicates the fractional change in the variance with respect to the Poisson standard deviation $\sigma$. Under this definition, when $\beta=0$ the distribution is Poissonian and if $\beta=1$ the standard deviation is twice that from a Poisson distribution.

The probability function of the negative binomial distribution is given by

[TABLE]

Here $\Gamma(x)=(x+1)!$ is the gamma function. The parameters $r$ and $p$ are determined by the first moment $\left<N\right>$ and second moment $\sigma^{2}$ of the distribution,

[TABLE]

Thus, we can control the width of the negative binomial distribution through the parameter $\beta$ and compute the value of $\sigma_{\rm NB}^{2}$.

4 Generating HOD mock catalogues

We now describe the procedure followed to build HOD mock galaxy catalogues using the HODs of the SAM samples. Section 4.1 presents the three methods we use to populate haloes with galaxies. In Section 4.2 we specify the treatment of the scatter in the HOD of satellites. Section 4.3 explains how we impose a standard Navarro-Frenk-White (NFW) density profile for satellites. Section 4.4 presents the impact of assembly bias in the SAM samples and explains why it must be removed from the SAM in order to compare with the HOD mock catalogues. Finally in Section 4.5 we discuss the treatment of the radial distribution of satellite galaxies within haloes.

4.1 The HOD models used to build mocks

We test three different HOD schemes of increasing complexity. This helps us to understand the level of complexity needed to obtain accurate clustering predictions. Each model uses occupation functions obtained from linear interpolations of the HOD values in each bin, rather than fitting a parametric form to the SAM HOD. The distribution of galaxies can be nearest-integer (centrals only) and Poisson or negative binomial (satellites).

4.1.1 1-HOD

The 1-HOD model builds mock catalogues using the HOD of all galaxies from the SAM sample (black solid lines in Fig. 2) including both centrals and satellites. The model assumes either a Poisson or negative binomial distribution for the occupation number. We adopt a Monte Carlo approach to obtain the final number of galaxies.

This approach does not distinguish between centrals and satellites. If the model predicts that $\rm N\geq 1$ we assume that this halo hosts a central and $\rm N_{sat}=N-1$. Because of this, the number of centrals and satellites in the 1-HOD mock catalogues can be notably different with respect to the SAM samples where there are haloes with satellites but no central. Moreover, the HODs of these two separate components in the mock catalogues are completely different with respect to the HODs of the SAM samples (see Appendix A). However, the total number of galaxies in these mock catalogues is essentially the same as in the SAM samples.

4.1.2 2-HOD

The 2-HOD model uses the HOD of centrals and satellites separately, i.e, the red and blue solid lines in Fig. 2, respectively. Thus, a particular distribution can be assumed for each component and the modelling is done independently for each one. For centrals, we use the nearest-integer distribution and for satellites the Poisson or negative binomial distribution.

This scheme predicts practically the same number of central and satellites as the SAM samples. Note that in a non-negligible number of realizations it is possible to get haloes without a central. This is more likely for haloes with masses for which $\rm\left\langle N_{cen}\right\rangle<1$ which is more frequently found in SFR-selected samples.

4.1.3 4-HOD

This model contains more information about the galaxy population of the SAM sample than the 2-HOD model. The 4-HOD requires us to store the number of haloes that host a central ($\rm N_{cen}$) and the number of haloes that do not host a central ($\rm N_{nocen}$) as a function of halo mass. Under this definition the total number of haloes in the volume is the sum of both quantities. Furthermore, the 4-HOD also needs knowledge of the number of satellites in haloes with a central ($\rm N_{sat\_cen}$) and without a central ($\rm N_{sat\_nocen}$). Thus, the total number of satellites is the sum of these two quantities. With these definitions, we build new HODs for satellites that take into account the population of centrals in the SAM samples. The SAM samples contain haloes with satellites but no centrals. This is more common in SFR selected samples. Indeed the HOD of centrals in these samples indicates that a large number of haloes do not host a central (see Fig 2), and the 4-HOD takes this feature into account.

We then define the satellite occupation functions conditioned on whether or not haloes host a central. With the four quantities explained above, we can define the conditional HODs,

[TABLE]

Fig. 3 shows the conditional HODs where the main differences are observed at low halo masses. Even though the ratio between these two HODs is close to unity, it is the galaxies hosted by these haloes ($\approx 10^{12}h^{-1}{\rm M_{\odot}}$) that dominate the amplitude of clustering. The conditional HODs are well fitted by a negative binomial distribution, including the HODs of the other number density samples. The 4-HOD method uses a Monte Carlo approach to decide if a halo hosts a central galaxy. Depending on this outcome, one of the two conditional HODs is then chosen to obtain the number of satellites.

4.2 Treatment of scatter in the HOD of satellites

A Poisson distribution is fully described by its first moment. In the case of satellites this is $\rm\langle N_{sat}\rangle$. If the distribution of the number of satellites follows instead a negative binomial distribution, an additional parameter $\beta$ is needed which specifies the increase in the scatter with respect to a Poisson distribution (see Eq. 1). We fix the $\beta$ value so that we reproduce as closely as possible the scatter of the HOD of satellites in a given SAM sample. Fig. 4 shows the scatter of the HOD of satellites in SAMs and 2-HOD mock catalogues for two illustrative $\beta$ values in a SFR selected sample. This shows that a small but non-zero $\beta$ is required to reproduce the HOD scatter of the SAM sample. The same is found for the other number density samples, and for the conditional HODs. We do not perform this analysis for the 1-HOD model as satellites are not treated independently in this case.

It is not possible to replicate the HOD scatter in SAMs more closely as this would require $\beta$ to be a function of mass. Instead, we assume a constant scatter for the HOD of satellites by using the same $\beta$ for all halo mass bins. The accuracy of the $\beta$ values used are judged by checking the quality of the resulting mocks via comparison of their 2PCFs with the clustering of the shuffled-NFW samples (see § 4.5 below for the definition of this catalogue). We show in Section 5, that when the scatter of the SAM and HOD mocks are matched up to $M_{h}\lesssim 10^{13.5}M_{\odot}h^{-1}$ ($\beta=0.05$ for SFR-selected samples), we obtain the most accurate clustering predictions. In contrast, using larger values for $\beta$ worsens the predictions (as does using $\beta=0$, which corresponds to Poisson scatter).

The satellite HOD is well described by the negative binomial distribution for a wide range of halo masses. Fig. 5 shows the satellite PDF in a particular mass bin for a stellar mass and a SFR selected sample. We show negative binomial distributions defined by $\beta=0.08$ and $\beta=0.05$. In order to compute the satellite distributions, we split satellites according to whether or not their haloes host a central galaxy, which is relevant for the 4-HOD model. The satellite distribution matches with the negative binomial when most of haloes in the bin are included. A similar close match is found when comparing with Poisson distributions ($\beta=0$). Note that in the SFR selection case most of the haloes do not host a central galaxy, as the HOD of centrals in that bin suggests. The opposite behaviour is observed when selecting by stellar mass.

4.3 The radial distribution of satellite galaxies in halos

The number of satellites in the mock catalogues is obtained from the adopted HOD model (see Sec 4.1). Their positions in haloes are set according the standard NFW density profile (Navarro et al., 1996) which requires two parameters, the concentration and scale radius. The former depends on halo mass and the latter is a function of the virial radius. For simplicity, we assume that all haloes in the simulation volume have the same concentration parameter $c=13.98$ which corresponds to the concentration of a halo at redshift $\rm z=0$ and mass $M_{h}=10^{12.5}{\rm M_{\odot}}h^{-1}$. We do not use a more realistic model for concentration as we are interested in comparing the HOD models rather than obtain a realistic redistribution of satellites. We impose that the maximum distance from a satellite to the halo center is two virial radii which depends on the halo mass. This defines the NFW mass profile used to obtain the satellite distances by a Monte Carlo approach. We modify the SAM output to impose a similar satellite distribution as described below (§ 4.5).

4.4 Removing assembly bias from the SAM output

In order to determine the best methodology to produce HOD mock catalogues, we aim to compare them with the clustering of the original SAM samples via their 2PCF. Before making this comparison it is necessary to remove assembly bias from the SAM samples.

The clustering of dark matter haloes depends on additional properties besides mass. For example, Gao et al. (2005) showed that the clustering of low mass haloes depends on their formation redshift and other works have found dependencies on concentration and subhalo occupation number (e.g. Wechsler et al., 2006) among other secondary properties. This additional contribution to the clustering is commonly known as assembly bias and potentially changes the galaxy clustering amplitude on large scales.

The standard HOD approach considers only halo mass as the variable regulating the galaxy population. SAMs include assembly bias because they follow the evolution of baryons in halo merger histories that are shaped by the large-scale environment in the N-body simulation. Namely, SAMs include a dependence on secondary halo properties as these affect the halo merger history and the evolution of galaxies that live within them. Thus, in order to compare the clustering between SAM samples and HOD mocks which use only halo mass as input, it is necessary to remove the assembly bias signal in the former samples.

Assembly bias can be eliminated from SAM samples through the shuffling technique introduced by Croton et al. (2007). This consists of randomly exchanging the galaxy populations between haloes of the same mass, thus removing any connection to the assembly history of the haloes. This procedure does not change the distances from satellites to their central galaxy in each halo. In clustering terms, the one-halo term of this “shuffled” catalogue is the same as the original SAM sample but its two-halo term is different because assembly bias is not present in the shuffled sample. If the SAM samples did not have assembly bias, we would measure the same 2PCFs for their shuffled samples as measured for the original output.

Fig. 6 shows the correlation functions of a SAM sample and its shuffled version, for both the stellar mass and the SFR selected samples. The assembly bias signature, shown in the middle panels, is evident in the clustering differences between these two catalogues at large separations. We also show the 2PCF of a modified shuffled sample that will be introduced below. The assembly bias signatures remain unchanged for the other samples, but they are noisier for the lowest number density samples as they contain fewer galaxies.

It can be seen that assembly bias increases the clustering for stellar mass selected samples, as was shown by Zehavi et al. (2018). SFR selected samples, on the other hand, show a decreased clustering amplitude. For the intermediate galaxy density sample, the assembly bias enhances the two-halo term by $\sim 12\%$ for the stellar mass selected sample and suppresses the amplitude in the SFR selection case by $\sim 4\%$. The enhance of clustering amplitude for the other stellar mass selections remains similar. For the SFR selections, we observe that suppression on clustering becomes weaker for higher density samples. Indeed, assembly bias can enhance the amplitude if the density of the sample if very high, as shown in (Contreras et al., 2019).

4.5 The shuffled-NFW target catalogue: changing the satellite distribution in the SAM

The radial profile of satellites in the G13 SAM deviates from the standard NFW profile of dark matter within halos because the SAM associates galaxies with subhalos (or a proxy, such as the most bound particle, in the case of subhalos which are no longer resolved). The radial profile of subhalos is different from that of the dark matter (see e.g. Angulo et al. 2009). The choice of which subhalos (and former subhalos) are associated with galaxies is driven by the galaxy formation model, which determines the luminosity of any galaxy associated with a subhalo and whether or not it has merged due to dynamical friction (only type 2 satellites, those that no longer have a resolved subhalo associated with them, are considered as candidates for galaxy mergers).

The final step before testing the accuracy of the HOD models is to modify the shuffled SAM catalogue to force the satellites in each halo to follow an NFW profile. We call the result the shuffled-NFW catalogue. Because satellite galaxies in the SAMs and shuffled samples do not follow an NFW profile, the one-halo term of their 2PCFs are different to the one-halo term of the shuffled-NFW sample, as shown in the bottom panel of Fig. 6. The shuffled-NFW catalogue does not contain assembly bias, and the satellites follow the same NFW profile as adopted in the HOD mocks. We note that the shuffled-NFW is not intended to be the “best” prediction of galaxy clustering but rather is the target sample for the reconstructions using the HOD mocks which has a controlled 1-halo clustering pattern to facilitate testing.

5 Testing the accuracy of the HOD models

5.1 Satellite radial distributions and clustering of HOD mocks

In Fig. 7 we show the satellite profiles in the SAM samples, for stellar mass and SFR selected samples separated into the contributions from type-1 satellites and from orphan galaxies. To examine the departure from NFW, we produce a SAM-NFW catalogue where satellites in SAM are forced to follow the same NFW profile as used in the HOD mocks. For this catalogue, we update the satellite positions in the SAM samples according the same NFW density profile used to produce the HOD mocks. Note that this SAM-NFW is different from the shuffled and shuffled-NFW catalogues mentioned above.

It can be seen that the NFW profile is different from the profile of type-1 satellites and orphans, particularly for the SFR selected sample. The profiles in the 2-HOD and 4-HOD mock catalogues are also shown, and they match with the NFW as expected from the construction of the HOD mock. Both models also reproduce the NFW profile in the other number density samples.

The masses of host haloes of 1-HOD satellites do not correspond with the masses in the original SFR selected samples (see Fig. 11). Thus, the virial radii of these haloes define NFW density profiles that are different from the profiles in the other models. This has an impact on the positions of satellites generating the striking difference with NFW in Fig. 7 for the SFR selection. The same occurs for the stellar mass selections but it is less extreme than the SFR case.

The HOD models predict different galaxy populations for the G13 SAM samples. Table 2 shows the satellite fraction of the SAM samples and HOD mocks built by the three different models, assuming a Poisson distribution for the HOD of satellites. Note that the 2-HOD and 4-HOD predict almost the same satellite fraction as in the SAM samples because of the separate HOD modelling of centrals and satellites.

We check the accuracy of each HOD model by comparing the HOD mocks with the shuffled-NFW sample via their 2PCFs. Fig. 8 shows the clustering of the shuffled-NFW and the HOD mocks built using the HOD models described in Sec. 4.1. These particular models assume a Poisson distribution for the HOD of satellites. It can be seen that the three schemes produce accurate clustering predictions on large scales. On small scales, the 2-HOD and 4-HOD models produce similar accurate results while the 1-HOD shows striking differences. These deviations come from the overprediction of the number of satellites in the stellar mass selected samples. For the SFR selection cases, the difference is due to the notably different occupation function of central and satellites in the 1-HOD mock (see Appendix A). The inaccuracy of the 1-HOD modelling is also present for the other number density samples too, whereas the 2-HOD and 4-HOD models produce similar quality results to the one shown here. As the 2-HOD and 4-HOD models are clearly the best, we drop the 1-HOD model henceforth.

5.2 Impact of the assumed HOD scatter

To study the impact of the scatter of the HOD on the clustering, we consider different dispersions for the negative binomial distribution in the construction of HOD mocks. Fig. 9 shows the 2PCF of HOD mocks using different $\beta$ values. For the stellar mass selected samples, the scatter of the HOD does not have a significant impact.

For the SFR selection, the amplitude of the clustering on small scales is very sensitive to the scatter in the number of satellites. We find that increasing $\beta$ changes the amplitude of the one-halo term. When we split the contribution to the clustering from low and high mass haloes, we observe that the scatter mainly impacts the one-halo term of low mass haloes. This feature is present in HOD mocks built using the 2-HOD and 4-HOD methods. This is reproduced by both HOD models, indicating that this is a feature particular to the SFR selected samples.

The most accurate clustering reconstructions for the G13 samples are obtained when we use the 2-HOD or 4-HOD to build mock catalogues assuming a negative binomial distribution for the HOD of satellites. Note that clustering predictions from both models do not show significant differences.

Fig 10 shows the particular results from the 4-HOD modelling for all the space density samples. For the G13 SFR selected samples, the 4-HOD (and the 2-HOD) modelling produces the best results when $\beta=0.05$, which corresponds to a distribution slightly wider than a Poisson distribution. For the case of stellar mass selected samples, the best reproduction is obtained with $\beta=0.08$. Using instead the Poisson distribution (i.e, $\beta=0$) produces worse results for both selections particularly in the one halo regime.

For SFR selections, when using $\beta=0.05$ and $\beta=0$, the departures from the shuffled-NFW catalogues are below $\sim 8\%$ and $\sim 15\%$, respectively. It can be seen that the dispersion of the 2PCFs becomes important in the lowest number density sample. However, the assumption of the negative binomial distribution still produces better results, especially in the transition from the one- to the two-halo term.

For the stellar mass selection cases, the impact on clustering when using different $\beta$ values is much less significant. Indeed, the weak relation between clustering and HOD scatter, shown in Fig. 9, suggests that it is not necessary to include additional scatter in the construction of HOD mock catalogues for stellar mass selections. To compare with the Poisson distribution, we show also the clustering prediction for the stellar mass samples using $\beta=0.08$, for which we obtain the best result.

Satellites in the G13 SAM are well described by a non-Poisson distribution. This is consistent with the HOD of subhaloes found in Boylan-Kolchin et al. (2010). The recipes that build mock catalogues of SFR selected samples using the HOD approach must undertake an analysis of the scatter of the HOD of satellites as it impacts the clustering. This analysis will provide the best $\beta$ to construct a HOD mock of a particular sample. For stellar mass samples, the HOD scatter has a weak impact on clustering, so the same analysis is not necessary in the context of HOD mock catalogues.

6 Conclusion

The next generation of surveys will measure the clustering of the galaxy distribution over a wide range of redshifts. Mock catalogues have proven to be important tools in preparation for this because of their multiple applications including error analysis, data interpretation and survey planning. SAMs are a physical approach to obtain such mocks, but sometimes their direct application to a simulation is not possible due to the limited resolution of the halo merger trees (Angulo et al., 2014) or the trees may not be available, as in the case of the Euclid flagship simulation (Potter et al., 2017). Even if the trees were available at the required resolution, the sheer number of halos in a giga-parsec side N-body simulation may preclude a direct calculation with a SAM.

The HOD model provides a simple yet efficient way to construct mock catalogues. The modelling consists of using a probability distribution to obtain the number of galaxies hosted in a halo of a particular mass. This simple method allows us to create large sets of mock catalogues for huge cosmological simulations. This is useful in the context of ELGs, as they are targets in current and coming surveys.

To determine the level of complexity needed to produce accurate mock catalogues, we test different HOD models. The 1-HOD uses the HOD of all galaxies making no distinction between centrals and satellites while the 2-HOD uses the HOD of these two components separately. The 4-HOD stores additional information about whether or not haloes host a central, and it constructs conditional HODs for satellites taking into account this information.

Because SAMs include assembly bias by construction, and in their simplest form HOD mocks do not, we remove the assembly bias from the G13 SAM samples by shuffling the galaxy populations among haloes of the same mass creating the shuffled catalogue. This allows us to make a direct comparison between the clustering of our mocks and the SAMs from which we extract the HOD measurements. For example, We find that, for the intermediate galaxy density sample in the G13 SAM, the assembly bias affects the 2-halo term of the 2PCF of stellar mass selected galaxies increasing the amplitude by $\sim 12\%$. For the SFR selected galaxies, in contrast, the assembly bias suppresses the clustering by $\sim 4\%$. We also impose the standard NFW profile for satellites in the shuffled catalogue as is done for satellites in the HOD mock catalogues. We can then check the accuracy of the HOD models through a comparison between the 2PCFs of the HOD mocks and the shuffled-NFW catalogues.

The 2-HOD and 4-HOD produce the best mock catalogues as their 2PCFs are in close agreement with the clustering of the shuffled-NFW sample. We obtain the best results using a negative binomial distribution for the (conditional) HOD (see Eq. 1); in previous works this was commonly considered to be a Poisson distribution. This is consistent with the subhalo HOD found in Boylan-Kolchin et al. (2010) using the Millennium-II simulation. Furthermore, we found that the assumption of this non-Poisonnian HOD changes the galaxy clustering. Previously, Jiang & van den Bosch (2017) found a similar result using subhaloes from the Bolshoi, the MultiDark N-body simulations, and the SAM presented in Jiang & van den Bosch (2016).

The scatter of the HOD of satellites in G13 is reproduced by a negative binomial distribution up to halo masses of ${\rm M_{h}\lesssim 10^{13.5}M_{\odot}}h^{-1}$. The galaxies in this halo mass range dominate the amplitude of the 2PCF. We quantify the departure from the Poisson distribution with the parameter $\beta$ (see Eq. 3). We obtain the best clustering predictions for SFR selected samples using $\beta=0.05$ and $\beta=0.08$ for stellar mass selected samples. These correspond to negative binomial distributions slightly wider than Poisson. Because of the specific modelling of different SAMs, we expect that the best $\beta$ values for each sample are model-dependent. For stellar mass selected samples, we find that the HOD scatter has a weak impact on clustering, making the addition of this additional parameter unnecessary in the context of mock catalogues.

The analysis of the HOD of satellites is important because the width of the distribution (determined by the $\beta$ parameter) has a large impact on the one-halo term of the 2PCF of mock catalogues that emulate SFR-selected sample and ELG samples. If we consider the Poisson distribution for the HOD of satellites ($\beta=0$) the 2PCF of the mock catalogues is underestimated with respect to the clustering of the shuffled-NFW. In contrast, using the negative binomial distribution increases the amplitude of clustering in the one-halo regime. If we assume a value of $\beta$ larger than the one present in the distribution of number of satellites, the clustering on small scales is further overestimated. We highlight the importance to perform a careful analysis of the satellite HOD if the HOD framework is used to produce mock catalogues, for ELGs or star forming galaxies, following a particular model or observation.

Acknowledgements

This work was made possible by the efforts of Gerard Lemson and colleagues at the German Astronomical Virtual Observatory in setting up the Millennium Simulation database in Garching. EJ, SC, NP & IZ acknowledge the hospitality of the ICC at Durham University. EJ acknowledges support from “Centro de Astronomía y Tecnologías Afines” BASAL 170002. NP acknowledges support from Fondecyt Regular 1191813. IZ acknowledges support by NSF grant AST-1612085. This project has received funding from the European Union’s Horizon 2020 Research and Innovation Programme under the Marie Skłodowska-Curie grant agreement No 734374. The calculations for this paper were performed on the Geryon computer at the Center for Astro-Engineering UC, part of the BASAL PFB-06, which received additional funding from QUIMAL 130008 and Fondequip AIC-57 for upgrades.

Appendix A 1-HOD occupation functions

The 1-HOD model uses the HOD of all SAM galaxies in a sample to produce mock catalogues. This does not necessarily reproduce the HOD of central and satellites separately.

The 1-HOD assumes a distribution for the full HOD, which we take to be Poisson or a negative binomial. However, centrals follow the nearest integer distribution, whereas satellites follow a Poisson or negative binomial distribution. As this distinction is not made, the predicted HODs of centrals and satellites in the resulting mock are different with respect to the original SAM samples. Fig. 11 shows that 1-HOD tends to put satellites in very low mass haloes. Moreover, the occupation function of centrals in the SFR selection is overpredicted over a wide halo mass range.