Measuring Bias via the Consistency Relations of the Large Scale Structure

TL;DR

This paper explores consistency relations in large scale structure, demonstrating how they can be used to measure bias and growth rate through observable correlations, validated by N-body simulations.

Contribution

It introduces a novel method to measure large scale halo bias and growth rate using consistency relations, independent of theoretical uncertainties.

Findings

BAO in bispectrum suppressed relative to power spectrum in squeezed limit

Consistency relations validated by N-body simulations

Method enables bias and growth rate measurement from observables

Abstract

Consistency Relations (CR) for the Large Scale Structure are exact equalities between correlation functions of different order. These relations descend from the equivalence principle and hold for primordial perturbations generated by single-field models of inflation. They are not affected by nonlinearities and hold also for biased tracers and in redshift space. We show that Baryonic Acoustic Oscillations (BAO) in the bispectrum (BS) in the squeezed limit are suppressed with respect to those in the power spectrum (PS) by a coefficient that depends on the BS configuration and on the bias parameter (and, in redshift space, also on the growth rate). We test these relations using large volume N-body simulations and show that they provide a novel way to measure large scale halo bias and, potentially, the growth rate. Since bias is obtained by comparing two directly observable quantities, the method is free from theoretical uncertainties both on the computational scheme and on the underlying cosmological model.

Tables

Table 1. Table 1: Best fit values for b m subscript 𝑏 𝑚 b_{m} (expected value, b m = 1 subscript 𝑏 𝑚 1 b_{m}=1 ) at z = 0 𝑧 0 z=0 for different values of the maximum allowed q 𝑞 q .

$q=0.019$ h/Mpc
$\mu$	$b_m$	${\sigma }_b$	${\tilde{\chi }}^2$
0.9	1.04	0.05	1.28
0.7	0.94	0.05	0.56
0.5	1.0	0.1	0.78
0.3	1.35	0.47	0.65

Table 2. Table 2: Bias values for halos with M > M m ​ i ​ n = 10 13 ​ M ⊙ 𝑀 subscript 𝑀 𝑚 𝑖 𝑛 superscript 10 13 subscript 𝑀 direct-product M>M_{min}=10^{13}M_{\odot} at z = 0 𝑧 0 z=0 . The expected value from Eq. ( 6 ) is b h = 1.46 ± 0.03 subscript 𝑏 ℎ plus-or-minus 1.46 0.03 b_{h}=1.46\pm 0.03 .

$q=0.019$ h/Mpc
$\mu$	$b_h$	${\sigma }_b$	${\tilde{\chi }}^2$
0.9	1.40	0.09	0.8
0.7	1.3	0.1	0.6
0.5	1.7	0.4	1.1
0.3	2.5	1.9	0.6

Table 3. Table 3: Bias values for halos with M > M m ​ i ​ n = 10 13 ​ M ⊙ 𝑀 subscript 𝑀 𝑚 𝑖 𝑛 superscript 10 13 subscript 𝑀 direct-product M>M_{min}=10^{13}M_{\odot} at z = 0 𝑧 0 z=0 , extracted from PT fits including one, two, and three bias parameters. Results for two different values for k max subscript 𝑘 max k_{\rm max} are reported, with k max low = 0.08 superscript subscript 𝑘 max low 0.08 k_{\rm max}^{\rm low}=0.08 h Mpc -1 and k max high = 0.3 superscript subscript 𝑘 max high 0.3 k_{\rm max}^{\rm high}=0.3 h Mpc -1 . The expected value for b h subscript 𝑏 ℎ b_{h} defined in Eq. ( 6 ) is b h = 1.46 ± 0.03 subscript 𝑏 ℎ plus-or-minus 1.46 0.03 b_{h}=1.46\pm 0.03 .

bias model	$b_1\left(k_{\max }^{\mathrm{low}}\right)$	$b_1\left(k_{\max }^{\mathrm{high}}\right)$	$b_2\left(k_{\max }^{\mathrm{low}}\right)$	$b_2\left(k_{\max }^{\mathrm{high}}\right)$	$b_s\left(k_{\max }^{\mathrm{low}}\right)$	$b_s\left(k_{\max }^{\mathrm{high}}\right)$
$b_1$	1.347	1.389	-	-	-	-
$b_1,b_2$	1.336	1.357	0.0355	0.1332	-	-
$b_1,b_2,b_s$	1.393	1.497	-0.1021	-0.3468	-0.1371	-0.3447

Equations

$$\nabla {\varphi }^{\prime }({\mathbf{x}}^{\prime },\tau )=\nabla \varphi (\mathbf{x},\tau )-\mathcal{H}\dot{\mathbf{d}}(\tau )-\ddot{\mathbf{d}}(\tau ),$$

$$\tilde{\mathbf{d}}(\mathbf{q},\tau )={\int }^{\tau }d{\tau }^{\prime }\mathbf{v}(q,{\tau }^{\prime })=\frac{\mathbf{v}(q,\tau )}{\mathcal{H}f}=i\frac{\mathbf{q}}{q^2}{\delta }_m(\mathbf{q},\tau ),$$

$$B_{\alpha \beta \gamma }(q,k_{+},k_{-};{\tau }_{\alpha },{\tau }_{\beta },{\tau }_{\gamma })\equiv$$

$$⟨{\delta }_{\alpha }(\mathbf{q};{\tau }_{\alpha }){\delta }_{\beta }(-{\mathbf{k}}_{+};{\tau }_{\beta }){\delta }_{\gamma }({\mathbf{k}}_{-};{\tau }_{\gamma }){⟩}^{\prime },$$

$$B_{\alpha \beta \gamma }(q,k_{+},k_{-};{\tau }_{\alpha },{\tau }_{\beta },{\tau }_{\gamma })\simeq \frac{\mathbf{k}\cdot \mathbf{q}}{q^2}{P}_{\alpha m}(q;{\tau }_{\alpha },{\tau }_{\alpha })\left[\frac{D({\tau }_{\beta })}{D({\tau }_{\alpha })}{P}_{\beta \gamma }(k_{-};{\tau }_{\beta },{\tau }_{\gamma })-\frac{D({\tau }_{\gamma })}{D({\tau }_{\alpha })}{P}_{\beta \gamma }(k_{+};{\tau }_{\beta },{\tau }_{\gamma })\right]+O\left({\left(\frac{q}{k}\right)}^0\right),$$

$$P_{\alpha \beta }(k;{\tau }_{\alpha },{\tau }_{\beta })\equiv ⟨{\delta }_{\alpha }(\mathbf{k};{\tau }_{\alpha }){\delta }_{\beta }(-\mathbf{k};{\tau }_{\beta }){⟩}^{\prime },$$

$$b_{\alpha }(q;{\tau }_{\alpha })\equiv \frac{P_{\alpha \alpha }(q;{\tau }_{\alpha },{\tau }_{\alpha })}{P_{\alpha m}(q;{\tau }_{\alpha },{\tau }_{\alpha })}.$$

$$\underset{q/k\to 0}{\lim }\frac{B_{\alpha \alpha \alpha }(q,k_{+},k_{-})}{P_{\alpha \alpha }(q)P_{\alpha \alpha }(k)}=$$

$$-\frac{{\mu }^2}{b_{\alpha }(q)}\frac{d\log P_{\alpha \alpha }(k)}{d\log k}+O\left({\left(\frac{q}{k}\right)}^0\right),$$

$$\underset{q/k\to 0}{\lim }\frac{B_{mmm}^{\mathrm{SPT}}\left(q,k_{-},k_{+}\right)}{P_m^0(q)P_m^0(k)}=-{\mu }^2\frac{d\log P_m^0(k)}{d\log k}$$

$$+\frac{13+8{\mu }^2}{7}+\mathrm{O}\left(\frac{q}{k}\right),$$

$$P_m^0(k)=P_m^{nw}(k)(1+A(k)\sin (k{r}_s)),$$

$$\displaystyle\!\!\!\!\!\!\!\!\!-\mu^{2}\Big{[}\left(\cos\left(kr_{s}\right)+\alpha(k)\sin\left(kr_{s}\right)\right)\,\frac{kr_{s}A(k)}{1+A(k)\sin\left(kr_{s}\right)}$$

$$\displaystyle\qquad\qquad\qquad\qquad\qquad\qquad+\frac{d\log P^{nw}_{m}(k)}{d\log k}\Big{]}\,,$$

$$\alpha (k)\equiv \frac{1}{k{r}_s}\frac{d\log A(k)}{d\log k},$$

$$p(k,q_{\max },{\mu }^2)=\sum _{i=-2}^n{a}_i(q_{\max },{\mu }^2)k^i.$$

$${\chi }^2=\sum _{i=1}^N\frac{\left(\frac{B}{PP}{|}_i-p_{B_i}\right)-\frac{1}{b_{\alpha }}\left(-{\mu }^2\frac{d\ln P}{d\ln k}{|}_i-p_{P_i}\right)}{{\sigma }_i^2},$$

$$\frac{\mu +{\mu }_k{\mu }_{q}f}{b_h+{\mu }_q^{2}f}\mu ,$$

$$\mathbf{x}+{\mathbf{d}}_{\alpha }(\mathbf{x}).$$

$$B_{\alpha \beta \gamma }(q,k_{+},k_{-})\equiv ⟨{\delta }_{\alpha }(\mathbf{q}){\delta }_{\beta }(-{\mathbf{k}}_{+}){\delta }_{\gamma }({\mathbf{k}}_{-}){⟩}^{\prime }.$$

$${\delta }_{\alpha }(\mathbf{x})\equiv {\tilde{\delta }}_{\alpha }\left(\mathbf{x}-{\mathbf{d}}_{\alpha }(\mathbf{x})\right)\simeq {\tilde{\delta }}_{\alpha }\left(\mathbf{x}\right)-{\mathbf{d}}_{\alpha }(\mathbf{x})\cdot \nabla {\tilde{\delta }}_{\alpha }\left(\mathbf{x}\right),$$

$${\delta }_{\alpha }(\mathbf{p})\simeq {\tilde{\delta }}_{\alpha }(\mathbf{p})-i\int \frac{d^3{q}^{\prime }}{(2\pi {)}^3}\left(\mathbf{p}-{\mathbf{q}}^{\prime }\right)\cdot {\mathbf{d}}_{\alpha }({\mathbf{q}}^{\prime }){\tilde{\delta }}_{\alpha }(\mathbf{p}-{\mathbf{q}}^{\prime }).$$

$$⟨{\tilde{\delta }}_{\alpha }(\mathbf{q}){\mathbf{d}}_{\beta }({\mathbf{q}}^{\prime }){\tilde{\delta }}_{\beta }(-{\mathbf{k}}_{+}-{\mathbf{q}}^{\prime }){\tilde{\delta }}_{\gamma }({\mathbf{k}}_{-}){⟩}^{\prime },$$

$$\simeq ⟨{\tilde{\delta }}_{\alpha }(\mathbf{q}){\mathbf{d}}_{\beta }({\mathbf{q}}^{\prime })⟩⟨{\tilde{\delta }}_{\beta }(-{\mathbf{k}}_{+}-{\mathbf{q}}^{\prime }){\tilde{\delta }}_{\gamma }({\mathbf{k}}_{-}){⟩}^{\prime }$$

$$(2\pi {)}^3{\delta }_D(\mathbf{q}+{\mathbf{q}}^{\prime })⟨{\tilde{\delta }}_{\alpha }(\mathbf{q}){\mathbf{d}}_{\beta }(-\mathbf{q}){⟩}^{\prime }⟨{\tilde{\delta }}_{\beta }(-{\mathbf{k}}_{-}){\tilde{\delta }}_{\gamma }({\mathbf{k}}_{-}){⟩}^{\prime },$$

$$\underset{q/k\to 0}{\lim }{B}_{\alpha \beta \gamma }(q,k_{+},k_{-})\simeq ⟨{\tilde{\delta }}_{\alpha }(\mathbf{q}){\tilde{\delta }}_{\beta }(-{\mathbf{k}}_{+}){\tilde{\delta }}_{\gamma }({\mathbf{k}}_{-}){⟩}^{\prime }+$$

$$\displaystyle\!\!\!\!\!\!\!\!\!i{\bf{k}}\cdot\Big{[}\langle\tilde{\delta}_{\alpha}({\bf{q}}){\bf d}_{\beta}(-{\bf{q}})\rangle^{\prime}P_{\beta\gamma}(k_{-})-\langle\tilde{\delta}_{\alpha}({\bf{q}}){\bf d}_{\gamma}(-{\bf{q}})\rangle^{\prime}P_{\beta\gamma}(k_{+})\Big{]}$$

$$+\cdots ,$$

$${\mathbf{d}}_m(\mathbf{q})=i\frac{\mathbf{q}}{q^2}{\delta }_m(\mathbf{q}),$$

$$\underset{q/k\to 0}{\lim }\frac{B_{\alpha \beta \gamma }(q,k_{+},k_{-})}{P_{\alpha \alpha }(q)P_{\beta \gamma }(k)}=$$

$$\frac{\mu }{b_{\alpha }(q)}\frac{k}{q}\frac{\left(P_{\beta \gamma }(k_{-})-P_{\beta \gamma }(k_{+})\right)}{P_{\beta \gamma }(k)}+O\left({\left(\frac{q}{k}\right)}^0\right)=$$

$$-\frac{{\mu }^2}{b_{\alpha }(q)}\frac{d\log P_{\beta \gamma }(k)}{d\log k}+O\left({\left(\frac{q}{k}\right)}^0\right).$$

$$⟨{\delta }_{\alpha }(-\mathbf{R}){\delta }_{\alpha }\left(\frac{\mathbf{r}}{2}\right){\delta }_{\alpha }\left(-\frac{\mathbf{r}}{2}\right)⟩=⟨{\bar{\delta }}_{\alpha }(-\mathbf{R}+\mathbf{d}(-\mathbf{R})){\bar{\delta }}_{\alpha }\left(\frac{\mathbf{r}}{2}+\mathbf{d}\left(\frac{\mathbf{r}}{2}\right)\right){\bar{\delta }}_{\alpha }\left(-\frac{\mathbf{r}}{2}+\mathbf{d}\left(-\frac{\mathbf{r}}{2}\right)\right)⟩$$

$$=r^j\frac{\partial {\xi }_{\alpha }(r)}{\partial r^i}⟨{\delta }_{\alpha }(-\mathbf{R}){\frac{\partial d^i(\mathbf{r})}{\partial r^j}|}_{\mathbf{r}=0}⟩+\cdots =\frac{2}{3{\mathcal{H}}^2}\frac{r^i{r}^j}{r^2}\frac{\partial {\xi }_{\alpha }(r)}{\partial \ln r}⟨{\delta }_{\alpha }(-\mathbf{R}){\frac{{\partial }^2\Phi (\mathbf{r})}{\partial r^i\partial r^j}|}_{\mathbf{r}=0}⟩+\cdots$$

Full text

††preprint: YITP-19-70

I Introduction

The Large Scale Structure of the Universe (LSS) is governed by nonlinear effects of different nature: the evolution of the dark matter (DM) field, redshift space distortions (RSD), and the bias of the field for the considered tracers (galaxies, halos…) with respect to the DM one. All these effects limit the application of analytical techniques to rather large scales, thus excluding large part of the data from actual analyses. It is therefore remarkable that fully nonlinear statements can be made, in the form of “consistency relations” (CR) Peloso and Pietroni (2013); Kehagias and Riotto (2013). These are statements about the effect of perturbations at large scales on small scales ones, expressed in terms of relations between correlation functions of different order.

The content of the CR’s is essentially kinematic: it describes the effect of a long wavelength displacement field on short distance fluctuations. In linear perturbation theory, assuming adiabatic and gaussian initial conditions, the displacement field is directly related to a long wavelength matter density fluctuation, with a coefficient growing as $1/q$ as $q\to 0$, where $q$ is the long wavelength wavenumber. Since a uniform displacement cannot have any consequence on equal-times correlation functions Scoccimarro and Frieman (1996), CR’s exhibit a charachteristic $1/q$ growth only when considering unequal time correlators. On the other hand, for equal-time correlators the behavior of correlation functions as prescribed by CR’s is in general ‘obscured’ by other terms, whose form is not dictated by the CR’s and are formally of the same order in $q$. Therefore the practical use of CR’s appears to be of a limited extent, as unequal-time correlation functions are not measurable (notice that “unequal-times” does not mean “unequal redshifts”, as the latter involve density fields in regions inside the past lightcone, which is not the case for the former).

In this paper we show that BAO’s provide a way to disentangle the CR-protected terms from the unprotected ones also for equal-times correlators. Indeed, the BAO oscillation feature is amplified in the CR protected terms, which can then be separated from the unprotected ones, which are mostly smooth. This opens the way to the possibility of using the bispectrum CR’s on real data and to extract the cosmological information imprinted on the CR coefficients. The main information we can extract is the large scale cross correlator between the velocity field and the density field for a given tracer. Therefore, by the CR’s we can extract large scale bias for that tracer in a completely model independent, and computational scheme independent, way. In this paper we present this approach by applying it on large volume N-body simulations, validating the use of CR’s on the matter and halo density fields.

The paper is organized as follows. In Sect. II we derive the equal time CR’s and check them both in linear perturbation theory and against simulations. In Sect. III we apply the CR’s to extract the bias of halos at different redshift, while in Sect. IV we compare our approach with the standard perturbative one, and give our conclusions. In Appendix A we give details on the derivations of the equal time CR’s considered in this paper.

II Consistency relations and BAO’s

The CR’s are based on two ingredients. At the dynamical level, the EP, which states that a change in the phase space comoving coordinates from (${\bf{x}}$, ${\bf{p}}$) to (${\bf{x}}^{\prime}$, ${\bf{p}}^{\prime}$), with ${\bf{x}}^{\prime}={\bf{x}}+{\bf{d}}(\tau)$ and ${\bf{p}}^{\prime}={\bf{p}}+am\dot{\bf{d}}(\tau)$, can always be absorbed by a change in the gravitational force from $\nabla\phi({\bf{x}},\tau)$ to

[TABLE]

where $\tau$ is conformal time, ${\bf{d}}(\tau)$ is an arbitrary uniform but time-dependent displacement, dots denote derivatives wrt $\tau$, and ${\cal H}=\dot{a}/a$. We stress that this is an invariance of the Vlasov equation, which describes the phase space evolution beyond the fluid approximation commonly advocated in analytical approaches, such as Perturbation Theory (PT). Therefore, the resulting CR’s hold not only at all PT orders but also beyond that, including all possible non-perturbative effects such as shell-crossing and multistreaming Pietroni (2018).

The second ingredient leading to CR’s comes from relating the displacement ${\bf{d}}(\tau)$ to actual long wavelength velocity modes of the Universe we live in. The connection is done, in Fourier space, by considering a wavenumber dependent displacement,

[TABLE]

where $\delta_{m}({\bf{q}},\tau)$ is the DM overdensity field, and we have used linear PT, assuming it holds small $q$ limit. The effect of a long wavelength velocity field on short scale perturbations is encoded in the *squeezed limit *of the BS

[TABLE]

where ${\bf{k}}_{\pm}={\bf{k}}\pm\frac{{\bf{q}}}{2}$, $q=|{\bf{q}}|$, $k_{\pm}=|{\bf{k}}_{\pm}|$, and the prime indicates that the expectation value has been divided by a $(2\pi)^{3}\delta_{D}(0)$ factor. $\delta_{\alpha,\beta,\gamma}$ indicate the density contrasts for different tracers (e.g. DM, baryons, a given galaxy type, …), evaluated at times $\tau_{\alpha,\beta,\gamma}$, respectively.

Assuming gaussian and adiabatic initial conditions in the growing mode, the BS in the squeezed limit, $q\ll k$, receives a contribution from the long wavelength displacement modes (2), given by Peloso and Pietroni (2014); Baldauf et al. (2015a),

[TABLE]

where the power spectra are defined as

[TABLE]

$D(\tau)$ is the linear matter growth factor and we have assumed the linear behavior at the soft scale $q$, $P_{\alpha m}(q;\tau_{\alpha},\tau_{\beta})=P_{\alpha m}(q;\tau_{\alpha},\tau_{\alpha})D(\tau_{\beta})/D(\tau_{\alpha})$. On the other hand, as we have already emphasised, the dynamics at the hard scale, $k$ is completely nonlinear. The key point is that the structure of the first term at the RHS is protected against any kind of, perturbative and nonperturbative, nonlinear effect. By contrast, the form of the remaining terms, indicated as $O((q/k)^{0})$, is not protected and will be modified in a less and less controllable way at increasing $k$ vaules and decreasing redshifts.

Considering equal times and equal species ($\alpha=\beta=\gamma$), and defining the bias parameter as

[TABLE]

we can rewrite the squeezed limit BS in terms of the logarithmic derivative of the PS, as

[TABLE]

where $\mu\equiv\hat{\bf{k}}\cdot\hat{\bf{q}}$, and we have omitted the time dependence. Notice that, as expected, and verified recently in N-body siumulations Esposito et al. (2019), the equal-time CR’s contains no $1/q$ pole. Finally, we stress that we define bias as the physical quantity in Eq. (6), and not as a parameter in a given bias model, as for instance $b_{1}$ in the usual bias expansion $\delta_{h}=b_{1}\,\delta_{m}+b_{2}\,\delta_{m}^{2}+b_{s}s^{2}\cdots$, (where $s^{2}$ is the “tidal bias” parameter, see e.g., Baldauf et al. (2012)), see the discussion in Sect. IV.

II.1 Check in Perturbation Theory

As a first test of the CR, we compute the matter BS in the squeezed limit at lowest order in SPT. It is given by

[TABLE]

where $P^{0}_{m}\equiv P^{0}_{mm}$ denotes the linear matter PS. Notice that the second term at the RHS is scale independent, as it is proportional to $(P^{0}_{m}(k_{+})+P^{0}_{m}(k_{-}))/P^{0}_{m}(k)\to 2$. On the other hand, the first term, although subdominant, can be isolated from the rest thanks to its scale dependence, which is induced mainly by BAO oscillations. To see this, one writes the PS as

[TABLE]

where $r_{s}$ is the comoving sound horizon and both $P^{nw}_{m}(k)$ and $A(k)$ are featurless functions. Performing the logarithmic derivative, we get the first term at the RHS of (8) as

[TABLE]

where we defined

[TABLE]

