Measuring non-local Lagrangian peak bias

Abstract

We investigate non-local Lagrangian bias contributions involving gradients of the linear density field, for which we have predictions from the excursion set peak formalism. We begin by writing down a bias expansion which includes all the bias terms, including the non-local ones. Having checked that the model furnishes a reasonable fit to the halo mass function, we develop a one-point cross-correlation technique to measure bias factors associated with χ²-distributed quantities. We validate the method with numerical realizations of peaks of Gaussian random fields before we apply it to N-body simulations. We focus on the lowest (quadratic) order non-local contributions |$-2\chi _{10}(\mathrm{\boldsymbol k}_1\cdot \mathrm{\boldsymbol k}_2)$| and |$\chi _{01}[3(\mathrm{\boldsymbol k}_1\cdot \mathrm{\boldsymbol k}_2)^2-k_1^2 k_2^2]$|⁠, where |$\mathrm{\boldsymbol k}_1$|⁠, |$\mathrm{\boldsymbol k}_2$| are wave modes. We can reproduce our measurement of χ₁₀ if we allow for an offset between the Lagrangian halo centre-of-mass and the peak position. The sign and magnitude of χ₁₀ is consistent with Lagrangian haloes sitting near linear density maxima. The resulting contribution to the halo bias can safely be ignored for M = 10¹³ M_⊙ h⁻¹, but could become relevant at larger halo masses. For the second non-local bias χ₀₁ however, we measure a much larger magnitude than predicted by our model. We speculate that some of this discrepancy might originate from non-local Lagrangian contributions induced by non-spherical collapse.

cosmology: theory, dark matter, large-scale structure of Universe

1 INTRODUCTION

Understanding the clustering of dark matter haloes has been a topic of active research for many years. A number of analytic approaches have been developed to tackle this issue such as the peak model (Bardeen et al. 1986, hereafter BBKS), the excursion set framework (Bond et al. 1991) or perturbation theory (see e.g. Bernardeau et al. 2002, for a review). Heuristic arguments like the peak-background split (Kaiser 1984) and approximations like local bias (Fry & Gaztanaga 1993) have been very helpful for modelling the clustering of dark matter haloes. Nevertheless, improvements in computational power and numerical algorithms as well as the advent of large-scale galaxy surveys have considerably increased the need for an accurate description of halo clustering. Until recently, however, it was unclear how the peak approach, which is thus far the only framework in which biased tracers form a discrete point set, relates to the more widespread excursion set theory, local bias approximation or peak-background split argument.

Working out this connection has been the subject of several recent papers. Desjacques (2013), building on earlier work by Desjacques et al. (2010), showed that correlation functions of discrete density peaks can be computed using an effective (i.e. which does not involve measurable counts-in-cells quantities) generalized bias expansion in which all the bias parameters, including those of the non-local terms,¹ can be computed from a peak-background split. In parallel, Paranjape & Sheth (2012) demonstrated how the peak formalism, which deals with statistics of density maxima at a fixed smoothing scale, can be combined with excursion set theory, whose basic building block is the density contrast at various filtering scales. Similar ideas can already be found in the early work of Bond (1989). Paranjape, Sheth & Desjacques (2013a, hereafter PSD) subsequently computed the mass function and linear bias of haloes within this excursion set peak (ESP) approach and showed that it agrees very well with simulation data.

The focus of this work is on the second-order non-local bias terms predicted by the ESP approach. These generate corrections to the Fourier peak bias of the form |$-2\chi _{10}(\mathrm{\boldsymbol k}_1\cdot \mathrm{\boldsymbol k}_2)$| and |$\chi _{01}[3(\mathrm{\boldsymbol k}_1\cdot \mathrm{\boldsymbol k}_2)^2-k_1^2 k_2^2]$| (Desjacques 2013). What makes them quite interesting is the fact that there are related to χ² rather than normally distributed variables. Here, we will show how one can measure their amplitude in the bias of dark matter haloes without computing any correlation function. Of course, this technique can also be applied to measure non-local Lagrangian bias contributions induced by e.g. the tidal shear, but this will be the subject of future work.

This paper is organized as follows. In the first part, we will advocate a slight modification of the original ESP formulation of PSD in order to easily write down the corresponding effective bias expansion (Section 2). Next, we will explain how the cross-correlation technique proposed by Musso, Paranjape & Sheth (2012), which has already been successfully applied to the bias factors associated with the density field (PSD; Paranjape et al. 2013b), can be extended to measure the second-order non-local bias factors χ₁₀ and χ₀₁ that weight the two quadratic, non-local bias contributions (Section 3). Finally, we will validate our method with peaks of Gaussian random fields before measuring χ₁₀ and χ₀₁ for dark matter haloes (Section 4). We conclude in Section 5.

2 EXCURSION SET PEAKS

In this section, we apply the excursion set approach to the peak model in the case of a moving barrier to get a prediction of the halo mass function which we compare to simulations. We then get expressions for bias parameters, generalizing results in Desjacques (2013) and Desjacques, Gong & Riotto (2013). We also point out a few changes to PSD. We show that, as far as the mass function is concerned, these modifications do not make much difference (only few per cent, in agreement with what PSD found), but they affect first- and second-order bias parameters, as new terms arise.

2.1 Notation

We will adopt the following notation for the variance of the smoothed density field (linearly extrapolated to present day) and its derivatives:

\begin{equation} \sigma ^2_{j\alpha } = \frac{1}{2\pi ^2}\int _0^\infty \!\! \mathrm{d}k\, P(k) k^{2(j+1)} W_\alpha ^2(kR_\alpha ), \end{equation}

(1)

where P(k) is the power spectrum of the mass density field, W_α(kR_α) and the subscript α = G or T will denote Gaussian or tophat filtering, respectively. Moreover, R_α is the Lagrangian smoothing scale (which may depend on the choice of kernel). Denoting δ_T and δ_G the linear density field smoothed with a tophat and Gaussian filter, respectively, we introduce the variables

\begin{eqnarray} \nu (\mathrm{\boldsymbol x}) &=& \displaystyle\frac{1}{\sigma _{0{\rm T}}}\delta _{\rm T}(\mathrm{\boldsymbol x}) \nonumber \\ u(\mathrm{\boldsymbol x}) &=& -\frac{1}{\sigma _{2{\rm G}}}\nabla ^2\delta _{\rm G}(\mathrm{\boldsymbol x}) \nonumber \\ \mu (\mathrm{\boldsymbol x}) &=& -\frac{{\rm d}\delta _{\rm T}}{{\rm d}R_{\rm T}}(\mathrm{\boldsymbol x}). \end{eqnarray}

(2)

Note that, while ν and u have unit variance, μ is not normalized. We will use the notation |$\langle \mu ^2 \rangle = \Delta ^2_0$| in what follows.

Cross-correlations among these three variables are useful and will be denoted as

\begin{equation} \langle \nu u \rangle = \gamma _1 = \frac{\sigma _{1X}^2}{\sigma _{0T}\sigma _{2{\rm G}}} \end{equation}

(3)

\begin{equation} \langle \nu \mu \rangle = \gamma _{\nu \mu } = \frac{1}{\sigma _{0T}}\int _0^\infty \!\! \frac{{\rm d}k}{2\pi ^2}\,P(k) k^2 W_{\rm T}(kR_{\rm T})\frac{{\rm d}W_{\rm T}(kR_{\rm T})}{{\rm d}R_{\rm T}} \end{equation}

(4)

\begin{equation} \langle u\mu \rangle = \gamma _{u\mu } = \frac{1}{\sigma _{2{\rm G}}}\int _0^\infty \!\! \frac{{\rm d}k}{2\pi ^2}\, P(k) k^4 W_{\rm G}(kR_{\rm G})\frac{{\rm d}W_{\rm T}(kR_{\rm T})}{{\rm d}R_{\rm T}} . \end{equation}

(5)

The first-order, mixed spectral moment σ_1X is

\begin{equation} \sigma ^2_{1X} = \frac{1}{2\pi ^2}\int {\rm d}k P(k) k^4 W_{\rm T}(kR_{\rm T})W_{\rm G}(kR_{\rm G}) , \end{equation}

(6)

i.e. one filter is tophat and the other Gaussian.

2.2 First-crossing and moving barrier

2.2.1 Summary of previous results

Let us first summarize the basic ideas behind the ESP approach introduced by Paranjape & Sheth (2012) and further developed in PSD and Desjacques et al. (2013).

