The spatial and velocity bias of linear density peaks and protohaloes in the cold dark matter cosmology
Abstract
We use high resolution Nbody simulations to investigate the Lagrangian bias of cold dark matter haloes within the cold dark matter cosmology. Our analysis focuses on ‘protohaloes’, which we identify in the simulation initial conditions with the subsets of particles belonging to individual redshiftzero haloes. We then calculate the numberdensity and velocitydivergence fields of protohaloes and estimate their auto spectral densities. We also measure the corresponding cross spectral densities with the linear matter distribution. We use our results to test a Lagrangianbias model presented by Desjacques and Sheth which is based on the assumption that haloes form out of local density maxima of a specific height. Our comparison validates the predicted functional form for the scaledependence of the bias for both the density and velocity fields. We also show that the bias coefficients are accurately predicted for the velocity divergence. On the contrary, the theoretical values for the density bias parameters do not accurately match the numerical results as a function of halo mass. This is likely due to the simplistic assumptions that relate virialized haloes to density peaks of a given height in the model. We also detect appreciable stochasticity for the Lagrangian density bias, even on very large scales. These are not included in the model at leading order but correspond to higher order corrections.
keywords:
methods:analytical – numerical – galaxies: haloes – cosmology: theory – largescale structure of Universe.1 Introduction
Galaxy redshift surveys are powerful probes of cosmology. The main observables able to constrain cosmological parameters are the overall shape of the galaxy power spectrum at wavenumbers and the baryonic acoustic oscillations within it. These are treated as proxies for the matter power spectrum for which we can make robust theoretical predictions. Galaxies, however, are biased tracers of the cosmic mass distribution and many features appearing in their power spectrum depend on how a specific observational sample was selected. To reconstruct the matter power spectrum we thus need an accurate bias model whose free coefficients should be used as nuisance parameters and marginalized over. In the era of precision cosmology, where measurements of the matter power spectrum with per cent accuracy are required, this task is particularly demanding.
Bias models can be divided into two broad classes. Eulerian biasing schemes relate the galaxy density contrast, , to the matter density distribution, , evaluated at the same time (but not necessarily at the same spatial location). After smoothing the fields on large scales, so that is typically much smaller than unity, one can write (Fry & Gaztañaga, 1993)
where all fields are evaluated at the same time and the details of the bias model are specified by the kernel functions, . If one further assumes that biasing is local (i.e. that all kernels can be written as products of Dirac delta distributions), this reduces to
(2) 
where now the bias coefficients are real numbers.
In Lagrangian bias models, on the other hand, one considers the regions in the initial conditions that will collapse to form galaxies (or their hosting darkmatter haloes) at time and writes their density contrast, , in terms of the linear density contrast, . Largescale expansions analogous to eqs. (1) and (2) can also be written in this case. As a second step, one must determine the final position, , of a fluid element initially located at , and compute out of . This LagrangiantoEulerian mapping (LEM) accounts for gravitationally induced motions that determine the final position of the objects.
Local Eulerian and Lagrangian bias schemes are not equivalent; they generate a different shape for the galaxy bispectrum (Catelan et al., 2000) and are not compatible within the framework of perturbation theory (Matsubara, 2011). In fact, Catelan et al. (1998) have shown that a local Lagrangian biasing scheme generates a nonlinear, nonlocal and stochastic bias in Eulerian space. On the other hand, nonlocal Eulerian and Lagrangian schemes are equivalent and can be seen as different mathematical representations of the same physical process (Matsubara, 2011).
Due to its simplicity, the local Eulerian model is by far the most widely used in practical applications, such as perturbative calculations. However, it is purely phenomenological and does not have a strong theoretical motivation. Detailed comparison with numerical simulations has also evidenced its limited validity (e.g. Roth & Porciani 2011, Pollack, Smith & Porciani 2011). Physical models of bias are generally given in the Lagrangian framework as conditions on galaxy (halo) formation are more easily imposed onto the linear density field using some model for the collapse of density perturbations. Mo & White (1996, hereafter MW) used a PressSchechterlike argument (Press & Schechter, 1974) to compute the bias coefficients of a local Lagrangian scheme as a function of halo mass. The same authors also showed how these parameters can be combined to calculate the bias coefficients of a local Eulerian scheme assuming that largescale density perturbations follow the spherical collapse model. The effect of nonlinear shear on the LEM was discussed by Catelan et al. (1998) using the Zel’dovich approximation (Zel’dovich, 1970). For halo masses (where is the characteristic mass for collapse at redshift ) the MW formula for is in good agreement with the predictions of Nbody simulations (Mo, Jing & White, 1997). At lower masses, however, the agreement rapidly deteriorates (Jing, 1998). Porciani, Catelan & Lacey (1999) and Jing (1999) showed that the discrepancy between the Nbody simulations and the analytical predictions is already present in Lagrangian space and should thus be attributed to the limitations of the PressSchechter formalism rather than to the approximated treatment of the LEM.
The Lagrangian bias emerging from the extended PressSchechter model (Bond et al., 1991) was first derived by Porciani et al. (1998), rediscussed in Scannapieco & Barkana (2002) and tested against simulations by Scannapieco & Thacker (2005). This approach follows correlated trajectories of at different Lagrangian locations as a function of the smoothing scale and looks for correlations in the firstcrossing scales of a density threshold.
According to the peakbackgroundsplit argument (Bardeen et al. 1986; Cole & Kaiser 1989), longwavelength density fluctuations modulate halo formation by modifying the collapse time of localized shortwavelength perturbations. This makes it possible to generalise the calculation of the Lagrangian MW bias coefficients to any model for the halo mass function (Sheth & Tormen 1999, Tinker et al. 2010, Giannantonio & Porciani 2010) and to improve the agreement with Nbody simulations (Seljak & Warren 2004, Tinker et al. 2005, Gao et al. 2005, Pillepich et al. 2010).
Kaiser (1984) explained the strong clustering of Abell clusters by assuming that they originate from the regions above a density threshold in the (suitably smoothed) linear density field. Following this line of reasoning, it is common to assume that darkmatter haloes form out of linear density peaks, as an alternative to the PressSchechter approach. Tests against Nbody simulations have shown that this is a well justified assumption, especially for massive haloes (Ludlow & Porciani 2011, see also Frenk et al. 1988). The statistical properties of the local maxima in a Gaussian random field have been extensively studied by Bardeen et al. (1986) (see also Peacock & Heavens 1985 and Hoffman & Shaham 1985). Mo, Jing & White (1997) introduced peaks theory in the MW formalism, while Matsubara (1999) evaluated the level of stochasticity in the Lagrangian clustering of density extrema.
Recently, Desjacques (2008) and Desjacques & Sheth (2010, hereafter DS) showed that the correlation of density peaks in real and redshift space can be interpreted in terms of a simple Lagrangian biasing scheme. Due to the peak constraint, the effective Lagrangian peak density at a given point not only depends on the local value of the mass density but also on its Laplacian. The firstorder peak bias depends on the mass and height of the peaks, and on the matter power spectrum. For high peaks, this reduces to the results by Matsubara (1999). Moreover, although peaks move with the matter at their positions, DS infer the existence of a statistical velocity bias due to the fact that local maxima can only exist at special locations. This also leads to a bias between the linear velocity spectra of peaks and matter which is predictable in quantitative terms.
In spite of the fact that the DS model provides the initial conditions for sophisticated models of the (Eulerian) halo distribution where the LEM is based either on resummed perturbation theory (Elia et al., 2011) or on the Zel’dovich approximation (Desjacques et al., 2010), its predictions for the Lagrangian clustering and velocities of the regions that will form collapsed structures have never been thoroughly tested against numerical simulations. This paper provides such a test, which is necessary if we are to use advanced bias models to extract useful information on the cosmological parameters through a comparison with observations.
The structure of the paper is as follows. In Section 2 we review the bias model first presented in DS. The details of our Nbody simulations and analysis techniques are outlined in Section 3, and a comparison between the model and numerical results is presented in Section 4. Finally, we summarize our main results in Section 5.
2 The DS model
In this section we summarize, for completeness, the peaks model described in Section 2 of DS. In order to do so, we first introduce and define some quantities relevant for peak statistics.
The spectral moments of the matter power spectrum are defined as
(3) 
where is the linear matter power spectrum at redshift , and is a smoothing kernel of characteristic length . In terms of the moments we define the spectral parameters:
(4) 
There are two common choices for : the Gaussian filter
(5) 
and the tophat filter
(6) 
The mass contained within the comoving length in these cases is
(7) 
where is the mean matter density of the Universe.
We will often characterize peaks in terms of their dimensionless peak height, , defined
(8) 
Here is the smoothed overdensity at the peak location linearly extrapolated to , and is the linear rms mass fluctuation in spheres of radius . DS linked density maxima of height to dark matter haloes of mass collapsing at redshift , assuming that coincides with the threshold for collapse, .
Bardeen et al. (1986) and DS computed the crosscorrelation between peaks and the underlying density field which also corresponds to the average density profile around density maxima. Similarly, Desjacques (2008) evaluated the leading order expressions (on large spatial separations) for the peak autocorrelation function and the lineofsight mean streaming for pairs of discrete local maxima of height . Desjacques (2008) and DS showed that the full expression for the crosscorrelation and the largescale asymptotic of the autocorrelation function are consistent with – and thus can be thought of as arising from – an effective bias relation. In the DS model, the number density and the velocity of peaks, and , are related to the dark matter density contrast and velocity fields, linearly extrapolated to , via:
(9) 
and
(10) 
The subscript “S” indicates that the fields are smoothed on the scale , and the bias parameters, and , are given by
(11) 
and
(12) 
Here is the mean curvature of the peaks, which can be approximated by (Bardeen et al., 1986)
(13) 
Note that coincides with the peak bias factor found by Bardeen et al. (1986) after neglecting the derivatives of the density correlation function.
Since, by definition, the gradient of the density field vanishes at peak locations, eq. (10) suggests that the peak and dark matter velocity fields must be coincident there. By construction, peaks move with the dark matter flow, yet the spatial bias induces a statistical velocity bias. We will consider the scaled velocity divergence, , rather than the velocity field. Here is the scale factor, the Hubble parameter, , with the linear growth factor. In these units, both and are dimensionless quantities. With these changes, eqs. (9) and (10) can be rewritten in Fourier space as
(14) 
and
(15) 
where we have defined
(16) 
In the limit of high peaks () it can be shown that the bias parameters obey the following asymptotic relations: and . This implies that the highest peaks are linearly biased tracers of the underlying matter field. This is consistent with the predictions of the peakbackground split (Mo, Jing & White, 1997) and DS showed that, indeed, is the appropriate generalization of the constant, largescale bias for low . Unlike the density bias factors, does not depend on .
In order to test this model against Nbody simulations, we will make use of the crossspectra between the peak and dark matter densities and velocities (denoted as and , respectively) and of the corresponding peak autospectra ( and ). From eqs. (14) and (15) we obtain
(17) 
(18) 
where and are the matter density and velocity divergence autospectra, respectively. We remind the reader that the expression for is only valid to first order in as , and that higher order corrections should be included to improve its accuracy (see e.g. Desjacques et al. 2010). On the contrary, the expression for the crossspectrum is exact, as shown in Bardeen et al. (1986) and in the Appendix A of DS. Eq. (17) has the same functional form as eq. (57) in Matsubara (1999) who studied the clustering of density extrema. The DS coefficients and match those in Matsubara (1999) only in the limit , for which nearly all extrema are density maxima.
3 Numerical Issues
In this section we provide a brief description of the main numerical issues relevant for this work. This includes a brief description of our numerical simulations in Section 3.1, our main analysis techniques in Section 3.2, and a characterization of halo collapse barriers in Section 3.3.
3.1 Nbody simulations
Our analysis focuses on two highresolution Nbody simulations of structure formation in the standard LCDM cosmology. The cosmological parameters for our runs were chosen to be consistent with the fifthyear WMAP data release (Komatsu et al., 2009). These are , , , , and . Each simulation was run with a lean version of the TreePM code Gadget2 (Springel, 2005), and followed the dark matter using collisionless particles. One simulation had a box sidelength of Mpc and a particle mass of ; the other used Mpc and had . The initial redshifts of the simulations were and for the larger and smaller box, respectively. Using these simulations we are able to probe a wide range of halo masses, spanning . These simulations were first studied in Pillepich et al. (2010), and later by Ludlow & Porciani (2011), and we refer the reader to those papers for further details.
Haloes were identified at using a friendsoffriends (FOF) algorithm (Davis et al., 1986) with a linking length of 0.2 times the mean interparticle distance. Protohaloes were identified by tracing backward to the initial redshift of all subsets of the particles belonging to FOF haloes. We use the centre of mass of each protohalo as a proxy for its spatial location; the massweighted linear velocity provides an estimate of the protohalo’s motion. As a test of the sensitivity of our results to the adopted halo finder, we also generated a spherical overdensity (SO) halo catalogue with an overdensity threshold of 200 times the critical density, . For a fixed halo mass, the two halofinders produce results consistent within per cent (in terms of all the bias coefficients) and so, in what follows, we will focus on results obtained for the FOF haloes, and only consider those containing at least 100 particles.
Haloes in each simulation are split into four separate mass bins in order to preserve their peculiar clustering properties. These bins are referred to as 1S to 4S for the small box, and 1L to 4L for the large one. To asses the impact of shot noise in the analysis of our smallbox simulation we consider an additional mass bin, labeled bin0S, which includes all haloes with . The mass ranges and total number of haloes in each bin are given in Table 1. We note that bins with label ‘S’ refer to masses for which darkmatter haloes are not expected to be in onetoone correspondence with linear peaks (Ludlow & Porciani, 2011).
3.2 Analysis
We construct protohalo density and momentum fields using cloudincell grid assignment on a mesh. Velocity fields are obtained by taking the ratio of the momentum and density fields, as described in Scoccimarro (2004). In the case of haloes, these distributions are smoothed to preclude the existence of empty cells; the smoothing scales used are Mpc for the large box and Mpc for the small one. These values are chosen in order to mask the effects of the grid, but we have explicitly verified that our results are not significantly affected by them. All power spectra have been computed using a fast Fourier transform technique.
The discreteness of darkmatter particles and haloes gives rise to a shotnoise component in the spectra. For the density fields, the estimated power spectrum, , includes a shotnoise term which is inversely proportional to the number density of objects (assuming Poisson sampling):
(19) 
Shot noise is therefore negligible for the matter spectra but may be significant for that of the haloes. The issue is more severe in Lagrangian space, because fluctuations in the initial conditions are small. We will consider two alternative estimates of the protohalo bias; one is determined from the shotnoise corrected autospectrum,
(20) 
and the other from the crossspectrum,
(21) 
Here the subscript “h” indicates the halo fields, and “m” the matter field (the analogous fields for peaks are indicated with the subscript “p”.) The relation between the two is
(22) 
where is the linear correlation coefficient, defined as
(23) 
These relations tell us that the two definitions of the bias are equivalent only if the bias is purely deterministic in Fourier space, i.e. . At leading order, in the model presented by DS, and (neglecting the filter function). Any stochasticity (represented by the higher order corrections in ) will degrade the correlation, yielding different estimates for the bias. Because of this, we will consider both the cross and autospectra to check for a potential stochastic element of the bias. ^{1}^{1}1Although the effective firstorder bias model by DS is deterministic in Fourier space, it is stochastic in configuration space (Matsubara 1999, Desjacques & Sheth 2010).
As for the density, we can define two estimates for the velocity bias, that we denote and ,
with a correlation coefficient .
There is no conclusive way to subtract the shot noise for the velocity divergence spectrum.
Hence we used our simulations to gain some insight into this issue. We introduced an artificial shot noise in the matter spectrum
by randomly drawing a fraction of the particles and found that asymptotically for large . Therefore we fitted the amplitude factor to obtain shotnoise corrected power spectra. Note however that other terms could be important at smaller wavenumbers; in this work we only consider data for which .
The DS model describes the biasing of linear density peaks that are expected to form haloes of a given mass at a specified redshift according to some collapse model (which determines the value for ). However, we analyse the density and velocity fields for the actual FOF and SO protohaloes. These are the quantities of physical interest for studying galaxy clustering in terms of halo occupation models (e.g. Cooray & Sheth 2002). Ludlow & Porciani (2011) showed that the vast majority of haloes in our Nbody simulations can be unambiguously associated with linear peaks in the initial conditions when smoothed on the mass scale of the halo. For example, per cent of all haloes can be matched with similarmass peaks in , and per cent for haloes with . Note, however, that the correspondence between the DS peaks and protohaloes will not be perfect.^{2}^{2}2The measured abundance of haloes in the simulations is of the same order of magnitude as (albeit a bit higher than) the number of peaks with computed as in Bardeen et al. (1986). On a given mass scale, some peaks with overdensities at the collapse threshold will evolve into substructures contained within larger virialized haloes (the socalled cloudincloud problem), or to haloes of significantly different mass. On top of this, numerical simulations have shown that the collapse threshold for haloes of a given mass and redshift has a broad probabilistic distribution rather than a fixed value (Porciani, Dekel & Hoffman 2002; Robertson et al. 2009) and possibly also depends on the form of the adopted smoothing kernel. We consider some of these issues in the following section.
3.3 Barrier heights for tophat and Gaussian filters
The peak model described in Section 2 requires welldefined spectral moments. However, due to its sharp boundary in real space, the tophat filter decays very slowly in Fourier space so that the integral defining is divergent in a LCDM model. Since eqs. (11), (12) and (13) depend on , DS instead adopted a Gaussian filter and assumed . Here, corresponds to the critical density for the collapse of a spherical tophat perturbation in an otherwise unperturbed EdS universe, and this does not necessarily apply to peaks in a smoothed Gaussian random field. Another issue is that the validity of the simple spherical collapse model is questionable, at best; the probability of a protohalo or a peak being spherical is null, since it would require the three eigenvalues of the tidal tensor being equal. Sheth, Mo & Tormen (2001) showed that the barrier height in the more general ellipsoidal collapse model (Bond & Myers, 1996) can be approximated by
(24) 
where is taken from the spherical collapse model, and and are determined from fits to the model results. The presence of the dispersion, , in eq. (24) results in a massdependent barrier height: lower mass haloes require, on average, higher overdensities for collapse since they must hold themselves together against larger tidal forces. It should be emphasized that eq. (24) describes the mean barrier height; the scatter about the mean can be approximated by (Robertson et al., 2009). These values are valid only for the tophat filter.
What is the appropriate barrier height corresponding to a Gaussian filter? In Figure 1 we show the linear overdensities measured at protohalo centres of mass after smoothing with a Gaussian filter on the halo mass scale. The shaded regions show the halo data, and connected circles the medians of the distribution. Open points correspond to results from our Mpc box simulation, and solid points to our Mpc box run. Squares show the median for the same sample of haloes, but after smoothing the linear density field with a tophat filter instead (note the good agreement with the result of Sheth, Mo & Tormen (2001), shown as a dotdashed line in Figure 1). All values have been linearly extrapolated to .
Bin  Mass range  haloes  
()  ()  ()  ()  ( Mpc)  
SG  EG  SG  EG  
0S  0.31  0.23  0.63  0.61  0.92  0.62  
1S  0.30  0.21  0.26  0.25  0.45  0.42  
2S  0.30  0.22  0.33  0.32  0.55  0.46  
3S  0.31  0.23  0.49  0.48  0.76  0.56  
4S  0.31  0.24  0.92  0.90  1.24  0.74  
1L  0.05  0.08  11.3  11.5  8.75  2.26  
2L  0.18  0.003  14.7  15.2  10.8  2.54  
3L  0.40  0.12  20.5  21.5  14.0  2.96  
4L  0.92  0.44  34.5  37.1  21.2  3.77 
It is clear that adopting a Gaussian filter entails lower values of at all masses^{3}^{3}3At any mass scale, , the smoothing length of a Gaussian filter exceeds that of a tophat filter by a factor of .. Nonetheless, our results are still accurately described by eq. (24), albeit with slightly different values for the model parameters: , and . We show this explicitly in Figure 1 using a dashed line. The onesigma dispersion about the mean trend is well described by a simple powerlaw, , as seen in the upper panel of the plot.
4 Results
In this section we test the predictions of the DS model against the measured density and velocity bias of dark matter protohaloes. In order to do so, we calculate the (modelpredicted) bias parameters , and for each mass bin in both of our simulations. Since the values of and depend on peak height (and hence smoothing scale) we adopt two models for the collapse barrier: one assumes the spherical collapse model (hereafter SG) and the other the ellipsoidal collapse model (hereafter EG). A Gaussian filter is used in both cases. Table 1 lists the bias parameters predicted by the DS model for these choices of collapse barrier.
Bin  ()  