which takes values around $10\,\%$ in the $k$ range of interest. So, the squeezed BS contains an oscillating component whose amplitude is enhanced by $kr_{s}\simeq 2\pi k/(0.05\,\mathrm{h\,Mpc^{-1}})$ with respect to Eq. (9) and whose phase is shifted by $\sim\pi/2-\alpha(k)$, and a smooth component given by the second line of Eq (10). The result is shown in Fig. 1, for $q=0.02\,\mathrm{h\,Mpc^{-1}}$ and $q=0.05\,\mathrm{h\,Mpc^{-1}}$ and for three different values for $\mu$. We see that the CR (Eq. (8)) reproduces the BS as long as the squeezed limit ($q/k\ll 1$) holds. When $q/k$ is not small, extra scale dependent terms come into play. However, the amplitudes of the oscillations in the RHS and in the LHS are still related by the CR, as we will show quantitatively in the next section.

II.2 Check in Simulations

Once nonlinear effects are included, the CR (7) ensures that the first line in Eq. (8) is modified just by changing the linear PS with the nonlinear one, while the second line is changed in an uncontrolled way. However, due to parity invariance, the BS is symmetric under $k_{+}\leftrightarrow k_{-}$, or, equivalently, $\mu\to-\mu$, which implies that the leading non-protected term in the squeezed limit should be proportional to $(P_{m}(k_{+})+P_{m}(k_{-}))/P_{m}(k)$, that is, still featureless, although with an unknown coefficient.

To check this explicitly in fully nonlinear dynamics, we use 10 realizations of $N$-body simulations with $N=2048^{3}$ mass elements performed in periodic cubes with the side length of $4\,h^{-1}\mathrm{Gpc}$. The error bars presented in this paper are obtained from the scattering among the 10 realizations.

The mass distribution evolved with the Gadget2 code Springel (2005) and the halo catalogs extracted by Rockstar algorithm Behroozi et al. (2013) at $z=0$ will be presented in what follows (see the “fiducial cosmology” in Nishimichi et al. (2018) for other cosmological/numerical parameters). We measure the bispectrum using the quick FFT-based algorithm presented in Baldauf et al. (2015b) with the aliasing correction following Sefusatti et al. (2016). We first store the data in bins of $(k,q,\mu)$ and check the results after summing up $B(q,k,\mu)/[P(q)P(k)]$ for different $q$ bins up to some $q_{\mathrm{max}}$ weighting by the number of triangles.

We compare the oscillating part of the squeezed BS to that of the logarithmic derivative of the PS, and check if their amplitudes are related by the $-\mu^{2}/b_{\alpha}(q)$ factor of Eq. (7). In order to extract the oscillating part, we will subtract smooth functions (“$poly$” in the figures) of the form

[TABLE]

The negative powers in $k$, at each fixed $q_{\mathrm{max}}$, account for subleading terms in the squeezed limit, of order up to $q_{\mathrm{max}}^{2}/k^{2}$, while the positive powers account for the extra scale dependence induced for instance by the $d\log P^{nw}_{m}/d\log k$ contribution. The coefficients $a_{i}(q_{\mathrm{max}},\mu^{2})$ are fixed such that $p(k,q_{\mathrm{max}},\mu^{2})$ is the best fit to the chosen data (BS or logarithmic derivative); then the function is subtracted from the data to obtain only the oscillatory part, which is not captured by the fit for small enough $n$. We truncate our fit at $n=2$, since with this value a satisfactory reduced $\chi^{2}$ (see below) is obtained.

To quantify the goodness of the CR we measure the bias $b_{\alpha}(q_{max})$ in each $\mu$ bin, by minimizing the $\chi^{2}$ function

[TABLE]

where all the $q_{max}$, $k$ and $\mu$ dependencies are omitted. We denote with $p_{B_{i}}$ and $p_{P_{i}}$ the fitting curves relative to, respectively, the LHS and the RHS of the CR. $\sigma_{i}$ are the errors for the BS in the $i-th$ $k$-bin, since the errors on the power spectra are much smaller. We use $N=26$ linearly spaced bins from $k_{min}=0.05$ h Mpc*-1* to $k_{max}=0.30$ h Mpc*-1*.

First, we performed the test on matter, for which the expected value is $b_{m}(q_{\mathrm{max}})=1$. In Fig. 2 we show the oscillating part of the BS (red points with error bars, estimated from the scatter among the 10 realizations) and of the logarithmic derivative of the PS multiplied by $-\mu^{2}/b_{m}$ (blue lines), for different $\mu$-bins (rows) and for two different values for $q_{max}$ (columns). In tab. 1 we give the best fit values for $b_{m}$ in the different $\mu$-bins and the relative $1-\sigma$ errors and reduced $\chi^{2}$. As expected, going to higher $q_{max}$ the quality of the fits gets worse, both because the squeezed limit is farther, and because the error bars are smaller since the number of BS configuration increases. The best fit values are always compatible with $b_{m}=1$ for q_{max}~{}\mbox{\raisebox{-2.58334pt}{\stackrel{{\scriptstyle<}}{{\sim}}}}~{}0.019\;\mathrm{h/Mpc}.

III Extracting Halo bias

Having tested the CR on the matter field, we then used them to measure the bias of a given halo population, $b_{h}(q_{\mathrm{max}})$, via the $-\mu^{2}/b_{h}(q_{\mathrm{max}})$ modulation of the BAO oscillations, and compare the results with those obtained via the definition in Eq. (6). We considered halos of masses $M>M_{min}=10^{13}M_{\odot}$ at $z=0$. The expected bias, as measured from Eq. (6), is $b_{h}=1.46\pm 0.03$.

The results are shown in Fig. 3, where, comparing the blue lines with the green-dashed ones, we see the effect of the halo bias in reducing the amplitude of the BAO oscillations with respect to the ones present in the logarithmic derivative of the PS. As shown in Fig. 4 and Tab. 2, the extracted values for $b_{h}$ are compatible with the expected one within the error bars, also for $q_{max}=0.061$ h Mpc*-1*, which gave bad fits in the DM case.

We have repeated the analysis for simulation data at redshift $z=1$, finding qualitatively similar results (see Fig. 5). The only nonlinear effect limiting the effectiveness of our method is the damping of the BAO wiggles, which is stronger at lower redshifts.

IV discussion

Extracting bias from a PT-based approach differs substantially from the approach discussed in this Letter. The PT prediction (at tree level as well as at higher orders) requires assuming: 1) a cosmological model, namely, the linear matter power spectrum, its normalization $\sigma_{8}$ and the growth factor, and, 2) a bias model.

On the other hand, we are not assuming any model, neither for the cosmology nor for bias. Moreover, we are not relying on a computational scheme either, being it PT expansion or numerical simulations.

The crucial point is that what we call bias is not a model parameter, but a physical quantity. Thanks to the consistency relation, the cross correlator $P_{\alpha m}(q)$, and then $b_{\alpha}(q)$, can be measured by comparing two directly measurable quantities, namely, the amplitude of the BAO oscillations in the logarithmic derivative of the nonlinear power spectrum and that in the bispectrum, normalized as in Eq. (7). The assumptions beyond our result are just the EP, adiabatic initial conditions, and the linear continuity equation at the very large scale q, which ensure that all tracers, and in particular matter and halos, fall with the same velocity field. This, from a theoretical point of view, singles out our proposal as being of a qualitatively different nature, rooting its robustness in being intrinsically nonperturbative and based on very general physical assumptions.

Moreover, in a PT approach, even at tree level, three bias parameters enter, namely, the linear bias $b_{1}$ and the second order bias parameters, $b_{2}$ and the “tidal bias” $b_{s}$ (see e.g., Baldauf et al. (2012)). If we try to fit the measured bispectrum with this expression, these parameters are degenerate with each other, with the normalization of the power spectrum, and possibly with other cosmological parameters.

Together with the subtleties in the treatment of the scale-dependences of the bias parameters, the end result of the fitting with the tree-level SPT template is clearly not as robust, from a theoretical point of view, as what we are proposing in this work.

To see this explicitly, we produced three plots, see Fig. 6, in which we show the halo bispectrum at z=0 with different values of $\mu$ (in different colors) at a fixed q (=0.04 h Mpc*-1*) as a function of k. The tree level calculation, but with different bias models (one for each figure, as indicated by the title) is fit to the data points, considering two different values of $k_{\rm max}$, below which the fitting was performed, $k_{\rm max}^{\rm low}=0.08$ h Mpc*-1* (solid lines), and $k_{\rm max}^{\rm high}=0.3$ h Mpc*-1* (dotted lines). We can see that, clearly, the tree-level PT does not provide a good fit over the BAO scale. We can see it quantitatively in Fig. 7, which shows the reduced chi-squared for the three bias models as a function of $k_{\rm max}$. Indeed, it is close to unity only at very small $k_{\rm max}$.

For completeness, we summarize below the best fit bias parameters at the two values of $k_{\rm max}$ in Tab. 3.

As we can see, the value of $b_{1}$, and, to a larger extent, those of $b_{2}$ and $b_{s}$, are affected by the choice of the bias model as well as the maximum wavenumber.

Moreover, we must keep in mind that these fits are performed by fixing the cosmological parameters of the PT computation to the true ones used in the simulation. In a real data analysis a scan over these parameters is needed, inducing an higher level of degeneracy.

Considering redshift space, the advantage becomes even more clear. It is known that, besides for very large scales, redshift space distortions, in particular, Fingers of God effects are particularly different to model, requiring some measure of empirical parameterization. Our approach is free from this problem, as the redshift space version of the CR assumes the Kaiser relation only where it is reliable, namely, at the large scale q, while no modelling is required at the short scale k.

Moreover, in redshift space, the possibility of measuring, besides bias, the growth function $f=d\log D/d\log a$ opens. Indeed, following Creminelli et al. (2014) , one realizes that the $\mu^{2}/b_{h}$ prefactor of Eq. (7) now becomes

[TABLE]

where $\mu_{k}=\hat{\bf{k}}\cdot\hat{\bf{z}}$, $\mu_{q}=\hat{\bf{q}}\cdot\hat{\bf{z}}$ are the angles between the two modes and line of sight, here taken to be the $z$ axis. Therefore, by considering different orientations of the triangles, or different multipoles, or again combining different tracers, $f$ and $b_{h}$ can be, at least in principle, measured independently.

We leave this analysis, as well as the quantitative assessment of the feasibility of these measurements in realistic future surveys, and the possible extension to higher order correlators for a forthcoming investigation.

Acknowledgements.

We acknowledge support by the Munich Institute for Astro- and Particle Physics (MIAPP) of the DFG cluster of excellence “Origin and Structure of the Universe”. MP acknowledges support from the EU Horizon 2020 programme under the Marie Sklodowska-Curie grants Invisible- sPlus RISE No. 690575, Elusives ITN No. 674896. TN was supported in part by JSPS KAKENHI Grant Numbers JP17K14273 and 19H00677, and by Japan Science and Technology Agency CREST JPMHCR1414. Numerical simulations and analyses were carried out on Cray XC50 at the Center for Computational Astrophysics, National Astronomical Observatory of Japan.

Appendix A Derivation of the consistency relation

We want to isolate the effect of nearly uniform (but time-dependent) displacements fields,

[TABLE]

on the bispectrum,

[TABLE]

In order to do that, we shift back to the initial positions, defining new density fields,

[TABLE]

which, in Fourier space reads (we use the same symbol for fields in real and in Fourier space),

[TABLE]

Inserting it in (16) we get two contributions containing the expectation values

[TABLE]

and the one obtained by this via the replacements $\beta\leftrightarrow\gamma$ and ${\bf{k}}_{+}\leftrightarrow-{\bf{k}}_{-}$. The CR is obtained by considering the $q\ll k$ limit, and assuming that ${\bf d}_{\beta}({\bf{q}}^{\prime})$ has support for $q^{\prime}\ll k$, and that it correlates with $\tilde{\delta}_{\alpha}({\bf{q}})$, that is, (19) can be written as

[TABLE]

and therefore the bispectrum is given, in this limit, by

[TABLE]

where dots indicate higher orders in the displacement fields.

Now, the EP and adiabatic initial conditions ensure that, at large scales, all species fall with the same acceleration under the effect of the potential generated by the total matter field, that is, ${\bf d}_{\beta}({\bf{q}}),{\bf d}_{\gamma}({\bf{q}})\to{\bf d}_{m}({\bf{q}})$, as $q\to 0$. Moreover, assuming that at large scales linear theory holds, we have

[TABLE]

and therefore the squeezed limit bispectrum can be written as

[TABLE]

where $\mu\equiv({\bf{q}}\cdot{\bf{k}})/(qk)$, and we have defined the bias parameter as in Eq. 6. If the fields $\delta_{\beta}$ and $\delta_{\gamma}$ were evaluated at different times, then the cancellation of the $1/q$ divergences would not be exact Peloso and Pietroni (2013), see Eq. (4). Here there is no such divergence, still the first term in the last line of Eq. (7) is protected by the CR, that is, the coefficient $\mu^{2}/b_{\alpha}(q)$ is not modified by any nonlinear effect, including the non-perturbative ones due to the breakdown of the fluid approximation at small scales (large $k$’s).

It is instructive to derive the configuration space counterpart of Eq. (23). A perfectly uniform displacement field cannot give any contribution to the equal times bispectrum. Therefore, the effect can depend only on the gradient of the large scale displacement/velocity field. More precisely, the $i-th$ component of the displacement field can affect the clustering on short scales along the $j-th$ direction only via its $\partial_{j}d^{i}({\bf{x}})$ component, leading to a contribution to the configuration space three point function proportional to

[TABLE]

where $\xi_{\alpha}(r)$ is the correlation function and $\Phi({\bf{r}})$ the gravitational potential, which has been related to the displacement ${\bf{d}}$ by means of linear PT.