The excursion set approach states that a region of mass M has virialized when the overdensity δ(R), where R ∼ M^1/3 is the filtering scale associated with the perturbation, reaches the spherical collapse threshold δ_c provided that, for any R^′ > R, the inequality δ(R) < δ_c holds. This last condition formally implies an infinite set of constraints (one at each smoothing scale). However, as was shown in Musso & Sheth (2012), the requirement δ(R + ΔR) < δ_c with ΔR ≪ 1 furnishes a very good approximation. This follows from the fact that the trajectory described by δ(R) as a function of R is highly correlated for large radii. As a result, if δ crosses δ_c at R, then it is almost certainly below the threshold at any larger radius.

This first-crossing condition can be combined with the peak constraint, so that peaks on a given smoothing scale are counted only if the inequality above is satisfied. In this case, the effective peak bias expansion introduced in Desjacques (2013) is modified through the presence of a new variable μ (equation 2) which, as was shown in Desjacques et al. (2013), reflects the dependence of bias to the first-crossing condition.

2.2.2 Modifications to PSD

We made a couple of modifications to the approach of PSD, which we will now describe in more detail.

First, PSD used the fact that μ ≡ u when Gaussian filtering is also applied to the density field, so that the first-crossing condition can be accounted for with the variable u only. When δ is smoothed with a tophat filter however, one should in principle deal explicitly with μ and, therefore, consider the trivariate normal distribution |${\cal N}(\nu ,u,\mu )$|⁠. We will proceed this way.

Secondly, Sheth, Mo & Tormen (2001) argued that, owing to the triaxiality of collapse, the critical density for collapse is not constant and equal to δ_c = 1.68, but rather distributed around a mean value which increases with decreasing halo mass. Analyses of N-body simulations have confirmed this prediction and showed the scatter around the mean barrier is always significant (Dalal et al. 2008; Robertson et al. 2009; Elia, Ludlow & Porciani 2012). Since the stochasticity induced by triaxial collapse is somewhat cumbersome to implement in analytic models of halo collapse (see e.g. Hahn & Paranjape 2014, for a tentative implementation with the peak constraint), we will consider a simple approximation calibrated with numerical simulations (note that it differs from the diffusing barrier approach of Maggiore & Riotto 2010). Namely the square-root stochastic barrier

\begin{equation} B= \delta _{\rm c} + \beta \sigma _0 , \end{equation}

(7)

wherein the stochastic variable β closely follows a log-normal distribution, furnishes a good description of the critical collapse threshold as a function of halo mass (Robertson et al. 2009). In PSD, this result was interpreted as follows: each halo ‘sees’ a moving barrier B = δ_c + βσ₀ with a value of β drawn from a log-normal distribution. Therefore, the first-crossing condition becomes

\begin{equation} B < \delta < B + \left(B^{\prime } + \mu \right) \Delta R , \end{equation}

(8)

where the prime designates a derivative w.r.t. the filtering scale. Here, however, we will assume that each halo ‘sees’ a constant (flat) barrier, whose height varies from halo to halo. Therefore, we will implement the first-crossing condition simply as

\begin{equation} B < \delta < B + \mu \, \Delta R . \end{equation}

(9)

Consequently, the variable μ will satisfy the constraint μ > 0 rather than μ > −B^′.

With the aforementioned modifications, the ESP multiplicity function reads

\begin{eqnarray} f_{\rm ESP}(\nu _c) &=& \left(\displaystyle\frac{V}{V_*}\right) \frac{1}{\gamma _{\nu \mu }\nu _c}\int _0^\infty \mathrm{d}\beta \, p(\beta ) \nonumber \\ && \times \int _0^{\infty }{\rm d}\mu \,\mu \int _0^{\infty } \mathrm{d}u\, f(u,\alpha =1)\, \mathcal {N}(\nu _c+\beta ,u,\mu ),\nonumber \\ \end{eqnarray}

(10)

where V is the Lagrangian volume associated with the tophat smoothing filter, V_* is the characteristic volume of peaks, p(β) is a log-normal distribution for which we take 〈β〉 = 0.5 and Var(β) = 0.25 as in PSD and f(u, α) is the slightly modified form (see Desjacques et al. 2010) of the original curvature function of Bardeen et al. (1986); see Appendix A. We can now apply Bayes’ theorem and write |$\mathcal {N}(\nu ,u,\mu )=\mathcal {N}(\nu ,u)\mathcal {N}(\mu |\nu ,u)$|⁠. The integral over μ,

\begin{equation} \int _0^\infty \!\!\mathrm{d}\mu \,\mu \, {\cal N}(\mu |\nu ,u), \end{equation}

(11)

is the same as in Musso & Sheth (2012) and, therefore, is equal to

\begin{equation} \bar{\mu } \left[\frac{1+{\rm erf}(\bar{\mu }/\sqrt{2}\Sigma )}{2} +\frac{\Sigma }{\sqrt{2\pi }\bar{\mu }}{\rm e}^{-\bar{\mu }^2/2\Sigma ^2}\right], \end{equation}

(12)

where

\begin{equation} \bar{\mu } = u\left(\frac{\gamma _{u\mu }-\gamma _1\gamma _{\nu \mu }}{1-\gamma _1^2}\right) +(\nu +\beta )\left(\frac{\gamma _{\nu \mu }-\gamma _1\gamma _{u\mu }}{1-\gamma _1^2}\right) \end{equation}

(13)

\begin{equation} \Sigma ^2 = \Delta _0^2-\frac{\gamma _{\nu \mu }^2-2\gamma _1\gamma _{\nu \mu }\gamma _{u\mu } +\gamma _{u\mu }^2}{1-\gamma _1^2}. \end{equation}

(14)

Substituting this expression into equation (10) and performing numerically the integrals over u and β, we obtain an analytic prediction for the halo mass function without any free parameter. Our ESP mass function differs at most by 2–3 per cent over the mass range 10¹¹–10¹⁵ M_⊙ h⁻¹ from that obtained with the prescription of PSD. Likewise, the linear and quadratic local bias parameters are barely affected by our modifications.

2.3 Comparison with numerical simulations

To test the validity of our approach, we compare the ESP mass function with that of haloes extracted from N-body simulations. For this purpose, we ran a series of N-body simulations evolving 1024³ particles in periodic cubic boxes of size 1500 and 250 h⁻¹ Mpc. The particle mass thus is 2.37 × 10¹¹ and 1.10 × 10⁹ M_⊙ h⁻¹, respectively. The transfer function was computed with CAMB (Lewis, Challinor & Lasenby 2000) assuming parameter values consistent with those inferred by WMAP7 (Komatsu et al. 2011): a flat ΛCDM cosmology with h = 0.704, Ω_m = 0.272, Ω_b = 0.0455, n_s = 0.967 and a normalization amplitude σ₈ = 0.81. Initial conditions were laid down at redshift z = 99 with an initial particle displacement computed at second order in Lagrangian perturbation theory with 2LPTic (Crocce, Pueblas & Scoccimarro 2006). The simulations were run using the N-body code GADGET-2 (Springel 2005) while the haloes were identified with the spherical overdensity halo finder AHF (Knollmann & Knebe 2009) assuming an overdensity threshold Δ_c = 200 constant throughout redshift.

In Fig. 1, we compare the simulated halo mass function to the ESP prediction at redshift z = 0 and 1. The latter can be straightforwardly obtained from the multiplicity function f_ESP(ν_c) as

\begin{eqnarray} \frac{{\rm d} \bar{n}_{\rm h}}{{\rm d}\,{\rm ln}M} &=& \displaystyle\frac{\bar{\rho }}{M} \nu _c f_{\rm ESP}(\nu _c,R_{\rm s}) \frac{{\rm d}\log \nu _c}{{\rm d}\log M} \nonumber \\ &=& -3 R_{\rm T}\left(\frac{\gamma _{\nu \mu }\nu _c}{\sigma _{0T}}\right) V^{-1}f_{\rm ESP}(\nu _c), \end{eqnarray}

(15)

where we used the fact that |$\gamma _{\nu \mu }= \sigma _{0T}^{\prime }$| to obtain the second equality. The ESP prediction agrees with the simulations at the 10 per cent level or better from 10¹⁴ M_⊙ h⁻¹ down to a halo mass 10¹¹ M_⊙ h⁻¹, where the correspondence between virialized haloes and initial density peaks should be rather vague. The abundance of very rare clusters with M > 10¹⁴ M_⊙ h⁻¹ is difficult to predict because of exponential sensitivity to δ_c. In this respect, it might be more appropriate to work with a critical linear density δ_c ≈ 1.60 if haloes are defined with a fixed non-linear threshold Δ_c = 200 relative to the mean density (see, e.g. Barkana 2004; Valageas 2009, for a discussion).

$Halo mass function measured from N-body simulation at redshift z = 0 (left-hand panel) and z = 1 (right-hand panel) with different box sizes as indicated in the figures. The error bars are Poisson. The data are compared to the theoretical prediction (equation 15) based on the ESP formalism and the fitting formula of Tinker et al. (2008). We also show the fractional deviation of the Tinker et al. (2008) and the measured halo mass function relative to our theoretical prediction.$

Figure 1.

Halo mass function measured from N-body simulation at redshift z = 0 (left-hand panel) and z = 1 (right-hand panel) with different box sizes as indicated in the figures. The error bars are Poisson. The data are compared to the theoretical prediction (equation 15) based on the ESP formalism and the fitting formula of Tinker et al. (2008). We also show the fractional deviation of the Tinker et al. (2008) and the measured halo mass function relative to our theoretical prediction.

Open in new tab Download slide

2.4 Bias parameters

The bias factors of ESP can be computed using the same formulas as in Desjacques (2013). With the additional variable μ, the ‘localized’ number density (in the terminology of Matsubara 2012) can be written as (Desjacques et al. 2013)

\begin{equation} n_{\rm ESP}({\boldsymbol w})=-\left(\frac{\mu }{\gamma _{\nu \mu }\nu _c}\right)\theta _H(\mu )\,n_{\rm pk}(\mathrm{\boldsymbol y}), \end{equation}

(16)

where n_pk is the localized number density of BBKS peaks and |${\boldsymbol w}=(\nu ,\eta _i,\zeta _{ij},\mu )\equiv ({\boldsymbol y}, \mu )$| is an 11-dimensional vector containing all the independent variables of the problem. Therefore,

\begin{eqnarray} \sigma _{0T}^i\sigma _{2{\rm G}}^j b_{ijk}&=& \displaystyle\frac{1}{\bar{n}_{\rm ESP}}\int \!\! {\rm d}^{11}\mathrm{\boldsymbol w}\, n_{\rm ESP}(\mathrm{\boldsymbol w}) H_{ijk}(\nu ,u,\mu ) P_1(\mathrm{\boldsymbol w}) \nonumber \\ \sigma _{1{\rm G}}^{2k}\chi _{k0}&=& \frac{(-1)^k}{\bar{n}_{\rm ESP}}\int \!\! {\rm d}^{11}\mathrm{\boldsymbol w}\, n_{\rm ESP}(\mathrm{\boldsymbol w}) L^{(1/2)}_k\!\!\left(\frac{3\eta ^2}{2}\right)P_1(\mathrm{\boldsymbol w}) \nonumber \\ \sigma _{2{\rm G}}^{2k}\chi _{0k}&=& \frac{(-1)^k}{\bar{n}_{\rm ESP}}\int \!\! {\rm d}^{11}{\boldsymbol w}\, n_{\rm ESP}(\mathrm{\boldsymbol w}) L^{(3/2)}_k\!\!\left(\frac{5\zeta ^2}{2}\right)P_1(\mathrm{\boldsymbol w}). \end{eqnarray}

(17)

Here, |$P_1(\mathrm{\boldsymbol w})$| is the one-point probability density

\begin{eqnarray} P_1(\mathrm{\boldsymbol w}){\rm d}^{11}\mathrm{\boldsymbol w}&= \mathcal {N}(\nu ,u,\mu )\, {\rm d}\nu \,{\rm d}u\,{\rm d}\mu \times \chi ^2_3(3\eta ^2)\,\mathrm{d}(3\eta ^2) \nonumber \\ & \qquad \times \chi ^2_5(5\zeta ^2)\,{\rm d}(5\zeta ^2)\times P({\rm angles}), \end{eqnarray}

(18)

where H_ijk(ν, u, μ) are trivariate Hermite polynomials and |$\chi _k^2(x)$| is a χ²-distribution with k degrees of freedom (d.o.f.). The probability density P (angles; which was missing² in Desjacques 2013) represents the probability distribution of the five remaining d.o.f. Since they are all angular variables, they do not generate bias factors because the peak (and halo) overabundance can only depend on scalar quantities (e.g. Catelan, Matarrese & Porciani 1998; McDonald & Roy 2009).

The behaviour of the bias factors b_ij0 and χ_kl as a function of halo mass is similar to that seen in fig. 1 of Desjacques (2013). The bias factors b_ijk with k ≥ 1 weight the contributions of μ^k terms to the clustering of ESP that are proportional to derivatives of the tophat filter w.r.t. the filtering scale R_T. Similar contributions appear in the clustering of thresholded regions (Matsubara 2012; Ferraro et al. 2013) since their definition also involve a first-crossing condition.

The effective bias expansion takes the form (Desjacques 2013; Desjacques et al. 2013)

\begin{eqnarray} \delta _{\rm pk}(\mathrm{\boldsymbol x}) &=& \sigma _{0T} b_{100}\nu (\mathrm{\boldsymbol x})+\sigma _{2{\rm G}} b_{010} u(\mathrm{\boldsymbol x}) + b_{001}\mu (\mathrm{\boldsymbol x}) \nonumber \\ &&\quad +\; \displaystyle\frac{1}{2}\sigma _{0T}^2 b_{200}\nu ^2(\mathrm{\boldsymbol x}) +\sigma _{0T}\sigma _{2{\rm G}} b_{110}\nu (\mathrm{\boldsymbol x}) u(\mathrm{\boldsymbol x}) \nonumber \\ &&\quad +\;\frac{1}{2}\sigma _{2{\rm G}}^2 b_{020} u^2(\mathrm{\boldsymbol x}) + \frac{1}{2}b_{002}\mu ^2(\mathrm{\boldsymbol x}) \nonumber \\ &&\quad +\;\sigma _{1{\rm G}}^2\chi _{10}\eta ^2(\mathrm{\boldsymbol x}) +\sigma _{2{\rm G}}^2\chi _{01}\zeta ^2(\mathrm{\boldsymbol x}) \nonumber \\ &&\quad +\; \sigma _{0T} b_{101}\nu (\mathrm{\boldsymbol x})\mu (\mathrm{\boldsymbol x}) +\sigma _{2{\rm G}} b_{011}u(\mathrm{\boldsymbol x})\mu (\mathrm{\boldsymbol x}) + \cdots . \end{eqnarray}

(19)

Here, the rule of thumb is that one should ignore all the terms involving zero-lag moments in the computation of |$\langle \delta _{\rm pk}(\mathrm{\boldsymbol x}_1),\ldots ,\delta _{\rm pk}(\mathrm{\boldsymbol x}_N)\rangle$| in order to correctly predict the N-point correlation function, as demonstrated explicitly in Desjacques (2013). The appearance of rotationally invariant quantities is, again, only dictated by the scalar nature of the peak overabundance. The variables of interest here are

\begin{eqnarray} \eta ^2(\mathrm{\boldsymbol x}) &=& \displaystyle\frac{1}{\sigma _{1{\rm G}}^2}\left(\nabla \delta \right)^2\!(\mathrm{\boldsymbol x}) \nonumber \\ \zeta ^2(\mathrm{\boldsymbol x}) &=& \frac{3}{2\sigma _{2{\rm G}}^2} {\rm tr}\!\biggl [\Bigl (\mathrm{\partial} _i\mathrm{\partial} _j\delta -\frac{1}{3}\delta _{ij}\nabla ^2\delta \Bigr )^2\biggr ]\!(\mathrm{\boldsymbol x}), \end{eqnarray}

(20)

so that |$3\eta ^2(\mathrm{\boldsymbol x})$| and |$5\zeta ^2(\mathrm{\boldsymbol x})$| are χ²-distributed with 3 and 5 d.o.f., respectively.

3 BIASES FROM CROSS-CORRELATION: EXTENSION TO χ² VARIABLES

In this section, we will demonstrate that the bias factors χ_ij can be measured with a one-point statistics. We will test our method on density peaks of a Gaussian random field before applying it to dark matter haloes.

3.1 Bias factors b_ijk: Hermite polynomials

Musso et al. (2012) showed that the bias factors of discrete tracers (relative to the mass density δ) can be computed from one-point measurements rather than computationally more expensive n-point correlations. Their idea was implemented by PSD and Paranjape et al. (2013b) to haloes extracted from N-body simulations in order to test the predictions of the ESP formalism. Namely haloes were traced back to their ‘protohalo’ patch (since one is interested in measuring Lagrangian biases) in the initial conditions, the linear density field was smoothed on some ‘large-scale’ R_l and the quantity H_n(ν_l = δ_l/σ_0l) was computed (for n = 1, 2 only) at the location of each protohalo. The average of H_n(ν_l) over all protohaloes reads

\begin{equation} \frac{1}{N}\sum _{i=1}^N H_n(\nu _{\rm l}) = \int _{-\infty }^{+\infty }\!\!{\rm d}\nu _{\rm l}\,{\cal N}(\nu _{\rm l}) \langle 1+\delta _{\rm h} \vert \nu _{\rm l} \rangle H_n(\nu _{\rm l}), \end{equation}

(21)

where δ_h is the overdensity of protohaloes. This expression assumes that the first-crossing condition can be implemented through a constraint of the form equation (9), so that P(ν_l) is well approximated by a Gaussian (Musso et al. 2012). For the ESP considered here, this ensemble average reads

\begin{eqnarray} &&\!\!\!\!&&{\displaystyle\frac{1}{\bar{n}_{\rm ESP}}\int _{-\infty }^{+\infty }\!\!{\rm d}\nu _{\rm l}\,{\cal N}(\nu _{\rm l}) \langle n_{\rm ESP} \vert \nu _{\rm l} \rangle H_n(\nu _{\rm l})} \nonumber \\ &&\quad= \frac{1}{\bar{n}_{\rm ESP}}\int \!\!{\rm d}^{11}\mathrm{\boldsymbol w}\,n_{\rm ESP}(\mathrm{\boldsymbol w})\,\left(-\epsilon _\nu \right)^n \nonumber \\ &&\qquad \times \left(\frac{\mathrm{\partial} }{\mathrm{\partial} \nu }+\frac{\epsilon _u}{\epsilon _\nu }\frac{\mathrm{\partial} }{\mathrm{\partial} u} +\frac{\epsilon _\mu }{\epsilon _\nu }\frac{\mathrm{\partial} }{\mathrm{\partial} \mu }\right)^n P_1(\mathrm{\boldsymbol w}). \end{eqnarray}

(22)

Here, ϵ_X denotes the cross-correlation between ν_l and the variables X = (ν, u, μ) defined at the halo smoothing scale. The right-hand side reduces to a sum of nth-order bias factors b_ijk weighted by products of ϵ_ν, ϵ_u and ϵ_μ. Relations between bias factors of a given order (which arise owing to their close connection with Hermite polynomials, see e.g. Musso et al. 2012) can then be used to extract a measurement of each b_ijk.

Before we generalize this approach to the chi-squared bias factors χ_ij, we emphasize that, in this cross-correlation approach, the smoothing scale R_l can take any value as long as it is distinct from the halo smoothing scale. Paranjape et al. (2013b) chose R_l ≫ R_s in the spirit of the peak-background split but this requirement is, in fact, not necessary as long as the correlation between the two scales is taken into account. In any case, we will stick with the notation R_l for convenience.

3.2 Bias factors χ_ij: Laguerre polynomials

The approach presented above can be generalized to χ² distributions. The main difference is the appearance of Laguerre polynomials |$L_n^{(\alpha )}$|⁠. Consider for instance the χ²-quantity 3η² smoothed at the scale R_l, i.e. |$3\eta _{\rm l}^2$|⁠. In analogy with equation (21), the ensemble average of |$L_n^{(1/2)}(3\eta _{\rm l}^2)$| at the peak positions is

\begin{eqnarray} \frac{1}{N}\sum _{i=1}^N L_n^{(1/2)}\!\left(\displaystyle\frac{3\eta _{\rm l}^2}{2}\right) &=& \int _0^\infty \!\!{\rm d}(3\eta _{\rm l}^2)\, \chi _3^2(3\eta _{\rm l}^2) \nonumber \\ &&\times \left\langle 1+\delta _{\rm h} \vert 3\eta _{\rm l}^2 \right\rangle L_n^{(1/2)}\!\!\left(\frac{3\eta _{\rm l}^2}{2}\right). \end{eqnarray}

(23)

The conditional average |$ \langle 1 + \delta _{\rm h} \vert 3\eta _{\rm l}^2 \rangle$| reads

\begin{eqnarray} \left\langle 1 + \delta _{\rm h} \vert 3\eta _{\rm l}^2 \right\rangle &=& \displaystyle\frac{1}{\bar{n}_{\rm ESP}}\int \!\! {\rm d}^{11}\mathrm{\boldsymbol w}\, n_{\rm ESP}(\mathrm{\boldsymbol w}) P_1\left (\mathrm{\boldsymbol w} \vert 3\eta _{\rm l}^2\right ) \nonumber \\ &=&\frac{1}{\bar{n}_{\rm ESP}}\int \!\!{\rm d}u \,{\rm d}\nu \,{\rm d}\mu \, {\cal N}(\nu ,u,\mu ) \nonumber \\ &&\qquad \times \int \!\!{\rm d}(3\eta ^2)\,\chi _3^2(3\eta ^2|3\eta _{\rm l}^2) \int \!\!{\rm d}(5\zeta ^2)\,\chi _5^2(5\zeta ^2) \nonumber \\ &&\qquad \times \int \!\!{\rm d}({\rm angles})\,P({\rm angles})\,n_{\rm ESP}(\mathrm{\boldsymbol w}). \end{eqnarray}

(24)

We substitute this relation into equation (23) and begin with the integration over the variable |$3\eta _{\rm l}^2$|⁠.

We use the following relation (which can be inferred from equation 7.414 of Gradshteyn & Ryzhik 1994):

\begin{equation} \int _0^\infty \!\!{\rm d}x\,{\rm e}^{-x} x^{j+\alpha }L_n^{(\alpha )}\!(x) =\frac{(-1)^n}{n!}\frac{j!\,\Gamma (j+\alpha +1)}{(j-n)!}. \end{equation}

(25)

With the aid of this result and on expanding the conditional χ²-distribution |$\chi _3^2(3\eta ^2|3\eta _{\rm l}^2)$| in Laguerre polynomials (see Appendix B for details), we obtain

\begin{eqnarray} &&\!\!\!\!&&{\int _0^\infty {\rm d}(3\eta _{\rm l}^2)\,\chi _3^2(3\eta _{\rm l}^2)\, L_n^{(1/2)}\!\left(\displaystyle\frac{3\eta _{\rm l}^2}{2}\right) \, \chi _3^2(3\eta ^2|3\eta _{\rm l}^2)} \nonumber \\ &&=\frac{(-1)^n}{n!}\frac{1}{\Gamma (3/2)} \left(\frac{3\eta ^2}{2}\right)^\alpha \frac{{\rm e}^{-3\eta ^2/2(1-\epsilon ^2)}}{2\left(1-\epsilon ^2\right)^{\alpha +1}} \nonumber \\ &&\qquad \times \sum _{j=0}^\infty \frac{j!}{(j-n)!}\left(\frac{-\epsilon ^2}{1-\epsilon ^2}\right)^j L_j^{(1/2)}\!\left[\frac{3\eta ^2}{2(1-\epsilon ^2)}\right]. \end{eqnarray}

(26)

For simplicity, let us consider the cases n = 0, 1 solely. For n = 0, the sum simplifies to

\begin{eqnarray} &&\!\!\!\!&&{\sum _{j=0}^\infty \left(\frac{-\epsilon ^2}{1-\epsilon ^2}\right)^j L_j^{(1/2)}\!\left[\frac{3\eta ^2}{2(1-\epsilon ^2)}\right]}\nonumber \\ &&\quad = \left(1-\epsilon ^2\right)^{3/2} \exp \!\left[\left(\frac{\epsilon ^2}{1-\epsilon ^2}\right)\frac{3\eta ^2}{2}\right], \end{eqnarray}

(27)

and the integral equation (26) (⁠|$L_0^{(1/2)}(3\eta _{\rm l}^2/2)\equiv 1$|⁠) is trivially equal to |$\chi _3^2(3\eta ^2)$| (as it should be, since we are essentially marginalizing over |$3\eta _{\rm l}^2$|⁠).

For n ≥ 1, the sum can be evaluated upon taking suitable derivatives of the right-hand side of equation (27), which indeed is a generating function for the Laguerre polynomials |$L_n^{(1/2)}$|⁠. For n = 1, a little algebra leads to