FOF  SO  FOF  SO  FOF  SO  FOF  SO  
0S  0.25  0.24  0.31  0.30  0.86  0.81  0.14  0.11 
1L  0.19  0.22  0.27  0.37  8.8  7.9  12.6  12.0 
2L  0.29  0.31  0.40  0.39  10.3  9.3  9.5  8.8 
3L  0.43  0.47  0.52  0.54  13.9  13.1  6.7  5.9 
4L  0.73  0.75  0.73  0.73  21.8  20.9  6.4  5.4 
4.1 Density spectra
In Figure 2 we plot and extracted from our simulations,^{4}^{4}4The values directly measured from the simulations are rescaled by the growth factor to match the theoretical estimates, where is linearly extrapolated to . Note that the actual Lagrangian bias is a factor larger than the values reported in the figure. along with the model prediction (eq. (17)). Although the model is not expected to work for very small scales, we nonetheless show the results for each box up to (we remind that is the smoothing scale needed for the velocity field). Overall, the model expression for the initial peak bias is able to describe the simulation results reasonably well. This can be seen from the solid lines in Figure 2, where and have been treated as free parameters and fitted to both the cross and the auto spectra. However, it is clear that the values for the fit parameters are different in the two cases (see Table 2), implying or, equivalently, . In particular, for the range of masses in Bin 0S, is negative, i.e. . The predictions from the SG and EG barriers, as apparent from the dashed and dotted lines in the right panel in Figure 2, are also negative. Hence the bias model matches more closely rather than . This is completely expected, because eq. (17) is exact only for the crosscorrelation between peaks and matter while it neglects higherorder corrections for the peak autocorrelation. However, neither the SG nor the EG barrier provide the appropriate values for the coefficients. In particular, the EG barrier performs better for low masses (bin 0S), while the situation is reversed for the higher mass bins.
Figure 2 suggests that the stochasticity is more of a problem for haloes with , and on smaller scales. Since highmass haloes are highly correlated with density peaks in the initial conditions (Ludlow & Porciani, 2011) and peaks follow eq. (14), the relationship between and is likely to be more deterministic for massive haloes. Note, however, that the estimate of is affected by the shotnoise correction, which is large in our samples. Therefore we cannot draw definitive conclusions regarding stochasticity.
The performance of the model is characterized in another way in Figure 3, where we compare directly the mass dependence of the bias parameters and predicted by the SG and EG models to the bestfit values obtained using the parameterization of the crossspectrum in eq. (17). We do not show the corresponding results from the autospectrum because of the uncertain shotnoise subtraction. Figure 3 reveals that overall the values derived from the SG model are closer to the fitted ones. In contrast, there is almost no difference between the two models for . For both parameters the mass dependence follows a different trend when we compare the fits to the models. In particular, the prediction for and becomes more and more inaccurate with decreasing halo mass (up to a factor of for ). On the other hand, in the mass range for which the model was developed (), is overpredicted by a factor of , independent of halo mass. Figure 2 shows that the disagreement is even worse when and are obtained from . All of this demonstrates that eqs. (11) and (12) cannot accurately describe the mass dependence of the Lagrangian density bias, especially for . We will discuss possible reasons for this in Section 5, but first turn our attention to the velocity spectra.
4.2 Velocity spectra
Figure 4 shows and for bins 1S4S and 1L4L; the theoretical model, , is also plotted. The data have been “desmoothed”, i.e. they have been divided by a filter function with smoothing scale , the same used to smooth the velocity field in the first place. It is apparent that on large scales the two estimates of the bias coincide, i.e. , indicating a strong correlation between the fields. Overall, we can conclude that the velocity bias is deterministic to a good approximation.
There is an excellent agreement (better than per cent) between the simulation results and the model predictions for all mass bins at scales for the large box and for the small one. Therefore, the peak model provides a faithful description of the velocity bias.
5 Summary
We have investigated the Lagrangian bias of darkmatter haloes by testing the theoretical model proposed by DS against Nbody simulations of structure formation. The model assumes that haloes form from density peaks and predicts a scaledependent bias for both the density and velocity fields. Our main results can be summarized as follows.

