On the correlations of galaxy peculiar velocities and their covariance

The rate of change of a matter density perturbation $\delta_{m}$ with time may be related to the surrounding velocity field using the continuity equation, in terms of the growth rate of cosmic structure $f=d(\ln{\delta_{m}})/d(\ln{a})$, where $a$ is the cosmic scale factor. The components of the 3D velocity field, $v_{i}(\mathbf{x})$, at vector position $\mathbf{x}$ relative to an origin, are given in “linear theory” by,

where $\tilde{v}_{i}(\mathbf{k})$ and $\tilde{\delta}_{m}(\mathbf{k})$ are the Fourier amplitudes of these fields in terms of wavevector $\mathbf{k}$, and $H=\dot{a}/a$ is the Hubble parameter. Eq.1 assumes linear perturbation theory and irrotational velocity fields (for a derivation, see for example Adams & Blake (2017), Appendix C). The Fourier amplitudes of matter fluctuations are related to the underlying matter power spectrum, $P_{m}(k)=\langle|\tilde{\delta}_{m}(\mathbf{k})|^{2}\rangle\,V$, which may be predicted by cosmological theories, where $V$ is the volume of the normalising Fourier cuboid and the angled brackets $\langle...\rangle$ indicate an average over many statistical realisations. We can trace the matter overdensity field through the locations of discrete galaxies. In a linear bias model, the galaxy overdensity at a location is $\delta_{g}(\mathbf{x})=b\,\delta_{m}(\mathbf{x})$, where $b$ is the galaxy bias factor.

Information about the matter power spectrum and growth rate can be gleaned by measuring the auto-correlation and cross-correlation functions of the galaxy overdensity and velocity fields. We define the galaxy auto-correlation function, $\xi_{gg}(\mathbf{r})$, as a function of vector separation $\mathbf{r}$, in terms of the galaxy overdensity field as,

where its relation to the galaxy power spectrum, $P_{gg}(k)=b^{2}\,P_{m}(k)$, follows by expressing $\delta_{g}(\mathbf{x})$ in terms of its Fourier components. Neglecting redshift-space distortions, $P_{gg}(k)$ does not depend on the direction of $\mathbf{k}$, through isotropy. Given the vector nature of the velocity field, we define the velocity correlation tensor $\psi_{ij}(\mathbf{r})$ as,

where the second equality in terms of the velocity power spectrum, $P_{vv}(k)=\frac{a^{2}H^{2}f^{2}}{k^{2}}\,P_{m}(k)$ follows by substituting in Eq.1. We also define the cross-correlation between the galaxy overdensity and each velocity component, which we write as $\zeta_{i}(\mathbf{r})$,

where we define the galaxy-velocity cross-power spectrum, $P_{gv}(k)=\frac{ibaHf}{k}\,P_{m}(k)$.

We note an important practical detail that measurements of these correlation functions often trace velocities using galaxies, which are preferentially located in overdense regions, hence tracing the “momentum” or mass-weighted correlation function weighted by $1+\delta_{g}(\mathbf{x})$ (Park, 2000; Howlett, 2019). This effect introduces higher-order corrections which we neglect in our current study, in which we choose to focus on large scales.

In order to compare correlation function measurements with theoretical predictions, it is convenient to use summary statistics which are “isotropic” – that is, which depend only on the magnitude of separation $r=|\mathbf{r}|$, not direction $\hat{\mathbf{r}}$. In that case, we may conveniently angle-average measurements over different directions of the separation vector. Neglecting redshift-space distortions, which we will consider in Sec.2.4, the galaxy auto-correlation is an isotropic function of $\mathbf{r}$, which we can see by evaluating Eq.2 as,

where we have used the relation for the spherical Bessel function $j_{0}$,

However, $\psi_{ij}(\mathbf{r})$ and $\zeta_{i}(\mathbf{r})$ are not isotropic functions, and are therefore not useful as a comparison point between measurements and theory.

The velocity correlation tensor $\psi_{ij}(\mathbf{r})$ can be expressed in terms of two isotropic functions $\psi_{\parallel}(r)$ and $\psi_{\perp}(r)$ as (Gorski, 1988),

where $\delta^{K}_{ij}$ is the Kronecker delta function, and the functions $\psi_{\parallel}(r)$ and $\psi_{\perp}(r)$ are given by,

and,

We have not found a convenient derivation of Eq.7 in the literature, so we provide one in Appendix A.

Using these results, we can readily derive two isotropic summary statistics from the 3D velocity field (Gorski, 1988). First, we may define the correlation function between the inward velocities between two points, $v_{r}(\mathbf{x})=\mathbf{v}(\mathbf{x})\cdot\hat{\mathbf{r}}=\sum_{i}v_{i}(\mathbf{x})\,r_{i}/r$,