\begin{eqnarray} &&\!\!\!\!&&{\sum_{j=0}^\infty j \left(\frac{-\epsilon^2}{1-\epsilon^2}\right)^{j-1} L_j^{(1/2)}\!\left[\frac{3\eta^2}{2(1-\epsilon^2)}\right]}\nonumber \\ &&\quad= \left(1-\epsilon^2\right)^{5/2} L_1^{(1/2)}\!\left(\frac{3\eta^2}{2}\right) \exp\left[\left(\frac{\epsilon^2}{1-\epsilon^2}\right)\frac{3\eta^2}{2}\right]. \end{eqnarray}

(28)

Hence, equation (26) with n = 1 equals |$\epsilon ^2 L_1^{(1/2)}(3\eta ^2/2)\chi _3^2(3\eta ^2)$|⁠. Performing the remaining integrals over ν, u, μ and 5ζ² (the integral over the angles is trivially unity) and taking into account the ESP constraint through the multiplicative factor |$n_{\rm ESP}(\mathrm{\boldsymbol w})$|⁠, equation (23) simplifies to

\begin{eqnarray} \int _0^\infty \!\!{\rm d}(3\eta _{\rm l}^2)\, \chi _3^2(3\eta _{\rm l}^2)\, \left\langle 1+\delta _{\rm h} \vert 3\eta _{\rm l}^2 \right\rangle \, L_1^{(1/2)}\!\!\left(\frac{3\eta _{\rm l}^2}{2}\right) = - \epsilon ^2 \sigma _1^2 \chi _{10}.\nonumber \\ \end{eqnarray}

(29)

For the variable 3η², the cross-correlation coefficient ϵ is

\begin{equation} \epsilon ^2 \equiv \frac{\langle \eta ^2\eta _{\rm l}^2 \rangle - \langle \eta ^2 \rangle \langle \eta _{\rm l}^2 \rangle }{\sqrt{ \left( \langle \eta ^4 \rangle - \langle \eta ^2 \rangle ^2\right) \left( \langle \eta _{\rm l}^4 \rangle - \langle \eta _{\rm l}^2 \rangle ^2\right) }} =\left(\frac{\sigma _{1X}^2}{\sigma _{1{\rm s}}\sigma _{1{\rm l}}}\right)^2, \end{equation}

(30)

which we shall denote as ϵ₁ in what follows. Furthermore,

\begin{eqnarray} \sigma _{n\times }^2 = \frac{1}{2\pi ^2} \int _0^\infty \!\!{\rm d}k\,k^{2(n+1)}\,P(k) W_{\rm G}(k R_{\rm s}) W_{\rm G}(k R_{\rm l}) \end{eqnarray}

(31)

designates the splitting of filtering scales, i.e. one filter is on scale R_s while the second is on scale R_l. It should be noted that, unlike σ_1X defined in equation (6), both filtering kernels are Gaussian.

The derivation of the bias factors χ_0k associated with the quadratic variable ζ² proceeds analogously. In particular,

\begin{eqnarray} \int _0^\infty \!\!{\rm d}(5\zeta _{\rm l}^2)\, \chi _5^2(5\zeta _{\rm l}^2)\, \left \langle 1+\delta _{\rm h} \vert 5\zeta _{\rm l}^2 \right\rangle \, L_1^{(3/2)}\!\!\left(\frac{5\zeta _{\rm l}^2}{2}\right) = - \epsilon ^2 \sigma _2^2 \chi _{01}.\nonumber \\ \end{eqnarray}

(32)

Here, the cross-correlation coefficient is |$\epsilon =\sigma _{2X}^2/(\sigma _{2{\rm s}}\sigma _{2{\rm l}})\equiv \epsilon _2$|⁠. Note that, in both cases, the cross-correlation coefficient drops very rapidly as R_l moves away from R_s for realistic CDM power spectra. In addition, one could in principle choose R_l < R_s (if there is enough numerical resolution) to measure χ_ij.

4 TEST WITH NUMERICAL SIMULATIONS

In this section, we first validate our predictions based on peaks of Gaussian random fields with measurements extracted from random realizations of the Gaussian linear density field, and then move on to calculate χ₁₀ and χ₀₁ for M ≳ M_⋆ haloes, where M_⋆ is the characteristic mass of the haloes.

4.1 Peaks of Gaussian random fields

We generate random realizations of the Gaussian, linear density field with a power spectrum equal to that used to seed the N-body simulations described above. To take advantage of Fast Fourier Transforms (FFTs), we simulate the linear density field in periodic, cubic boxes of side 1000 h⁻¹ Mpc. The size of the mesh along each dimension is 1536. We smooth the density field on scale R_s = 5 h^{− 1} Mpc with a tophat filter and find the local maxima by comparing the density at each grid point with its 26 neighbouring values.

We then smooth the density field on the larger scales R_l = 10, 15 and 20 h⁻¹ Mpc with a Gaussian filter and compute

\begin{equation} \eta _{\rm l}^2 = \frac{1}{\sigma _{1{\rm l}}^2}(\nabla \delta _{\rm l})^2 \end{equation}

(33)

\begin{equation} \zeta _{\rm l}^2 = \frac{3}{2\sigma _{2{\rm l}}^2} {\rm tr}\!\biggl [\Bigl (\mathrm{\partial} _i\mathrm{\partial} _j\delta _{\rm l} -\frac{1}{3}\delta _{ij}\nabla ^2\delta _{\rm l}\Bigr )^2\biggr ]. \end{equation}

(34)

These density fields with derivatives sensitively depend on the smoothing scales used. To illustrate this, we show in Fig. 2 sections of ν_l, |$3 \eta _{\rm l}^2$| and |$5 \zeta _{\rm l}^2$|⁠. The sections, each of which of dimensions 200 × 200 h⁻² Mpc², were generated at z = 99 with the same random seed. The first row corresponds to R_s = 5 h^{− 1} Mpc, whereas the second row displays results on the filtering scale R_l = 10 h^{− 1} Mpc. We note that, for the normalized density field ν_l, an increase in the smoothing scale washes out the small-scale features, but the large-scale pattern remains. For the quadratic variable |$\eta _{\rm l}^2$| however, the resemblance between the features at the small and large filtering scale is tenuous. This is even worse for |$\zeta _{\rm l}^2$|⁠.

$Sections for νl, $3 \eta _{\rm l}^2$ and $5 \zeta _{\rm l}^2$ (from left to right). A filtering scale of Rl = 5 and 10 h−1 Mpc is used for the first and second row, respectively. Note that a tophat kernel is applied for νl, while a Gaussian window is used for $\eta _{\rm l}^2$ and $\zeta _{\rm l}^2$. In each panel, the dimension of the section is 200 × 200 h−2 Mpc2.$

Figure 2.

Sections for ν_l, |$3 \eta _{\rm l}^2$| and |$5 \zeta _{\rm l}^2$| (from left to right). A filtering scale of R_l = 5 and 10 h⁻¹ Mpc is used for the first and second row, respectively. Note that a tophat kernel is applied for ν_l, while a Gaussian window is used for |$\eta _{\rm l}^2$| and |$\zeta _{\rm l}^2$|⁠. In each panel, the dimension of the section is 200 × 200 h⁻² Mpc².

Open in new tab Download slide

Compared to ν_l, the fields |$\eta _{\rm l}^2$| and |$\zeta _{\rm l}^2$| have one and two additional derivatives which give rise to an effective window function whose isotropic part is given by

\begin{equation} W_{\rm eff}(k,R) = k^n {\rm e}^{ -( k R)^2 /2 }, \end{equation}

(35)

where n = 0, 1 and 2 for ν_l, |$\eta _{\rm l}^2$| and |$\zeta _{\rm l}^2$|⁠, respectively. For n = 0, the window becomes narrower as R_l increases, yet remains unity for wavenumbers k ≲ 1/R_l. W_eff reaches a maximum at |$\sqrt{n} /R$|⁠. Hence, for n = 1 and 2, W_eff selects predominantly wavemodes with k ∼ 1/R. Consequently, since in a Gaussian random field the wavemodes at different scales are uncorrelated, patterns in the fields |$\eta _{\rm l}^2$| and |$\zeta _{\rm l}^2$| can change drastically as R_l varies. This effect is expected to be most significant for n = 2, i.e. |$\zeta _{\rm l}^2$|⁠.