When averaged over a spherical Lagrangian volume containing the appropriate mass, the linear density contrast measured at protohalo centres depends sensitively on the choice of smoothing kernel. For a Gaussian kernel, for example, the resulting density contrasts are systematically lower than those computed using a tophat filter. This is because, at fixed mass, the smoothing length of a Gaussian kernel exceeds that of a tophat filter by a factor of about resulting in systematically lower density estimates. Nonetheless, the median barrier height computed with a Gaussian kernel can be accurately parameterized by the same fitting formula first advocated by Sheth, Mo & Tormen (2001) for the case of a tophat filter, albeit with different values for the numerical parameters. We use this result to approximate the collapse threshold for dark matter halo formation when adopting a Gaussian filter.

The functional forms for the density and velocity bias relations derived by DS  our eqs. (17) and (18)  accurately describe the results obtained from our simulations, provided the parameters of the model are allowed to vary with respect to the modelpredicted values. In both cases, the Lagrangian bias is characterized by a constant term that dominates on very large scales, and a scaledependent term proportional to .

Quantitatively, the velocity bias predicted by the DS model is able to reproduce the measured protohalo velocity bias in our simulations to better than per cent, provided we limit ourselves to quasilinear scales (). The predicted density bias (eqs. (11) and (12)), on the other hand, does not match the Lagrangian density bias extracted from the simulations. This is likely due to the more complex nature of the density bias, which depends additionally on peak height and on the exact definition of a protohalo. These results are independent of whether one adopts a barrier height consistent with either the spherical or ellipsoidal collapse model.

We have measured the mass dependence of the density bias, and , by fitting the density powerspectra obtained for haloes in several different mass bins in each of our two simulations. The most massive haloes identified in our simulations are the most strongly biased, and have characteristic overdensities that correlate with the underlying matter density. In contrast, the distribution for low mass haloes () exhibits an anticorrelation. A similar trend is also predicted by the DS model, although with noticeable differences (see Figure 3). We emphasize that the bestfit values for and obtained in our analysis apply for a LCDM cosmogony (with our adopted cosmological parameters) and only over the limited mass range probed by our simulations. Future work should consider how the Lagrangian density bias depends on the underlying cosmology.

Our comparison of the auto and crossspectra in Section 4 suggests the presence of a stochasticity in the Fourierspace density bias of dark matter protohaloes, but none for the velocity bias. This is in disagreement with the DS model for the bias at leading order, and seems to corroborate the need for higherorder terms in the expression for the autospectrum of the protohaloes, as predicted in Desjacques et al. (2010). However, we cannot draw firm conclusions regarding stochasticity in the bias estimates due to uncertainties in the shotnoise subtraction.
We have tested the Lagrangian bias model of DS against a pair of high resolution Nbody simulations of structure formation and found that it is able to accurately reproduce the velocity bias, but not that of the density. Our analysis focused on protohaloes, the high redshift progenitors of darkmatter haloes, whereas the model describes the biasing of density peaks. One possible explanation for the differences in the model predictions and simulation results comes form the differences between the expressions for the bias parameters: , and . The latter of the three is solely determined by the linear matter power spectrum and the smoothing scale corresponding to a given mass, . The density bias parameters, however, depend additionally on explicit properties of the peaks, such as their height, , mean curvature, , and on an assumed collapse threshold for their identification as haloes of a given mass. Figure 1 shows that there is a large halotohalo variation in at any given mass scale; characterizing the collapse threshold as a single massdependent value may, therefore, be too simplistic. This added complexity introduces a significant margin for error in the model’s estimates of the density bias. It remains to be seen whether more realistic models of the collapse barrier  such as those that account for the statistical scatter in the linear overdensities of protohaloes at a given mass  will improve the model’s predictive power.
Another possibility for the discrepancies stems from the fact that we are identifying linear density peaks with protohaloes in our simulation initial conditions. Although the majority of our dark matter protohaloes form in the vicinity of linear density peaks of the same characteristic mass, the fate of all peaks with the same mass and overdensity is unclear. Uncertainties associated with the identification of protohaloes within the linear density field may have adverse effects on the predictive power of the density bias model. A more detailed understanding of how protohaloes map onto linear density peaks, and vice versa, will, no doubt, provide valuable insight into the mechanisms behind halo biasing.
Acknowledgments
AE acknowledges financial support from the SFBTransregio 33 “The Dark Universe” by the Deutsche Forschungsgemeinschaft (DFG) and a scholarship from the BonnCologne Graduate School (BCGS). We acknowledge interesting discussion on halo finders with Steffen Knollmann and Andrea Macciò. We would like to thank Jeremy Tinker for providing us with his SO halo finder, and for help applying it to our simulations. We are very grateful to the referee, Vincent Desjacques, for useful suggestions and valuable insights.
References
 Bardeen et al. (1986) Bardeen J.M., Bond J.R., Kaiser N., Szalay A.S., 1986, ApJ, 304, 15
 Bond et al. (1991) Bond J.R., Cole S., Efstathiou G., Kaiser N., 1991, ApJ, 379, 440
 Bond & Myers (1996) Bond J.R., Myers S.T., 1996, ApJ Suppl., 103, 1
 Catelan et al. (1998) Catelan P., Lucchin F., Matarrese S., Porciani C., 1998, MNRAS, 297, 692
 Catelan et al. (2000) Catelan P., Porciani C., Kamionkowski M., 2000, MNRAS, 318, L39
 Cole & Kaiser (1989) Cole S., Kaiser N., 1989, MNRAS, 237, 1127
 Cooray & Sheth (2002) Cooray A., Sheth R.K., 2002, PhR, 372, 1
 Davis et al. (1986) Davis M., Efstathiou G., Frenk C. S., White S. D. M., 1985, ApJ, 292, 371
 Desjacques (2008) Desjacques V., 2008, Phys. Rev. D, 78, 103503
 Desjacques & Sheth (2010) Desjacques V., Sheth R.K., 2010, Phys. Rev. D, 81, 023526
 Desjacques et al. (2010) Desjacques V., Crocce M., Scoccimarro R., Sheth R.K., 2010, Phys. Rev. D, 82, 103529
 Elia et al. (2011) Elia A., Kulkarni S., Porciani C., Pietroni M., Matarrese S., 2011, MNRAS, 416, 1703
 Frenk et al. (1988) Frenk C. S., White S. D. M., Davis M., Efstathiou G., 1988, ApJ, 327, 507
 Fry & Gaztañaga (1993) Fry J.N., Gaztañaga E., 1993, ApJ, 413, 447
 Kaiser (1984) Kaiser N., 1984, ApJ, 284, L9
 Komatsu et al. (2009) Komatsu E., 2009, ApJ Suppl., 180, 330
 Gao et al. (2005) Gao L., White S.D.M., Jenkins A., Frenk C.S., Springel V., 2005, MNRAS, 363, 379
 Giannantonio & Porciani (2010) Giannantonio T., Porciani C., 2010, Phys. Rev. D, 81, 063530
 Hoffman & Shaham (1985) Hoffman Y., Shaham J., 1985, ApJ, 297, 16
 Jing (1998) Jing Y.P., 1998, ApJ, 503, L9
 Jing (1999) Jing Y.P., 1999, ApJ, 515, L45
 Ludlow & Porciani (2011) Ludlow A., Porciani C., 2011, MNRAS, 413, 1916
 Matsubara (1999) Matsubara T., 1999, ApJ, 525, 543
 Matsubara (2011) Matsubara T., 2011, Phys. Rev. D, 83, 083518
 Mo & White (1996) Mo H.J., White S.D.M., 1996, MNRAS, 282, 347
 Mo, Jing & White (1997) Mo H.J., Jing Y.P., White S.D.M., 1997, MNRAS, 284, 189
 Peacock & Heavens (1985) Peacock J.A., Heavens A.F., 1985, MNRAS, 217, 805
 Pillepich et al. (2010) Pillepich A., Porciani C., Hahn O., 2010, MNRAS, 402, 191
 Pollack, Smith & Porciani (2011) Pollack J.E., Smith R.E., Porciani C., 2011, arXiv 1109.3458
 Porciani, Catelan & Lacey (1999) Porciani C., Catelan P., Lacey C., 1999, ApJ, 513, L99
 Porciani, Dekel & Hoffman (2002) Porciani C., Dekel A., Hoffman Y., 2002, MNRAS, 332, 339
 Porciani et al. (1998) Porciani C., Matarrese S., Lacey C., Lucchin F., Catelan P., 1998, MNRAS, 298, 109
 Press & Schechter (1974) Press W. H., Schechter P., 1974, ApJ, 187, 425
 Robertson et al. (2009) Robertson B.E., Kravtsov A.V., Tinker J., Zentner A.R., 2009, ApJ, 696, 636
 Roth & Porciani (2011) Roth N., Porciani C., 2011, MNRAS, 415, 829
 Scannapieco & Barkana (2002) Scannapieco E., Barkana R., 2002, ApJ, 571, 585
 Scannapieco & Thacker (2005) Scannapieco E., Thacker R.J., 2005, ApJ, 619, 1
 Scoccimarro (2004) Scoccimarro R., 2004, Phys. Rev. D, 70, 083007
 Seljak & Warren (2004) Seljak U., Warren M.S., 2004, MNRAS, 355, 129
 Sheth & Tormen (1999) Sheth R.K., Tormen G., 1999, MNRAS, 308, 119
 Sheth, Mo & Tormen (2001) Sheth R.K., Mo H.J. & Tormen G., 2001, MNRAS, 323, 1
 Springel (2005) Springel V., 2005, MNRAS, 364, 1105
 Tinker et al. (2005) Tinker J.L., Weinberg D.H., Zheng Z., Zehavi I., 2005, ApJ, 631, 41
 Tinker et al. (2010) Tinker J., Robertson B.E., Kravtsov A.V., Klypin A., Warren M.S., Yepes G., Gottlöber S., 2010, ApJ, 724, 878
 Zel’dovich (1970) Zel’dovich Ya.B., 1970, A & A, 5, 84