where the final equality is derived by substituting in Eq.7. (We emphasise that throughout this paper, $v_{r}$ represents the approach velocity of one galaxy towards another along separation vector $\mathbf{r}$, not the radial velocity with respect to an observer that we introduce below.) Second, we can form the isotropic total velocity correlation function,

For the galaxy-velocity cross-correlation function, we can define an isotropic correlation between galaxy position and inward velocity (Fisher, 1995; Okumura et al., 2014; Adams & Blake, 2017),

where we have used the relation for the spherical Bessel function $j_{1}$,

When analysing observational data, we measure the line-of-sight velocity $u(\mathbf{x})=\mathbf{v}(\mathbf{x})\cdot\hat{\mathbf{x}}=\sum_{i}v_{i}(\mathbf{x})\,x_{i}/x$, where $\mathbf{x}$ is the vector position of the velocity tracer with respect to the observer. Since the projection of the 3D velocity to the line-of-sight velocity depends on the position relative to the observer, the theoretical expressions now depend on the distribution of the velocity tracers through space. We represent this by writing the correlation functions as estimators over the survey volume $V$ (this was implicitly present in our definitions of the 3D correlations in Sec.2.1 and Sec.2.2, but not salient because the statistics were independent of position $\mathbf{x}$).

We first consider the line-of-sight velocity auto-correlation function, which we estimate as,

where we write $\mathbf{y}=\mathbf{x}+\mathbf{r}$ in this section. Substituting in Eq.7, the expectation value of the line-of-sight velocity correlation function is,

where $\langle\cos{\theta_{12}}\rangle$ represents the average of the cosine of the angle subtended by the position vectors $\mathbf{x}$ and $\mathbf{y}$ at the observer, and $\langle\cos{\theta_{1}}\cos{\theta_{2}}\rangle$ is the average of the product of the cosines between the position and separation vectors. (We refer to Fig.1 of Turner et al. (2021) for a useful visualisation of the geometry.) Note that this expression encodes the exact curved-sky effects which arise from the variation of the line-of-sight direction with position. In a flat-sky approximation, $\cos{\theta_{12}}=1$ and $\cos{\theta_{1}}=\cos{\theta_{2}}$.

In order to gain insight into Eq.15 and the nature of the $\psi_{\perp}(r)$ and $\psi_{\parallel}(r)$ functions, it is useful to consider the correlation in line-of-sight velocities of a single pair of velocity tracers in the flat-sky approximation. If the separation vector of that pair is perpendicular to the line-of-sight, then $\cos{\theta_{i}}=0$ and $\langle u_{1}\,u_{2}\rangle=\psi_{\perp}(r)$. If the separation vector is parallel to the line-of-sight, then $\cos{\theta_{i}}=1$ and $\langle u_{1}\,u_{2}\rangle=\psi_{\parallel}(r)$. This analysis shows that the radial velocity auto-correlation function is dependent on the direction of the separation vector, and therefore angle-averaging $\hat{\xi}_{uu}(\mathbf{r})$ may lose information.

In order to extract more information from the line-of-sight velocity correlations, we may utilise the $\psi_{1}(r)$ and $\psi_{2}(r)$ statistics defined by Gorski et al. (1989), which apply pairwise weights for each statistic depending on the angles subtended with respect to the observer and separation vectors,

and,

where $N_{\psi_{1}}$ and $N_{\psi_{2}}$ are normalisation constants we identify below. Using Eq.7, the expectation value of Eq.16 is given by,

Defining for each separation bin the normalisation, $N_{\psi_{1}}=\frac{1}{\langle\cos^{2}{\theta_{12}}\rangle}$, and the geometry factor, $A(r)=\frac{\langle\cos{\theta_{1}}\cos{\theta_{2}}\cos{\theta_{12}}\rangle}{\langle\cos^{2}{\theta_{12}}\rangle}$, we find,

Likewise for $\psi_{2}$,

Defining the normalisation, $N_{\psi_{2}}=\frac{1}{\langle\cos{\theta_{12}}\,\cos{\theta_{1}}\,\cos{\theta_{2}}\rangle}$, and the geometry factor, $B(r)=\frac{\langle\cos^{2}{\theta_{1}}\cos^{2}{\theta_{2}}\rangle}{\langle\cos{\theta_{12}}\,\cos{\theta_{1}}\,\cos{\theta_{2}}\rangle}$, we find,

In the flat-sky approximation, $\psi_{1}(r)$ and $\psi_{2}(r)$ are linear combinations of the first two multipoles of the expansion of $\xi_{uu}(\mathbf{r})$ about the angle to the line-of-sight ($\psi_{1}=\xi^{0}_{uu}$, $\psi_{2}=\xi^{0}_{uu}+\frac{2}{5}\xi^{2}_{uu}$), which encode all the angular information.