For each local density maxima, we store the peak height ν as well as the value of |$\eta _{\rm l}^2$| and |$\zeta _{\rm l}^2$| at the peak position. The left-hand panel of Fig. 3 displays as histograms the resulting probability distribution |$P(3\eta _{\rm l}^2|{\rm pk})$| for three different values of R_l = 10, 15 and 20 h⁻¹ Mpc. The solid curves represent the theoretical prediction equation (B4) with x = 〈3η²|pk〉 = 0 and ϵ₁ = 0.71, 0.44 and 0.29 (from the smallest to largest R_l) as was measured from the random realizations. The dashed curve is the unconditional χ²-distribution with 3 d.o.f. The theory gives excellent agreement with the simulations. Note also that we did not find any evidence for a dependence on the peak height, as expected from the absence of a correlation between ν and |$\eta _{\rm l}^2$|⁠. The right-hand panel of Fig. 3 shows results for |$\zeta _{\rm l}^2$|⁠. Here, however, since the cross-correlation coefficient diminishes very quickly when R_l even slightly departs from R_s, we show result for R_l = 10 h^{− 1} Mpc only, which corresponds to ϵ₂ = 0.57. In addition, because one should expect a dependence of the shape of the density profile around peaks to the peak height, we consider three different ranges of ν as indicated in the figure. The solid curves indicate the theoretical prediction equation (B4) with ϵ₂ = 0.57 and x = 〈5ζ²|pk〉, where

\begin{equation} \langle 5\zeta ^2|{\rm pk}\rangle = -2\mathrm{\partial} _\alpha {\rm ln}\int _{\nu _{\rm min}}^{\nu _{\rm max}}{\rm d}\nu \, G_0^{(\alpha )}\!(\gamma _1,\gamma _1\nu ). \end{equation}

(36)

Here, |$G_0^{(\alpha )}$| is the integral of f(u, α) over all the allowed peak curvatures. The average 〈5ζ²|pk〉 increases with the peak height to reach 5 in the limit ν → ∞. The figure shows a clear deviation from the unconditional distribution |$\chi _5^2(5\zeta _{\rm l}^2)$| (shown as the dashed curve) and a dependence on ν consistent with theoretical predictions.

$Conditional probability distribution for the variables $3\eta _{\rm l}^2$ (left-hand panel) and $5\zeta _{\rm l}^2$ (right-hand panel) measured at the position of maxima of the linear density field smoothed with a Gaussian filter on scale R = 5 h−1 Mpc. Left-hand panel: histograms indicate the results for Rl = 10, 15 and 20 h−1 Mpc, which leads to ϵ1 = 0.71, 0.44 and 0.29 as quoted in the figure. Right-hand panel: histograms show the results for a fixed Rl = 10 h− 1 Mpc (which implies ϵ2 = 0.57) but several peak height intervals. In all cases, the solid curves are the theoretical prediction (see text) whereas the dashed (green) curves represents the unconditional distribution $\chi _k^2(y)$.$

Figure 3.

Conditional probability distribution for the variables |$3\eta _{\rm l}^2$| (left-hand panel) and |$5\zeta _{\rm l}^2$| (right-hand panel) measured at the position of maxima of the linear density field smoothed with a Gaussian filter on scale R = 5 h⁻¹ Mpc. Left-hand panel: histograms indicate the results for R_l = 10, 15 and 20 h⁻¹ Mpc, which leads to ϵ₁ = 0.71, 0.44 and 0.29 as quoted in the figure. Right-hand panel: histograms show the results for a fixed R_l = 10 h^{− 1} Mpc (which implies ϵ₂ = 0.57) but several peak height intervals. In all cases, the solid curves are the theoretical prediction (see text) whereas the dashed (green) curves represents the unconditional distribution |$\chi _k^2(y)$|⁠.

Open in new tab Download slide

4.2 Dark matter haloes

Having successfully tested the theory against numerical simulations of Gaussian peaks, we will now attempt to estimate the bias factors χ₁₀ and χ₀₁ associated with dark matter haloes. For this purpose, we first trace back all dark matter particles belonging to virialized haloes at redshift z = 0 to their initial position at z = 99. We then compute the centre-of-mass positions of these Lagrangian regions and assume that they define the locations of protohaloes. We can now proceed as for the Gaussian peaks and compute ν, |$\eta _{\rm l}^2$| and |$\zeta _{\rm l}^2$| at the position of protohaloes.

The quadratic bias factors χ₁₀ and χ₀₁ could be in principle computed analogously to Paranjape et al. (2013b), i.e. by stacking measurements of |$\eta _{\rm l}^2$| and |$\zeta _{\rm l}^2$| at the locations of protohaloes as

\begin{equation} \sigma _{1{\rm s}}^2 \hat{\chi }_{10} = -\frac{1}{N\epsilon _1^2}\sum _{i=1}^N L_1^{(1/2)}\!\left(\frac{3\eta _{\rm l}^2}{2}\right) \end{equation}

(37)

and

\begin{equation} \sigma _{2{\rm s}}^2 \hat{\chi }_{01} = -\frac{1}{N\epsilon _2^2}\sum _{i=1}^N L_1^{(3/2)}\!\left(\frac{5\zeta _{\rm l}^2}{2}\right). \end{equation}

(38)

Here, N is the number of haloes, s designates smoothing at the halo mass scale with a Gaussian filter W_G on scale R_G(R_T), whereas l designates Gaussian smoothing at the large-scale R_l. However, because the cross-correlation coefficient is fairly small unless R_l is very close to R_G, we decided to compute χ₁₀ and χ₀₁ by fitting the probability distribution |$P(3\eta _{\rm l}^2|{\rm halo})$| and |$P(5\zeta _{\rm l}^2|{\rm halo})$| with the conditional χ²-distribution |$\chi _k^2(y|x)$|⁠. Namely

\begin{eqnarray} \sigma _{1{\rm s}}^2\hat{\chi }_{10} &=& \displaystyle\frac{1}{2}\left(\langle 3\eta ^2|{\rm halo}\rangle -3\right) \nonumber \\ \sigma _{2{\rm s}}^2\hat{\chi }_{01} &=& \frac{1}{2}\left(\langle 5\zeta ^2|{\rm halo}\rangle -5\right), \end{eqnarray}

(39)

where 〈3η²|halo〉 and 〈5ζ²|halo〉 are the best-fitting values obtained for x. We used measurements obtained at the smoothing scale R_l = 10 h^{− 1} Mpc only to maximize the signal.

To predict the value of R_G given R_T, we followed PSD and assumed that R_G(R_T) can be computed through the requirement that 〈δ_G|δ_T〉 = δ_T. This yields a prediction for the value of the cross-correlation coefficients ϵ₁ and ϵ₂ as a function of halo mass, which we can use as an input to |$\chi _k^2(y|x)$| and only fit for x. However, we found that using the predicted ϵ₁ leads to unphysical (negative) values for x when one attempts to fit |$P(3\eta _{\rm l}^2|{\rm halo})$|⁠. Therefore, we decided to proceed as follows.

Estimate both ϵ₁ and x = 〈3η²|halo〉 by fitting the model |$\chi _3^2(y|x;\epsilon _1)$| to the measured |$P(3\eta _{\rm l}^2|{\rm halo})$|⁠.
Compute ϵ₂ assuming that the same R_G enters the spectral moments.
Estimate |$x=\langle 5\zeta _{\rm l}^2|{\rm halo}\rangle$| by fitting the theoretical model |$\chi _5^2(y|x;\epsilon _2)$| to the simulated |$P(5\zeta _{\rm l}^2|{\rm halo})$|⁠.

We considered data in the range |$0<3\eta _{\rm l}^2<8$| and |$0<5\zeta _{\rm l}^2<12$| and gave equal weight to all the measurements (assuming Poisson errors does not affect our results significantly). Table 1 summarizes the best-fitting values obtained for four different halo bins spanning the mass range 10¹³-10¹⁵ M_⊙ h^{− 1}, whereas the measured probability distributions together with the best-fitting models are shown in Fig. 4. The data are reasonably well described by a conditional χ²-distribution, but the fit is somewhat poorer when the cross-correlation coefficient is close to unity.

$Conditional probability distribution for $3\eta _{\rm l}^2$ (top panel) and $5\zeta _{\rm l}^2$ (bottom panel) measured at the centre-of-mass position of protohaloes. The filter is Gaussian with Rl = 10 h− 1 Mpc. The various curves show the best-fitting theoretical predictions for the halo mass bins considered here. Halo mass range is in unit of 1013 M⊙ h−1. Poisson errors are much smaller than the size of the data points and, therefore, do not show up in the figure.$

Figure 4.

