# Sample Variance of the Higher-Order Cumulants of Cosmic Density and Velocity Fields

## 1 Introduction

It is now commonly believed that the large-scale structure in the Universe observed today is explained by gravitational evolution of small initial density inhomogeneities (Peebles 1980). These initial fluctuations are usually assumed to obey random Gaussian distribution which is not only plausible from the central limit theorem but is also predicted by standard inflation models (Guth & Pi 1982, Hawking 1982, Starobinski 1982). Quantitative analysis of statistical measures of present-day cosmic fields are very important to confirm or disprove the structure formation scenario based on gravitational instability from primordially Gaussian fluctuations. Their higher-order cumulants provide useful tools for this purpose.

In addition to the anisotropy of the cosmic microwave background radiation (CMB), distribution of galaxies and the peculiar velocity field are basic measures for the statistical analysis of the large-scale inhomogeneities in the Universe. The former has been widely investigated and there are two ongoing large redshift surveys now, namely, the Anglo-Australian 2dF Survey (Colless 1998) and the Sloan Digital Sky Survey (Gunn & Weinberg 1995) which are expected to revolutionarily improve our knowledge of three-dimensional galaxy distribution.

We should notice, however, that what we can directly observe is the distribution of galaxies whereas what we can discuss from the first principle is the distribution of the underlying matter. In spite of the rapid increase of observational data, our understanding of the relation between distribution of galaxies and that of underlying gravitating matter, namely biasing (Kaiser 1984), is far from satisfactory. This hampers straight forward comparison between theories and observations.

In contrast to the number-density field of galaxies, the cosmic velocity field reflects the dynamical nature of underlying matter fluctuation and is basically independent of the poorly understood biasing relation at least on large scales (Dekel 1994). This is a fundamental merit of the cosmic velocity field. On the other hand, we must point out that the survey depth of the comic velocity field is currently limited only to $`L70h^1\mathrm{Mpc}`$ around us even in the case of recent Mark III Catalog of Galaxy Peculiar Velocities (Willick et al. 1997), which is much smaller than that of the redshift surveys<sup>1</sup><sup>1</sup>1Bernardeau et al. (1995) used $`L40h^1\mathrm{Mpc}`$ as a practical current limit of high-quality data of the peculiar velocity field.. Therefore uncertainties due to the finiteness of our survey volume, or the sample variance, must inevitably become large and are therefore very important in the analysis of the cosmic velocity field.

The higher-order cumulants of velocity divergence field have been extensively investigated in the framework of nonlinear perturbation theory, and are expected to work as useful quantities in observational cosmology, for example, to constrain the density parameter independent of the biasing (Bernardeau et al. 1995, 1997). In this Letter we investigate the sample variance of the higher-order cumulants of the velocity divergence field assuming that the initial fluctuation obeys isotropic random Gaussian distribution. Our formalism is also applicable to the density field and is similar to Srednicki (1993) who analyzed the skewness parameter of cosmic microwave background radiation.

In the previous Letter (Seto & Yokoyama 1998) we discussed the sample variance of the second-order moment or the variance of a component of the peculiar velocity and that of the linear density fluctuation, whose expectation values and sample variances are of the same order in perturbation. In the case of higher-order cumulants, their expectation values vanish in linear theory and are generated from nonlinear mode coupling in higher-order perturbation theory. Nevertheless their sample variance is nonvanishing even at the linear order. Thus the sample variance is expected to be much more important for higher-order cumulants. In this Letter we compare the expectation values of the lowest-order contribution of the higher-order cumulants of these fields and their sample variances predicted by linear theory.

It is true that there are other sources of errors in the observational analysis of these fields. Using Monte Carlo calculations we could basically take various effects into account at one time (e.g. Borgani et al. 1997). But the sample variance due to the finiteness of the survey region can be regarded as a fundamental limitation in the sense that this uncertainty is independent of how accurately we could measure the cosmic fields in a specific region in the Universe. In addition, to investigate the sample variance by means of a numerical simulation, we need a simulation box much larger than the (mock) survey region, as the sample variance is heavily weighted to Fourier modes which are comparable or greater than the survey depth. Considering these two factors, the simple analysis presented in this Letter is a useful and convenient approach to estimate a fundamental limitation in the observational determination of higher-order cumulants of cosmic fields.

## 2 Formulation

We denote the density contrast field by $`\delta (𝒙)`$ and the velocity divergence field by $`\theta (𝒙)H_0^1𝑽(𝒙)`$ where $`𝑽(𝒙)`$ is the peculiar velocity field and $`H_0`$ is the Hubble parameter. At the linear order, which is indicated by the suffix “lin” hereafter, we have the following relation.

$$\theta _{\mathrm{lin}}(𝒙)=f(\mathrm{\Omega }_0)\delta _{\mathrm{lin}}(𝒙),$$
(1)   
where the function $`f`$ is the logarithmic time derivative of the logarithm of the linear growth rate of the density contrast $`\delta _{\mathrm{lin}}(𝒙)`$ and is well fitted by $`f(\mathrm{\Omega }_0)\mathrm{\Omega }_0^{0.6}`$ with $`\mathrm{\Omega }_0`$ being the density parameter (Peebles 1980). We define the linear dispersion of these fields as

$$\sigma ^2\delta _{\mathrm{lin}}^2(𝒙),\sigma _\theta ^2\theta _{\mathrm{lin}}^2(𝒙),$$
(2)   
where $`X`$ represents to take an ensemble average of the a field $`X`$. From equation (1) the linear root-mean-square (RMS) fluctuation of the velocity divergence field $`\sigma _\theta `$ is written in terms of $`\sigma `$ and $`f(\mathrm{\Omega }_0)`$ as

$$\sigma _\theta =f(\mathrm{\Omega }_0)\sigma .$$
(3)   
These two quantities $`\sigma `$ and $`\sigma _\theta `$ work as the expansion parameters for perturbative analysis in this Letter.

In the observational study of the cosmic fields in the framework of perturbation theory a smoothing operation is crucially important to get rid of strong nonlinearities on small scales and noises due to the discreteness of galaxies which work as tracers of these fields. In this Letter we only discuss fields smoothed with a Gaussian filter defined by

$$W(𝒙)\frac{1}{\sqrt{(2\pi R_\mathrm{s}^2)^3}}\mathrm{exp}\left(\frac{𝒙^2}{2R_\mathrm{s}^2}\right).$$
(4)   
Here $`R_\mathrm{s}`$ is the smoothing radius but we omit its explicit dependence in most part of this Letter for notational simplicities.

Next we briefly summarize the expectation values of the third- and forth-order cumulants for two fields, $`\delta (𝒙)`$ and $`\theta (𝒙)`$. We introduce their first-nonvanishing contributions predicted by higher-order (nonlinear) Eulerian perturbation theory (e.g. Peebles 1980). First we assume that the power spectrum of density fluctuation has a power-law form characterized by a single power index $`n(3)`$ as

$$P(k)k^n.$$
(5)   
In this case the third-order cumulants or the skewness of $`\delta (𝒙)`$ and $`\theta (𝒙)`$ smoothed with the Gaussian filter (4) are evaluated perturbatively as follows.

$$\delta (𝒙)^3=S_3(n)\sigma ^4+O(\sigma ^6),\theta (𝒙)^3=S_{3\theta }(n)\sigma _\theta ^4+O(\sigma _\theta ^6),$$
(6)   
where the factors $`S_3(n)`$ and $`S_{3\theta }(n)`$ are of order unity and have been given by Matsubara (1994) and Łokas et al. (1995) in terms of the hypergeometric function $`F`$.

$`S_3(n)`$ $``$ $`3F({\displaystyle \frac{n+3}{2}},{\displaystyle \frac{n+3}{2}},{\displaystyle \frac{3}{2}},{\displaystyle \frac{1}{4}})\left(n+{\displaystyle \frac{8}{7}}\right)F({\displaystyle \frac{n+3}{2}},{\displaystyle \frac{n+3}{2}},{\displaystyle \frac{5}{2}},{\displaystyle \frac{1}{4}}),`$ (7)   
$`S_{3\theta }(n)`$ $``$ $`{\displaystyle \frac{1}{f(\mathrm{\Omega }_0)}}\left[3F({\displaystyle \frac{n+3}{2}},{\displaystyle \frac{n+3}{2}},{\displaystyle \frac{3}{2}},{\displaystyle \frac{1}{4}})\left(n+{\displaystyle \frac{16}{7}}\right)F({\displaystyle \frac{n+3}{2}},{\displaystyle \frac{n+3}{2}},{\displaystyle \frac{5}{2}},{\displaystyle \frac{1}{4}})\right].`$ (8)   
Here we neglect extremely weak dependence on cosmological parameters except for the function $`f(\mathrm{\Omega }_0)`$ in $`S_{3\theta }(n)`$. In principle, equations (7) and (8) are valid only for a pure power-law spectrum as equation (5) but we extrapolate them to more realistic power spectra with an effective power index $`n(R_\mathrm{s})`$ defined at the smoothing scale by the following equation (see Bernardeau et al. 1995):

$$n(R_\mathrm{s})=3\frac{d\mathrm{ln}\sigma ^2(R_\mathrm{s})}{d\mathrm{ln}R_\mathrm{s}}=3\frac{d\mathrm{ln}\sigma _\theta ^2(R_\mathrm{s})}{d\mathrm{ln}R_\mathrm{s}}.$$
(9)   
In the same manner the forth-order cumulants, or the kurtosis, of $`\delta (𝒙)`$ and $`\theta (𝒙)`$ are written perturbatively as follows.

$`\delta (𝒙)^43\sigma ^4`$ $`=`$ $`S_4(n)\sigma ^6+O(\sigma ^8),`$ (10)   
$`\theta (𝒙)^43\sigma _\theta ^4`$ $`=`$ $`S_{4\theta }(n)\sigma _\theta ^6+O(\sigma _\theta ^8).`$ (11)   
Łokas et al. (1995) derived analytic formulas for $`S_4(n)`$ and $`S_{4\theta }(n)`$ based on higher-order perturbation theory and evaluated them numerically for $`3n1`$. In the present analysis we use fitting formulas for $`S_4(n)`$ and $`S_{4\theta }(n)`$ given in their paper.

In observational cosmology, we usually estimate an ensemble average $`X(𝒙)`$ of a field $`X`$ by taking its volume average $`A(X,𝒱)`$ in an observed patch $`𝒱`$,

$$A(X,𝒱)\frac{1}{𝒱}_𝒱X(𝒙)d^3xX.$$
(12)   
We can commute the ensemble average with a volume integral above to obtain

$$A(X,𝒱)=X.$$
(13)   
Thus the ensemble average of the volume average $`A(X,𝒱)`$ is identical to the universal value $`X`$. However, the observed value $`A(X,𝒱)`$ in one specific patch $`𝒱`$ is expected to fluctuate around its mean $`X`$ because of the spatial correlation and inhomogeneity of the field $`X(𝒙)`$ beyond the patch $`𝒱`$. These fluctuations are nonvanishing even in linear theory and we define its RMS value $`E_{\mathrm{lin}}(X,𝒱)`$ as follows.

$$E_{\mathrm{lin}}(X,𝒱)\left\{A(X_{\mathrm{lin}},𝒱)X_{\mathrm{lin}}\right\}^2^{1/2}.$$
(14)   
Our basic strategy is to compare the magnitude of this linear sample variance $`E_{\mathrm{lin}}(X,𝒱)`$ with the expectation value $`X`$. This fluctuation should be smaller than the expectation value $`X`$; otherwise the particular value of $`A(X,𝒱)`$ obtained in one survey volume would lose its universality and one could not extract any cosmological information from it.

Let us now calculate the linear fluctuation $`E_{\mathrm{lin}}(X,𝒱)`$ for the skewness and the kurtosis of the velocity divergence field. Using the nature of the multivariate Gaussian variables, we obtain the sample variance of the skewness of $`\theta (𝒙)`$ as

$$E_{\mathrm{lin}}^2(\theta ^3,𝒱)=\frac{3\sigma _\theta ^6}{𝒱^2}_𝒱d^3x_𝒱d^3y\mathrm{\Xi }(r_{xy})\{3+2\mathrm{\Xi }(r_{xy})^2\},$$
(15)   
where we have denoted the separation between two points $`𝒙`$ and $`𝒚`$ by $`r_{xy}=|𝒙𝒚|`$ and defined the normalized linear two-point correlation function $`\mathrm{\Xi }(r)`$ as

$$\mathrm{\Xi }(r_{xy})=\frac{\theta _{\mathrm{lin}}(𝒙)\theta _{\mathrm{lin}}(𝒚)}{\sigma _\theta ^2}=\frac{\delta _{\mathrm{lin}}(𝒙)\delta _{\mathrm{lin}}(𝒚)}{\sigma ^2}=_0^{\mathrm{}}\frac{k^2dk}{2\pi ^2\sigma ^2}\frac{\mathrm{sin}kr}{kr}P(k)\mathrm{exp}(k^2R_\mathrm{s}^2).$$
(16)   
In the same manner the linear fluctuation for the forth-order and second-order cumulants are given as follows.

$`E_{\mathrm{lin}}^2(\theta ^43\sigma _\theta ^4,𝒱)`$ $`=`$ $`{\displaystyle \frac{24\sigma _\theta ^8}{𝒱^2}}{\displaystyle _𝒱}d^3x{\displaystyle _𝒱}d^3y\mathrm{\Xi }(r_{xy})^2\{3+\mathrm{\Xi }(r_{xy})^2\},`$ (17)   
$`E_{\mathrm{lin}}^2(\theta ^2,𝒱)`$ $`=`$ $`{\displaystyle \frac{2\sigma _\theta ^4}{𝒱^2}}{\displaystyle _𝒱}d^3x{\displaystyle _𝒱}d^3y\mathrm{\Xi }(r_{xy})^2.`$ (18)   
In the next section we calculate the ratio of $`E_{\mathrm{lin}}(X,𝒱)`$ to the lowest nonvanishing order of $`X`$. Here one should notice that the skewness and kurtosis given in equations (6), (10), and (11) are obtained from higher-order contributions in perturbation in contrast to $`E_{\mathrm{lin}}(X,𝒱)`$ obtained in linear theory. We summarize the order of the expansion parameter $`\sigma _\theta `$ for these ratios below.

$`{\displaystyle \frac{E(\theta ^2,𝒱)}{\theta ^2}}`$ $`=`$ $`O(1),`$ (19)   
$`{\displaystyle \frac{E(\theta ^3,𝒱)}{\theta ^3}}`$ $`=`$ $`O(\sigma _\theta ^1),`$ (20)   
$`{\displaystyle \frac{E(\theta ^43\sigma _\theta ^4,𝒱)}{\theta ^43\sigma _\theta ^4}}`$ $`=`$ $`O(\sigma _\theta ^2).`$ (21)   
So far we have mainly discussed linear fluctuation of the velocity divergence field $`\theta (𝒙)`$ but extension to the case of the density field $`\delta (𝒙)`$ is simple and straightforward (see eq.).

## 3 Results

In this section we calculate the ratio $`E_{\mathrm{lin}}(X,𝒱)/|X|`$ using specific cosmological models. For calculational simplicity we assume that our survey patch $`𝒱`$ is a sphere with radius $`L`$ and volume $`𝒱=(4\pi /3)L^3`$. In this case we can simplify the six-dimensional integral of equations (15), (17), and (18) to a three-dimensional one owing to the rotational symmetry (Seto & Yokoyama 1998).

We investigate two cold-dark-matter (CDM) models with different density parameter $`\mathrm{\Omega }_0`$, namely, $`\mathrm{\Omega }_0=0.3`$ (open model) and $`\mathrm{\Omega }_0=1.0`$ (Einstein de-Sitter model) both with vanishing cosmological constant<sup>2</sup><sup>2</sup>2Cosmological constant is almost irrelevant in our analysis, even if it is nonvanishing within the current observational limit. and the Hubble parameter $`h=H_0/(100\mathrm{k}\mathrm{m}/\mathrm{sec}/\mathrm{Mpc})=0.7`$. As for the initial matter fluctuation we use CDM power spectrum given in Efstathiou et al. (1992) as

$$P(k)=\frac{Bk}{\left\{1+\left[\alpha k+(\beta k)^{3/2}+(\gamma k)^2\right]^\mu \right\}^{2/\mu }},$$
(22)   
where $`\alpha =(6.4/\mathrm{\Gamma })h^1\mathrm{Mpc}`$, $`\beta =(3.0/\mathrm{\Gamma })h^1\mathrm{Mpc}`$, $`\gamma =(1.7/\mathrm{\Gamma })h^1\mathrm{Mpc}`$, $`\mu =1.13`$, and the normalization factor $`B=(96\pi ^2/5)\mathrm{\Omega }_0^{1.54}H_0^4(Q_{rms}/T_0)^2`$ with the current temperature of CMB $`T_0=2.73`$K and the quadruple fluctuation amplitude of it $`Q_{rms}=15.3\mu `$K from 4yr COBE data (Górski et al. 1996). We fix the shape parameter $`\mathrm{\Gamma }`$ by $`\mathrm{\Gamma }=h\mathrm{\Omega }_0`$.

As explained before, the cosmic velocity field is considered to be less contaminated by the poorly understood biasing effect but its survey depth is much smaller than that of the redshift surveys of galaxies. Therefore we mainly consider a typical observational situation of the cosmic velocity field and adopt a Gaussian filter with $`R_\mathrm{s}=12h^1\mathrm{Mpc}`$ following the POTENT analysis (Bertschinger & Dekel 1989, Dekel 1994).

Using the formulas given in the previous section, we plot the sample variance due to the smallness of the survey volume in Fig.1. The expansion parameter $`(\sigma ,\sigma _\theta )`$ is (0.37,0.18) for the open model $`\mathrm{\Omega }_0=0.3`$ and (0.38,0.38) for the Einstein de-Sitter model $`\mathrm{\Omega }_0=1.0`$ in our case. We have $`\theta (𝒙)\delta (𝒙)`$ at the linear order and thus the relative fluctuations of the second-order moments are identical for these two fields (thick solid lines in Fig.1). For the higher-order cumulants, nonlinear mode coupling arises in a different manner for $`\theta (𝒙)`$ and $`\delta (𝒙)`$ and their relative fluctuations are no longer identical.

As is seen in Fig.1, at the current survey depth $`L40h^1\mathrm{Mpc}`$ (Bernardeau et al. 1995), the sample variance of the skewness of the velocity divergence field $`\theta (𝒙)^3`$ remains as much as the expectation value itself. To reduce it smaller than $`30\%`$, we have to take the survey radius $`L`$ as deep as $`L180h^1\mathrm{Mpc}`$ for the open model and $`L150h^1\mathrm{Mpc}`$ for the Einstein de-Sitter model. These values are much larger than the current observational limit. Bernardeau (1995) proposed a method to estimate the density parameter using a relation $`\theta ^3/\theta ^2^2\mathrm{\Omega }_0^{0.6}`$. Analysis above show that even if we could take the survey radius $`L`$ as big as $`150h^1\mathrm{Mpc}`$, the error of such estimation of $`\mathrm{\Omega }_0`$ would remain as large as $`50\%`$.

## 4 Summary

In this Letter we have discussed the magnitude of the sample variance in the observational determination of the reduced higher-order cumulants of smoothed density and velocity divergence fields. We have compared the sample variance predicted by linear theory with the lowest-order nonvanishing contribution to the cumulants assuming that the primordial fluctuation was random Gaussian distributed. We have paid much attention to the velocity divergence field as (i) it is less contaminated by the biasing relation and extensively investigated in the framework of perturbation theory but (ii) velocity survey is currently limited to a relatively small region.

The skewness of the velocity divergence is an interesting quantity to characterize the non-Gaussianity induced by gravity and is expected to constrain the density parameter with small theoretical ambiguities as long as primordial fluctuations are Gaussian distributed. But according to the present analysis we cannot determine the skewness of this field with an error less than $`30\%`$ if our survey depth is not as deep as $`200h^1\mathrm{Mpc}`$.

In the previous Letter (Seto & Yokoyama 1998) we have shown that the peculiar velocity field suffers from much larger sample variance than the linear density field because the former depends on small wavenumber modes much more strongly. On the other hand, the peculiar velocity divergence field discussed here has the same spectral dependence as the density field in linear theory (eq. ). Hence the large relative sample variance we have encountered in the present Letter is entirely due to the fact that it is nonvanishing even in linear theory whereas the expectation values of the higher-order cumulants become nonvanishing only after nonlinear effects are taken into account.

This work was partially supported by the Japanese Grant in Aid for Science Research Fund of the Monbusho Nos. 3161 (N. S.) and 09740334 (J. Y.).