We now consider the cross-correlation function between the line-of-sight velocity and galaxy overdensity. Since its monopole is zero, we consider the dipole contribution, which Turner et al. (2021) defined as,

Using Eq.12, this function has expectation value,

Hence, with the identification $N_{\psi_{3}}=\frac{1}{\langle\cos^{2}{\theta_{12}}\rangle}$, we have $\langle\psi_{3}(r)\rangle=\xi_{gv}(r)$.

The isotropy of correlation functions is broken in observational datasets by the effects of redshift-space distortions (RSD), through which the apparent positions of galaxies, as determined by their redshifts, are shifted according to their line-of-sight velocities. The presence of RSD requires us to move to an isotropic basis of correlation function multipoles. RSD has no effect on the velocity auto-correlation function to first order (Dam et al., 2021), hence we only consider redshift-space modifications to the cross-correlation and the galaxy auto-correlation functions. The multipole estimators are, for the galaxy auto-correlation function,

and for the cross-correlation function,

(where we can identify $\xi^{1}_{gu}=\psi_{3}$).

To develop the relations between these correlation function multipoles and the underlying power spectra, it is convenient to express the dependence of the power spectra on the angle with respect to the line-of-sight using multipoles, such that the power spectra become position-dependent (e.g, Blake, 2019),

The expectation values of these estimators are given by, for the galaxy auto-correlation function multipoles (e.g., Satpathy et al., 2017),

and for the cross-correlation multipoles (Adams & Blake, 2020; Lai et al., 2023; Turner et al., 2023),

with coefficients,

where the matrix in Eq.29 is a Wigner 3j-symbol. For reference we give the linear-theory terms entering these predictions (Adams & Blake, 2020):

In our study we focus on the subset of correlation function multipoles which can currently be measured with high signal-to-noise: $\xi_{gg}^{0}$, $\xi_{gg}^{2}$ and $\xi_{gu}^{1}$.

In the following analysis, we hence consider two groups of correlation functions. First, “simulation” correlations, in which we utilise 3D velocities and analyse the three correlation functions $[\xi_{gg},\xi_{gv},\xi_{vv}]$. Second, “observation” correlations, in which we utilise line-of-sight velocities and redshift-space multipoles, and analyse the five correlation functions $[\xi^{0}_{gg},\xi^{2}_{gg},\xi^{1}_{gu},\psi_{1},\psi_{2}]$.

In order to evaluate the covariance of a correlation measurement in a realistic scenario, we need to adapt our correlation models to include the effects of the survey selection function, measurement noise, and any position-dependent optimal weights that are applied. We model two separate discrete datasets: a “density” sample which is used to trace the galaxy overdensity $\delta_{g}$, and a “velocity” sample which is used to probe the 3D or radial velocities. We assume that:

• •

  The density and velocity samples may have different galaxy number densities, which vary with location as $n_{g}(\mathbf{x})$ and $n_{v}(\mathbf{x})$.

• •

  Discrete density and velocity tracers are assigned position-dependent weights $w_{g}(\mathbf{x})$ and $w_{v}(\mathbf{x})$.

• •

  Measurement of the galaxy overdensity and velocity fields contains Poisson noise from discrete objects.

• •

  Tracer velocities contain measurement noise, which may be applied either to all velocity components (for a simulation), or in the line-of-sight direction (for observations).

In this context, the correlation between the (weighted, noisy) galaxy overdensity field at two positions $\mathbf{x}$ and $\mathbf{y}$ becomes (Feldman et al., 1994),

where $f_{g}(\mathbf{x})=w_{g}(\mathbf{x})\,n_{g}(\mathbf{x})\,V$ is the selection function of the density tracers. We model the uncertainty in the galaxy overdensity using Poisson noise, for which $\sigma_{g}(\mathbf{x})=w_{g}(\mathbf{x})\sqrt{n_{g}(\mathbf{x})\,V}$ (Feldman et al., 1994; Blake, 2019). We note here that through the inclusion of the normalising volume $V$, $f_{g}(\mathbf{x})$ and $\sigma_{g}(\mathbf{x})$ are dimensionless functions. The equivalent expression for the correlation of two velocity components $v_{i}(\mathbf{x})$ and $v_{j}(\mathbf{y})$ is,

Here, $f_{v}(\mathbf{x})=w_{v}(\mathbf{x})\,n_{v}(\mathbf{x})\,V$ is the selection function of the velocity tracers, and $\sigma_{v}(\mathbf{x})=\epsilon_{v}(\mathbf{x})\,w_{v}(\mathbf{x})\,\sqrt{n_{v}(\mathbf{x})\,V}$ is the velocity noise at a given location as a function of the individual-tracer velocity error $\epsilon_{v}(\mathbf{x})$, where this noise is applied in the radial direction from an observer. This last statement can be established by considering,