Conditional probability distribution for |$3\eta _{\rm l}^2$| (top panel) and |$5\zeta _{\rm l}^2$| (bottom panel) measured at the centre-of-mass position of protohaloes. The filter is Gaussian with R_l = 10 h^{− 1} Mpc. The various curves show the best-fitting theoretical predictions for the halo mass bins considered here. Halo mass range is in unit of 10¹³ M_⊙ h⁻¹. Poisson errors are much smaller than the size of the data points and, therefore, do not show up in the figure.

Open in new tab Download slide

Table 1.

Open in new tab

Best-fitting parameter values as a function of halo mass. The latter is in unit of 10¹³ M_⊙ h⁻¹. Note that we also list the values of ϵ₂ even though it is not directly fitted to the data (see text for details).

Halo mass	〈3η²\|halo〉	ϵ₁	〈5ζ²\|halo〉	(ϵ₂)
M > 30	0.71	0.80	2.98	(0.70)
10 < M < 30	1.24	0.66	4.49	(0.52)
3 < M < 10	1.62	0.54	5.82	(0.37)
1 < M < 3	1.94	0.49	6.12	(0.31)

Halo mass	〈3η²\|halo〉	ϵ₁	〈5ζ²\|halo〉	(ϵ₂)
M > 30	0.71	0.80	2.98	(0.70)
10 < M < 30	1.24	0.66	4.49	(0.52)
3 < M < 10	1.62	0.54	5.82	(0.37)
1 < M < 3	1.94	0.49	6.12	(0.31)

Table 1.

Open in new tab

Halo mass	〈3η²\|halo〉	ϵ₁	〈5ζ²\|halo〉	(ϵ₂)
M > 30	0.71	0.80	2.98	(0.70)
10 < M < 30	1.24	0.66	4.49	(0.52)
3 < M < 10	1.62	0.54	5.82	(0.37)
1 < M < 3	1.94	0.49	6.12	(0.31)

Halo mass	〈3η²\|halo〉	ϵ₁	〈5ζ²\|halo〉	(ϵ₂)
M > 30	0.71	0.80	2.98	(0.70)
10 < M < 30	1.24	0.66	4.49	(0.52)
3 < M < 10	1.62	0.54	5.82	(0.37)
1 < M < 3	1.94	0.49	6.12	(0.31)

The second-order bias factors χ₁₀ and χ₀₁ of the dark matter haloes at z = 0 can be readily computed from equation (39) using the best-fitting values of 〈3η²|halo〉 and 〈5ζ²|halo〉. The results are shown in Fig. 5 as the data points. Error bars indicate the scatter among the various realizations and, therefore, likely strongly underestimate the true uncertainty. The dashed curves indicate the predictions of the ESP formalism. The measurements, albeit of the same magnitude as the theoretical predictions, quite disagree with expectations based on our ESP approach, especially χ₀₁ which reverses sign as the halo mass drops below 10¹⁴ M_⊙ h⁻¹.

$The bias factors $\sigma _1^2\chi _{10}$ and $\sigma _2^2\chi _{01}$ of dark matter haloes identified in the N-body simulations at z = 0 are shown as filled (green) circle and (blue) triangle, respectively. Error bars indicate the scatter among six realizations. The horizontal dashed (green) line at −3/2 and the dashed (blue) curve are the corresponding ESP predictions. The dotted (blue) curve is $\sigma _2^2\chi _{01}$ in a model where haloes are allowed to collapse in filamentary-like structures. The solid curves are our final predictions, which take into account the offset between peak position and protohalo centre-of-mass (see text for details).$

Figure 5.

The bias factors |$\sigma _1^2\chi _{10}$| and |$\sigma _2^2\chi _{01}$| of dark matter haloes identified in the N-body simulations at z = 0 are shown as filled (green) circle and (blue) triangle, respectively. Error bars indicate the scatter among six realizations. The horizontal dashed (green) line at −3/2 and the dashed (blue) curve are the corresponding ESP predictions. The dotted (blue) curve is |$\sigma _2^2\chi _{01}$| in a model where haloes are allowed to collapse in filamentary-like structures. The solid curves are our final predictions, which take into account the offset between peak position and protohalo centre-of-mass (see text for details).

Open in new tab Download slide

4.3 Interpretation of the measurements

To begin with, we note that, if haloes were forming out of randomly distributed patches in the initial conditions, then both χ₁₀ and χ₀₁ would be zero since 〈3η²〉 = 3 and 〈5ζ²〉 = 5 for random field points.

The measured dimensionless bias factor |$\sigma _1^2\chi _{10}$| is always negative, which indicates that haloes collapse out of regions which have values of η² smaller than average. In our ESP approach, we assume that the centre-of-mass position of protohaloes exactly coincides with that of a local density peak, so that |$\sigma _1^2\chi _{10}\equiv -3/2$|⁠. However, simulations indicate that, while there is a good correspondence between protohaloes and linear density peaks, the centre-of-mass of the former is somewhat offset relative to the peak position (see e.g. Porciani, Dekel & Hoffman 2002; Ludlow & Porciani 2011). To model this effect, we note that, if the protohalo is at a distance R from a peak, then the average value of 3η² is |$ \langle 3\eta ^2 \rangle (R) = \epsilon _1^2(R) ( \langle 3\eta ^2|{\rm pk} \rangle -3)$| (in analogy with the fact that the average density at a distance R from a position where δ ≡ δ_c is 〈δ〉(R) = ξ_δ(R) δ_c). Assuming that the offset R follows a Gaussian distribution, the halo bias factor is

\begin{equation} \sigma _{1{\rm s}}^2\chi _{10} = -\frac{3}{2}\sqrt{\frac{2}{\pi }}\int _0^\infty \! \frac{{\rm d}R}{\sigma }\left(\frac{R}{\sigma }\right)^2 {\rm e}^{-R^2/2\sigma ^2}\epsilon _1^2(R). \end{equation}

(40)

The rms variance σ(M) of the offset distribution, which generally depends on the halo mass, can be constrained from our measurements of χ₁₀ for dark matter haloes. The best-fitting power-law function,

\begin{equation} \sigma (M) = 2.50 \left(\frac{M}{10^{13}\,{\rm M}_{\odot }\,h^{-1}}\right)^{0.063}\,{\rm {\it h}^{-1}\,Mpc}, \end{equation}

(41)

turns out to be a weak function of halo mass. In unit of the (tophat) Lagrangian halo radius, this translates into σ/R_T ≈ 0.79 and ≈0.36 for a halo mass M = 10¹³ and 10¹⁴ M_⊙ h⁻¹, respectively. The resulting theoretical prediction is shown as the solid curve in Fig. 5 and agrees reasonably well with our data. This crude approximation demonstrates that an offset between the protohalo centre-of-mass and the peak position can have a large impact on the inferred value of χ₁₀, since the latter is very sensitive to small-scale mass distribution.

Likewise, an offset between the protohalo centre-of-mass and the position of the linear density peak will also impact the measurement of χ₀₁, yet cannot explain the observed sign reversal. In this regard, one should first remember that density peaks become increasingly spherical as ν → ∞. Nevertheless, while their mean ellipticity 〈e〉 and prolateness 〈p〉 converge towards zero in this limit, 〈v〉 = 〈ue〉 approaches 1/5 at fixed u (see equation 7.7 of Bardeen et al. 1986). Hence, 〈ζ²〉 = 〈3v² + w²〉 does not tend towards zero but rather unity, like for random field points. Consequently, |$\sigma _2^2\chi _{01}\rightarrow 0$| in the limit ν → ∞. Secondly, at any finite ν, our ESP approach predicts that χ₀₁ be negative because we have assumed that protohaloes only form near a density peak (λ₃ > 0, where λ₁ ≥ λ₂ ≥ λ₃ are the eigenvalues of − ∂_i∂_jδ). However, N-body simulations strongly suggest that a fraction of the protohaloes collapse along the ridges or filaments connecting two density maxima, and that this fraction increases with decreasing halo mass (Ludlow & Porciani 2011). To qualitatively assess the impact of such primeval configurations on χ₀₁, we extend the integration domain in the plane (v, w) to include all the points with λ₂ > 0 and λ₃ < 0 (but still require that the curvature u be positive). This way we not only consider density peaks, but also extrema that correspond to filamentary configurations. The resulting curvature function f(u, α) can be cast into the compact form

\begin{eqnarray} f(u,\alpha ) &=& \displaystyle\frac{1}{\alpha ^4} \Biggl \lbrace \frac{{\rm e}^{-{5\alpha u^2\over 2}}}{\sqrt{10\pi }} \left(\alpha u^2-\frac{16}{5}\right) \nonumber \\ &&+\,\frac{{\rm e}^{-{5\alpha u^2\over 8}}}{\sqrt{10\pi }}\left(31\alpha u^2+\frac{32}{5}\right) +\frac{\sqrt{\alpha }}{2}\left(\alpha u^3-3u\right) \nonumber \\ && \times \,\left[{\rm erf}\left(\sqrt{\frac{5\alpha }{2}}\frac{u}{2}\right)+ {\rm erf}\left(\sqrt{\frac{5\alpha }{2}}u\right)-1\right]\Biggr \rbrace. \end{eqnarray}

(42)

The dotted curve in Fig. 5 shows |$\sigma _2^2\chi _{01}$| when the filamentary configurations are included. While it agrees with the original ESP prediction at large halo mass, it reverses sign around 10¹⁴ M_⊙ h⁻¹ because, as the peak height decreases, configurations with λ₃ < 0 or, equivalently, large values of ζ² become more probable. The solid curve takes into account, in addition to filamentary configurations, an offset between the protohalo and the peak position according to the simple prescription discussed above. This is our final prediction for |$\sigma _2^2\chi _{01}$|⁠. It is clearly at odds with the measurements, which strongly suggest that |$\sigma _2^2\chi _{01}$| can be very different from zero for M ≳ 10¹³ M_⊙ h⁻¹.

It is beyond the scope of this paper to work out a detailed description of the measurements. Using a value of R_G different than that obtained through the condition 〈δ_G|δ_T〉 = δ_T has a large impact on the mass function, suggesting that it will be difficult to get a good fit of both the mass function and the bias factors χ₁₀ and χ₀₁. Before concluding however, we note that, if the Lagrangian clustering of haloes also depends on |$s_2(\mathrm{\boldsymbol x})=s_{ij}(\mathrm{\boldsymbol x})s^{ij}(\mathrm{\boldsymbol x})$|⁠, where (in suitable units)

\begin{equation} s_{ij}(\mathrm{\boldsymbol x}) = \mathrm{\partial} _i\mathrm{\partial} _j\phi (\mathrm{\boldsymbol x}) - \frac{1}{3}\delta _{ij}\delta (\mathrm{\boldsymbol x}), \end{equation}

(43)

then we are not measuring χ₀₁ but some weighted and scale-dependent combination of both χ₀₁ and the Lagrangian bias γ₂ associated with |$s_2(\mathrm{\boldsymbol x})$|⁠. Recent numerical work indeed suggests that γ₂ might be non-zero for massive haloes (Baldauf et al. 2012; Chan, Scoccimarro & Sheth 2012; Sheth, Chan & Scoccimarro 2013). In this regards, our approach can furnish a useful cross-check of these results since it can provide a measurement of γ₂ which is independent of the bispectrum.

5 CONCLUSION

Dark matter haloes and galaxies are inherently biased relative to the mass density field, and this bias can manifest itself not only in n-point statistics such as the power spectrum or bispectrum, but also in simpler one-point statistics. In this work, we took advantage of this to ascertain the importance of certain non-local Lagrangian bias factors independently of a two-point measurement. We extended the cross-correlation technique of Musso et al. (2012) to χ²-distributed variables, focusing on the quadratic terms |$\eta ^2(\mathrm{\boldsymbol x})$| and |$\zeta ^2(\mathrm{\boldsymbol x})$| (see equation 20) which arise from the peak constraint and for which we have theoretical predictions. In principle, however, our approach could be applied to measure the Lagrangian bias factor associated with any χ²-distributed variable such as the tidal shear for instance. We validated our method with peaks of Gaussian random field before applying it to a catalogue of dark matter haloes with mass M > 10¹³ M_⊙ h⁻¹. Including an offset between the protohalo centre-of-mass and the peak position in the modelling (motivated by the analysis of Ludlow & Porciani 2011), we were able to reproduce our measurements of the non-local bias |$\sigma _1^2\chi _{10}$|⁠. Our result χ₁₀ < 0 is consistent with the findings of Ludlow & Porciani (2011), who demonstrated that protohaloes with M > 10¹³ M_⊙ h⁻¹ preferentially form near initial density peaks (χ₁₀ ≡ 0 for a random distribution). However, we were unable to explain the measurements of |$\sigma _2^2\chi _{01}$|⁠, even with the additional assumption that a fraction of the haloes collapse from filamentary-like structures rather than density peaks. We speculate that a dependence of the halo Lagrangian bias on |$s_2({\boldsymbol x})$| might be needed to explain this discrepancy.

The dependence on |$\eta ^2({\boldsymbol x})$| induces a correction |$-2\chi _{10}(\mathrm{\boldsymbol k}_1\cdot \mathrm{\boldsymbol k}_2)$| to the halo bias which, for collinear wavevectors |$\mathrm{\boldsymbol k}_1$| and |$\mathrm{\boldsymbol k}_2$| of wavenumber 0.1 h⁻¹ Mpc, is Δb ≈ 0.02 (0.05) and ≈0.30 (0.88) for haloes of mass M = 10¹³ and 10¹⁴ M_⊙ h⁻¹ at redshift z = 0 (1), respectively. Relative to the evolved, linear halo bias |$b_1^{\rm E}\equiv 1 + b_{100}$|⁠, the fractional correction is |$\Delta b/b_1^{\rm E}\sim$|2 per cent and ∼15 per cent for the same low and high halo mass in the redshift range 0 < z < 1. Hence, this correction can safely be ignored for M = 10¹³ h⁻¹ Mpc, but it could become relevant at larger halo masses.

We also refined the ESP approach of PSD so that clustering statistics can be straightforwardly computed from the (effective) bias expansion equation (19) (following the prescription detailed in Desjacques 2013). We checked that the predicted halo mass function, from which all the bias factors can be derived, agrees well with the numerical data. However, some of the model ingredients, especially the filtering of the density field, will have to be better understood if one wants to make predictions that are also accurate at small scales.

VD would like to thank the Perimeter Institute for Theoretical Physics and CCPP at New York University for their hospitality while some of this work was being completed there. MB, KCC and VD acknowledge support by the Swiss National Science Foundation.

To facilitate the comparison with other studies, we will call non-local terms all contributions to Lagrangian clustering that are not of the form |$\delta ^n(\mathrm{\boldsymbol x})$|⁠, where |$\delta (\mathrm{\boldsymbol x})$| is the linear mass density field.

We thank Marcello Musso for pointing this out to us.

REFERENCES

Baldauf

Seljak

Desjacques

McDonald

. ,

Phys. Rev. D

2012

, vol.

pg.

083540

Month:	Total Views:
January 2017	1
February 2017	1
July 2017	2
August 2017	1
September 2017	1
October 2017	3
December 2017	8
January 2018	9
February 2018	5
March 2018	10
April 2018	8
May 2018	5
June 2018	3
July 2018	14
August 2018	7
September 2018	6
October 2018	6
November 2018	3
December 2018	5
January 2019	4
February 2019	3
March 2019	9
April 2019	5
May 2019	6
June 2019	3
July 2019	6
August 2019	3
September 2019	7
October 2019	5
November 2019	8
January 2020	2
February 2020	1
March 2020	5
April 2020	5
May 2020	1
June 2020	1
July 2020	5
August 2020	6
September 2020	5
October 2020	1

Article Contents

Abstract

1 INTRODUCTION

2 EXCURSION SET PEAKS

2.1 Notation

2.2 First-crossing and moving barrier

2.2.1 Summary of previous results

2.2.2 Modifications to PSD

2.3 Comparison with numerical simulations

2.4 Bias parameters

3 BIASES FROM CROSS-CORRELATION: EXTENSION TO χ2 VARIABLES

3.1 Bias factors bijk: Hermite polynomials

3.2 Bias factors χij: Laguerre polynomials

4 TEST WITH NUMERICAL SIMULATIONS

4.1 Peaks of Gaussian random fields

4.2 Dark matter haloes

4.3 Interpretation of the measurements

5 CONCLUSION

REFERENCES

APPENDIX A: THE CURVATURE FUNCTION OF DENSITY PEAKS

APPENDIX B: BIVARIATE χ2 DISTRIBUTIONS

Citations

Views

Altmetric

Email alerts

Astrophysics Data System

Citing articles via

Latest

Most Read

Most Cited

This Feature Is Available To Subscribers Only

3 BIASES FROM CROSS-CORRELATION: EXTENSION TO χ² VARIABLES

3.1 Bias factors b_ijk: Hermite polynomials

3.2 Bias factors χ_ij: Laguerre polynomials

APPENDIX B: BIVARIATE χ² DISTRIBUTIONS