where the final step uses $\sum_{i}\frac{x_{i}^{2}}{x^{2}}=1$. If velocity noise is instead applied independently to each velocity component (for example, the non-linear contribution to the velocities in a simulation), the final term in Eq.33 would be $\sigma_{v}^{2}(\mathbf{x})\,\delta^{K}_{ij}\,\delta_{D}(\mathbf{x}-\mathbf{y})$. Finally, the expression for the correlation between the galaxy overdensity field $\delta_{g}(\mathbf{x})$ and a component of velocity $v_{i}(\mathbf{y})$ is,

In this case there is no noise term, because the noise contribution to $\delta_{g}$ and $v_{i}$ is uncorrelated.

To calculate the analytical covariance we will apply several approximations, stated as follows:

• •

  We will assume that the galaxy overdensity and velocity fields are Gaussian random fields.

• •

  We will neglect higher-order correlations; for example, the fact that velocity tracers preferentially sample overdense locations.

• •

  We will neglect the variation of the survey selection function and noise fields on the scale of the separation vector.

• •

  When correlating line-of-sight velocities we will assume the “local plane-parallel” approximation, that the position vectors of points on the scale of the separation are parallel (where this direction varies for different pairs).

These approximations become increasingly inaccurate for small separations, or when the separations are large compared to the radial distances. However, we will show in the remainder of this paper that these approximations produce analytical covariances with a significant range of validity.

Applying these approximations, we restate our estimators of the correlation function statistics in terms of the noisy, weighted fields, in the form that we will use to compute the analytical covariances. Starting with the galaxy correlation function we have,

where $N_{\xi_{gg}}$ is a normalisation constant. Using Eq.32, we can confirm that $\langle\hat{\xi}_{gg}\rangle=\xi_{gg}$ if $N_{\xi_{gg}}=\int\frac{d^{3}\mathbf{x}}{V}\,f_{g}(\mathbf{x})\,f_{g}(\mathbf{x}+\mathbf{r})\approx\int\frac{d^{3}\mathbf{x}}{V}f_{g}^{2}(\mathbf{x})$, where we have applied the approximation that the selection function does not vary on the separation scale. We can average the above estimator over directions by writing,

Generalising to the galaxy correlation function multipoles,

For the inward velocity auto-correlation function we can similarly write,

where using the same approximation as above, $N_{\xi_{vv}}\approx\int\frac{d^{3}\mathbf{x}}{V}f_{v}^{2}(\mathbf{x})$. For the cross-correlation between inward velocity and galaxy position,

where $N_{\xi_{gv}}=\int\frac{d^{3}\mathbf{x}}{V}f_{g}(\mathbf{x})\,f_{v}(\mathbf{x})$. We now turn to the correlations involving line-of-sight velocities. The estimators for the $\psi_{1}$, $\psi_{2}$ and $\psi_{3}$ statistics have already been provided as Eq.16, Eq.17 and Eq.22 above where, for evaluating the covariance, we additionally apply the local plane-parallel approximation $\hat{\mathbf{x}}\approx\hat{\mathbf{y}}$. For the cross-correlation multipole statistics we use,

These estimators form the basis of our analytical covariance evaluations.

We illustrate the calculation of the analytical covariance using the example of $\psi_{1}$. First, we express the covariance in terms of moments of the fields, using the estimator defined in Eq.16, noting in the following that $\mathbf{y}$ is now a general position vector and does not satisfy $\mathbf{y}=\mathbf{x}+\mathbf{r}$:

where the last step uses Wick’s theorem for a Gaussian random field. We now substitute in the correlation model from Eq.33, where we will just show the calculation of the first term for clarity, noting that the second term produces an identical result where $\mathbf{r}-\mathbf{s}$ is replaced by $\mathbf{r}+\mathbf{s}$:

where we have again applied the approximation that the selection function and noise do not vary on the scale of the separations, i.e. $f_{v}(\mathbf{x}+\mathbf{r})\approx f_{v}(\mathbf{x})$ and $\sigma_{v}(\mathbf{x}+\mathbf{r})\approx\sigma_{v}(\mathbf{x})$. Now, we convert the evaluation to Fourier space, by substituting in the expression for the velocity correlation tensor from Eq.3 and for the delta function, $\delta_{D}(\mathbf{r})=\int\frac{V\,d^{3}\mathbf{k}}{(2\pi)^{3}}e^{-i\mathbf{k}\cdot\mathbf{r}}$: