Open Access
Issue
A&A
Volume 663, July 2022
Article Number A124
Number of page(s) 22
Section Cosmology (including clusters of galaxies)
DOI https://doi.org/10.1051/0004-6361/202142526
Published online 19 July 2022

© F. Bernardeau 2022

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

1. Introduction

In the context of cosmological studies, the concept of counts-in-cells statistics has been put forward for a long time as a unique way to quantify the statistical properties of the cosmological fields (White 1979; Colombi et al. 1995; Balian & Schaeffer 1989; Bernardeau & Schaeffer 1999). It was then shown in particular that counts-in-cells statistics, which represents a discrete representation of the local density probability distribution function (PDF), could be directly related to the correlation hierarchy of the density field.

Interest in these types of observables was recently renewed for several reasons. The size of the surveys makes accurately measuring these quantities more realistic. This is already the case for surveys such as the Dark Energy Survey (DES collaboration; Abbott et al. 2018), the Kilo-Degree Survey (KIDS; Heymans et al. 2021), and the Hyper Suprime Cam (HSC; Hikage et al. 2019). The future promises even larger and more powerful surveys such as Euclid (Laureijs et al. 2011; Amendola et al. 2018) and the Rubin Observatory (Ivezić et al. 2019). Moreover, the theoretical foundations for these constructions (at least in the cosmological context) has been considerably strengthen with the realization that the large-deviation theory (LDT; for a general review, see Touchette 2011) could successfully be invoked, as shown in Bernardeau & Reimberg (2016). It clarifies the applicability of the theory to the cosmological density field and places previous works on a much more solid foundation (Valageas 2002; Bernardeau et al. 2014). The ability of density PDF to constrain cosmology was emphasized in Codis et al. (2016b) and completed in Friedrich et al. (2020) and in Uhlemann et al. (2020), who showed that these observable could efficiently constrain the neutrino mass or primordial non-Gaussianities. Finally, although the matter PDF is not a direct observable, as is matter density, it can be closely approached with the help of luminous tracer statistics (Repp & Szapudi 2020), more convincingly in weak-lensing fields, as advocated in numerous recent papers (Barthelemy et al. 2021; Bernardeau & Valageas 2000), or with combined approaches such as density-split statistics (Gruen et al. 2018; Friedrich et al. 2018; Brouwer et al. 2018), which proved to be particularly promising.

The construction of a full theory of these observable requires a detailed analysis of its global error budget, however, due to finite-size surveys, imperfect tracers, and so on. Some of these aspects have been explored in early studies such as Szapudi & Colombi (1996) and Szapudi et al. (1999), but a full theory is still lacking. The developments presented in this paper are made in this context. More precisely, the purpose of this study is to explore what determines the expression of the covariance of data vectors whose elements are local quantities, such as the density contrast and density profiles, in cosmological contexts, that is, in classical random fields with long-range correlations. Derivations were made furthermore assuming statistical homogeneity and isotropy. The domain of application encompasses both counts-in-cells statistics, basically 2D or 3D counts of density tracers, or proxies to projected densities such as mass maps for weak-lensing tomographic observations.

In order to gain insights into the different contributions and the effects that might contribute to the covariance, we rely on the use of the hierarchical models to derive results we think rather general. The immense advantage of using such models is that they naturally incorporate many of the features expected in density cosmological fields (e.g., the magnitude of the high-order correlation functions), and there are also models for which many exact results can be obtained in particular for counts-in-cells statistics. The goal of these constructions is to eventually infer precisely what the performance of PDF measurements would be on the determination of cosmological parameters, taking advantage of results such as those found in Boyle et al. (2021), which give the response function of these observable to various cosmological parameters

Section 2 is devoted to the presentation of the general framework. The subsequent section explores different contributions, from large-scale effects with the derivation of several bias functions to short-distance contributions. Results are derived in a framework as general as possible, including discrete noise associated with the use of a finite number of tracers. Section 4 presents the general hierarchical models, and more specifically, the Rayleigh-Levy flight model that we use as a toy model to evaluate the performances of approximate schemes. In Sect. 5, simplified models for the covariance matrix are presented and evaluated with the help of a set of numerical experiments. Section 6 summarizes the results that have be found and specifies their expected range of application.

The text is complemented by appendices that contain a large amount of material. They present the hierarchical models, their mathematical description, and the mean-field approximation that is used throughout for explicit computations. Appendix C is more specifically devoted to the minimal tree model and the construction of the exact mean-field covariance matrix.

2. General framework. Construction of covariance matrices

The purpose of this section is to show how the elements of the covariance matrix are related to the joint density PDFs within a given survey. We first formalize this relation in a general framework before we explore its consequence in the context we are interested in. We assume we are interested in the PDF of some local quantity, μ, that can be evaluated within a survey, thus defining a field μ(x) throughout the survey. The a priori typical example of this quantity is the density (see below for a more precise illustration of what this quantity could be). The value of μ is assumed to lie in some ensemble ℳ (that can be simply the real numbers), and the data vector we are interested in consists of the probabilities pi that μ lie within the subsets ℳi (which are a priori nonzero within ℳ). The one-point PDF of μ is then given by

(1)

if ℘(μ,x) dμ is the PDF of μ at location x. pi(x) is then assumed to be independent of x in the context we are interested in, for which statistical homogeneity is assumed. More formally, we can define the characteristic function χi(x), which takes the value 1 where μ(x) ∈ ℳi and 0 otherwise.

An estimation of pi would then be given by the volume fraction of the survey where μ(x) ∈ ℳi. We note this estimate as Pi 1,

(2)

which is then itself a random variable, the properties of which we are interested in. More precisely, we would like to derive an operational form for the likelihood function of a set of Pi variable. We limit our investigation here to the construction of the likelihood from the covariance matrix, assuming that the likelihood of Pi is close enough to a Gaussian distribution2.

The ensemble average of Pi is

(3)

We can further define a joint PDF of the same field, ℘(μ,x;μ′,x′), which is the joint PDF of μ and μ′ in locations x and x′. Defining pij(x, x′) as the joint ensemble average of ℘(μ,x;μ′,x′), we have

(4)

The elements of the covariance matrix of Pi are then formally

(5)

This gives the relation between the covariances and joint PDF. If pij(x, x′) depends only on the relative distance between x and x′, this expression can be recast in terms of the distribution of such distances, Ps(rd), in the form

(6)

The precise form of Ps(rd) depends on the detail of the survey. Explicit forms can be given in case of simple regular surveys such as square surveys3. In the context of statistically homogeneous and isotropic random fields, this latter expression is used. In particular, we wish to determine the configurations that contribute most to ij. They obviously depend on both the random processes we consider and on the definition of ℳi and ℳj. In order to be more specific, we assume in the following that μ is a local density assigned to be in bins (i) centered on ρi and with width dρi (assumed a priori to be arbitrarily small), so that

(7)

If necessary, local densities could be obtained after the field μ(x) has been convolved with a window function WR(x), associated with a scale R that is

(8)

It is then assumed R is small compared to the sample size in order to identify what the leading contributions to the joint PDFs might be. In practice, WR might also be a simple top-hat window function, but this is not necessarily the case. It could be more elaborated filters, such as compensated filters (of zero average), such as those introduced for cosmic shear analysis (Schneider 1996; Kaiser 1998; Bernardeau & Valageas 2000).

We furthermore allow the estimated densities ρi to be defined with respect the overall density of the sample ρs,

(9)

For instance, we could be interested in or which are frequently encountered situations in praxis. Then ρs is itself a random variable whose correlation with ρ(x) ought to be taken into account. We then need to explore the properties of either P(ρi, ρj; x, x′) or P(ρs, ρi, ρj; x, x′) from which functions of interest can be built, that is,

(10)

(11)

from which the covariance elements such as

(12)

can be derived and whose properties we wish to explore. We wish in particular to build a model of the likelihood function from such a covariance, requiring full knowledge of its eigenvalues and eigendirections. In this respect, it is implicit that the number of bins (i) to be used is finite. We nonetheless present at least in this first section the results in the continuous limit for ρi. It is finally to be noted that as stated before, we restrict our analysis to covariance matrices, but higher-order correllators might also be considered by generalizing the relation (5) to a higher number of variables.

3. PDF covariances in the context of cosmological models

3.1. Modeling the joint PDF

To make progress, we need to make further assumptions about the mathematical structure of the joint PDF. In the following, we assume in particular that joint PDFs can be obtained from their cumulant generating functions (CGF)4, φ(λi, λj; rd). The latter is defined as

(13)

and it is assumed that this relation can be inverted to give the joint PDF from Laplace inverse transformations,

(14)

where the integrations are made a priori along the imaginary axis. The CGFs themselves are closely related to the averaged correlation functions of the underlying field. In the following, we develop models for which these correlation functions can be computed precisely.

3.2. Large-distance contributions to the joint density PDF

We start by assuming that covariances are dominated by long-range correlation and not by proximity effects (e.g., densities taken in nearby regions). Whether this assumption is correct obviously depends on the particular model and setting we consider, as we detail below. There are large sets of models for which general expressions can be given in this regime. They are the so-called hierarchical models, originally introduced in Peebles (1980), discussed in more detail in Fry (1984a,b), Balian & Schaeffer (1989), Bernardeau & Schaeffer (1992), and further formalized in Bernardeau & Schaeffer (1999) as described below; it is also true in the quasilinear regime as originally pointed out in Bernardeau (1996) and derived in more detail in Codis et al. (2016a). In these regimes, we obtain the following functional form (see previous references and the detailed derivation in Appendix B):

(15)

where ξ(x, x′) is the two-point correlation function of the density field at positions x and x′, and φ0(λ) and φ1(λ) are specific functions of λ that depend on the details of the model.

Then, setting λs to zero, we can easily obtain the expression of the joint PDF at leading order in ξ(rd),

(16)

Here P(ρi) is the one-point density PDF (i.e., implicitly at scale R), and b(ρi) is the density-bias function (to be distinguished from the standard halo-bias function). It also depends on ρi (and on the scale R) so that in the previous expression, the dependence on ρi, ρj, and rd can be factorized out.

To be more precise, P(ρi) is given by the inverse Laplace transform of the CGF (see, e.g., Balian & Schaeffer 1989 and Bernardeau 2013 for a detailed derivation of this inversion),

(17)

where φ0(λ) is the CGF of the density taken at scale R (i.e., for the filter WR). The function b(ρi) is defined through a similar relation,

(18)

The function φ1(λ) can be explicitly computed in the context of perturbation theory calculations (Codis et al. 2016a). This is the case also for models in the so-called hierarchical models (see appendices). In the latter case, these calculations can be extended to higher order, as we describe below, providing ways to better assess the domain of validity of this expansion.

According to the previous relation, this implies that this form translates into the expression of the covariance coefficients for the density PDF. More precisely, we expect

(19)

where is the average value of the two-point correlation function ξ(rd) within the sample.

It is to be noted, however, that this is true if

  • The term in is indeed the leading contribution of the expansion (15). This is obviously not the case for samples with periodic boundary conditions, for which vanishes by construction;

  • The density is defined regardless of the density of the sample. Its expectation value therefore does not coincide with ρs for a given sample.

It can also be noted that in the Gaussian limit, we have b(ρi)=δi/ξ. Applying the relation (16) to the density within one cell and to the sample density ρs = 1 + δs leads then to the following expression for the conditional expression of density PDF,

(20)

This leads to the interpretation of the function b(ρi) as the response function of the density PDF to the sample density. This means that although the density-bias function cannot be derived from the density PDF alone, we should be able to derive it if we are in possession of an operational method to compute the density PDF for arbitrary cosmological parameters (in a way similar to the derivation of halo-bias function as pioneered in Mo & White 1996). Undoubtedly, this result is closely related to the so-called supersample effects (as described for the power spectra covariance in Takada & Hu 2013), that is, the impact of modes of scale comparable to or larger than the sample size. This is not necessarily their only contribution (because subdominant large-scale contributions can also contribute to the covariance), however, but likely to be the most important contribution, as described below.

The density-bias function obeys the following consistency relations:

(21)

(22)

as initially pointed out in Bernardeau & Schaeffer (1992).

3.3. Case of relative densities

The previous formula applies to the local densities, evaluated regardless of the sample density. It does not apply in particular to standard settings (e.g., densities measured out of galaxy counts) where the density is defined with respect to the mean density of the sample. To address this case in particular, we should consider

(23)

as the observable for which the covariance is to be computed. In this case, the formal derivation of the PDFs is presented in the appendix, and it leads to

(24)

(25)

We then use relation (15) to compute the form of this function. At this stage, it is to be noted that the expressions , ∫dx0ξ(x0, x1) and ξ12 all take the same averaged value when integrated over the sample. We note here this common value as . Inserting the resulting expressions of the CGF in both the expressions of and and expanding all terms at linear oder in , we obtain

(26)

This leads to the definition of the first sample-bias function,

(27)

which can be re-expressed in terms of the density-bias function defined in Eq. (18) and the derivative of with respect to

(28)

In this case, the covariance matrix elements are then expected to be given by

(29)

Remarkably, bs1(ρ) can entirely be expressed in terms of b(ρ). For the sake of completeness, we also consider the case of . In this case, it is easy to show that

(30)

(31)

Following the same approach as for the previous case, the leading-order expression in of the connected joint PDF is

(32)

It leads to the definition of the second sample-bias function,

(33)

so that

(34)

The three bias functions are therefore closely related. Although the density-bias function b(ρ) cannot be derived from the shape of P(ρ) alone, as mentioned before, the relations between b(ρ) and either bs1(ρ) and bs2(ρ) depend on the PDF alone. Furthermore, the two relative density bias functions obey the following consistency relations:

(35)

(36)

The second relation is at variance with the corresponding relation (22) for the density-bias function. It indicates that for typical values of ρP(ρ), the sample bias functions, bs#(ρ), are likely to be smaller than the density-bias function b(ρ).

3.4. Structure of the covariance matrix

The consequences of these formulae on the structure of the covariance matrix are illustrated below with the help of the Rayleigh-Levy flight model. Figure 3 compares the results from exact derivations of the covariance matrix with these prescriptions. The diagonal parts of the covariance matrices are well accounted for by these formulae. The root mean square of the measured local density PDF in particular exhibits the expected density dependence, at least for mild values of the density.

In all the formulae (19), (29), and (34), the expression of the covariance exhibits a simple structure, as it is factorizable in the two densities. This implies, for instance, that the reduced covariance matrix

(37)

has an extremely simple structure: it is given by the sign of the product of the bias functions (i.e., sign[b(ρi)b(ρj)], sign[bs1(ρi)bs1(ρj)], and sign[bs2(ρi)bs2(ρj)] for the three different measurement strategies). This leads to the butterfly-like structure in the plotted matrices, as illustrated in Fig. 4. This simple form betrays the fact that the density covariance is only poorly known. To be more specific, formulae (19), (29), and (34) give only a single eigendirection of the covariance matrix (namely b(ρi)P(ρi)) and the amplitude of a single eigenvalue associated with it. The numerical calculations suggest that it is the leading one when does not vanish, as illustrated on Fig. 7. These formulae do not offer any indication of the amplitude of the covariance in orthogonal directions, however. Taken at face value, they imply that the other eigenvalues all vanish, preventing the covariance matrix from being invertible. These formulae therefore cannot be used alone to model the covariance for practical purposes, and complementary contributions have to be derived from other (and a priori subdominant) effects.

3.5. Beyond leading-order effects

In the previous subsection, we identified the long-distance leading contributions. As mentioned before, this leads to only limited information of the covariance structure. This difficulty is even more acute for covariances evaluated in numerical experiments consisting of a collection of independent samples, each of them with periodic boundary conditions (this does not have to be so, but it is often the case in practice). By construction, the mean correlation function within the sample then vanishes, making the term we have computed identically zero. All these considerations indicate that further contributions need to be identified. The identification of the next-to-leading order effects in Eq. (16) is difficult to do a priori, however:

  • One natural next-to-leading contribution is obtained by taking into account second-order terms in ξ(d) in Eq. (16), that is, by considering doubles lines between cells in a diagrammatic representation. This would induce a term of about ξ(d)2, whose average never vanishes5. As shown in the appendix, these contributions can be formally derived in the context of the hierarchical models. This leads to correction terms that can be organized in a sum of factorized terms. Therefore, although it can indeed provide corrective terms to the covariance matrix, only a limited number of eigendirections can be generated.

  • Other contributions naturally come from proximity effects due to the fact that cells are finite, and could even overlap, which makes the expansion in ξ(d) ineffective. In a diagrammatic point of view, they are due to the fact that many more diagrams contributed when cells are too close. This has dramatic effects for overlapping cells. For hierarchical models, an approximate form can be used to help model these effects, which we use below.

  • Finally, effects due to the fact that discrete tracers are used in count-in-cells statistics might also play a role at short distances. They are also tentatively modeled below.

In the following, we propose some modeling of these effects and explore how they depend on the properties of the survey.

3.5.1. Joint PDF at short distances

There are no general forms for the joint PDF at close distance. The hierarchical models suggest the following form (derived from the saddle point approximation, which is valid for moderate values of and of the density contrast), however:

(38)

where ρm = (ρi + ρj)/2 and δρ = (ρi − ρj)/2, and where α is model-dependent parameter. In other words, the PDF of the difference between ρi and ρj can be described by a simple Gaussian with a known width driven by the expression of provided it is small compared to . We note that Δξ(d) obviously vanishes at d = 0, it then leads to a Kronecker δ function at zero separation as expected, and generically scales like d2 at short distances6. Interestingly, for the minimal tree model the form (38) is exact for α = 1 (see appendix). In general, this is also the expected form based on the saddle point approximation (valid when is small) for generic hierarchical tree models. The value of α can be identified from small-order cumulants,

(39)

where S3 is the reduced third-order cumulant,

(40)

This form is probably not very accurate in general. It can be used to model the impact of close distances to the covariance matrix, however, as shown below.

3.5.2. Poisson noise and minimum separation

A further contribution to this joint PDF can come from discrete effects that arise because the density is evaluated from the counting of discrete tracers (as explored in Szapudi & Colombi 1996 or more recently in Repp & Szapudi 2021). In this subsection, we assume that the density corresponds to the density obtained after application of a top-hat filter and that tracers are Poisson realizations of continuous fields (although it is possible to encounter sub- or over-Poissonnian noises; Friedrich et al. 2018). The use of other filters can be explored but would require specific developments that we do not pursue here. Within such hypotheses then, the joint distribution of counts-in-cells Ni is given by the convolution of the joint density PDF, ℘({ρi}), in the continuous limit convolved by Poisson counts-in-cells as

(41)

where PPoisson(N;) is more precisely the probability of having N tracers in a cell where the mean density of tracers is . In practice, i is given by nVi, where n is the number density of tracers and Vi is the volume of the cell Vi.

Discrete effects would then induce further scatter between the estimated values of ρi and ρj. The latter are given by Poisson noise induced by the nonoverlapping parts of the cells, as shown in Szapudi & Colombi (1996), further contributing to the scatter. The scatter in the difference in the number of points is

(42)

It can be incorporated as a contribution to the variance of the PDF of δρ of the form

(43)

where fe(d) is the fraction of the volume of each cell that does not overlap with the other as a function of the cell distance. For short distances (i.e., for about d ≲ R), it is in the 2D case given by

(44)

The expression (43) is then a priori to be added to the variance term that appears in Eq. (38) so that the total variance for the density difference reads

(45)

Equation (38) then fully encodes the fact that nearby cells are likely to have similar densities. This encodes, for instance, that nearby cells are within the same haloes. This contribution is expected to enhance the covariance terms. It shows that the amount of information is limited at small scales: there is therefore a minimum separation between cells smaller than which no gain in precision is expected of PDF measurements. The minimum distance depends on the bin size: dmin is the distance such that the densities in two cells separated by less than dmin are almost certainly in the sale density bin. dmin therefore depends on the bin width Δi. From Eq. (45), it is possible to infer this value. We desire

(46)

This suggests that a minimum distance between cells can be adopted, given by

(47)

The other upper limit comes on d from the expression of Δξ as a function of d. The latter depends on both the shape of the power spectrum and on the filter that is used. In general (e.g., for a Gaussian filter), Δξ(d) scales like d2/R2, where R is the filtering radius, with a coefficient cns that depends on the power spectrum index ns and is proportional to its amplitude. Top-hat filters have different analytical properties. We give here the formal expression of Δξ(d) at 2D for a power-law spectrum of index ns,

(48)

This is the situation we encounter below in the numerical experiments we perform. This leads to the following form:

(49)

It is to be noted that it can be in practice a rather short distance, shorter than the filtering scale R. For instance, for a bin width of 1/4, a variance of about unity, dminhalo is about R/5.

Equation (38), together with the expressions of the bias functions described above, is the main results of this paper. We illustrate below how they can be used to compute the covariance matrices.

4. Hierarchical models

In order to illustrate the previous findings, we make use of toy models for which explicit computations can be made.

4.1. General formalism

Hierarchical models are a class of non-Gaussian fields whose correlation functions follow the so-called hierarchical ansatz,

(50)

where the sum is made over all possible trees that join the p points (diagram without loops), and the tree value is obtained by the product of a fixed weight (that depends only on the tree topology) and the product of the two-point correlation functions for all pairs that are connected together in the given tree. This construction ensures that the average p-point function, , scales like the , where is the averaged two-point function. More precisely, there are Sp parameters such that

(51)

The precise value of the Sp parameters depend on the Qp parameters and on the averages of the product of ξ(rij) functions. A very good approximation is to assume that the average of the products of this function is given by the product of these averages. Then the Sp coefficients depend solely on Qp,

(52)

4.2. The (minimal) tree model

The tree models are based on a further assumption on the Qp parameters. It is basically assumed that tree expressions can be computed locally7, that is,

(53)

where νp is a weight attributed to all vertices with p incoming lines (ν0 = ν1 = 1 for completion). In this formalism, the vertex generating function is generally introduced,

(54)

The minimal tree model is a model in which ν2 alone does not vanish. In the minimal model8, its value is fixed and is given by ν2 = 1/2, so that

(55)

Together with the Gaussian case (which corresponds to ζ(τ)=1 + τ), this is the only case for which we are sure that it can be effectively built (in the sense that other models may be unphysical).

In this model, it is possible to build the cumulant generating function of the local density. For the one-point case, assuming the mean-field approximation, the CGF is given by

(56)

with

(57)

This is not the result of large deviation principle calculations, but of mere combinatorics, although it leads to the same formal transformation between the CGF and the vertex-generating function. In case of the minimal model, Eq. (57) takes a simple form that can be easily solved. We finally have

(58)

The one-point PDF of the density can then be computed explicitly (see appendix),

(59)

as can the density-bias function,

(60)

where is the averaged two-point correlation function within the cell.

4.3. Rayleigh-Levy flight model

The minimal tree model can be implemented with Rayleigh-Levy random walks (or rather Rayleigh-Levy flights, as described in Peebles 1980). This is a Markov random walk where the PDF of the step length ℓ follows a power law,

(61)

with a regularizing cutoff at small separation, and where α satisfies

(62)

The sample points are then all the step points reached by the walker.

More precisely, the cumulative distribution function of step of length ℓ is

(63)

(64)

where ℓ0 is a small-scale parameter. The two- and higher-order correlation functions can then be explicitly computed. Starting with a first point at position r0, the density of the subsequent point (first descendant) at position r is given by

(65)

(66)

In the following, the dimension of space is denoted D. The density of the descendants, assuming there are an infinity of them, of a point at position r0 is then given by a series of convolutions,

(67)

with subsequent convolutions and where the integral is done in the whole space. Defining the Fourier transform of f1(r) as ψ(k),

(68)

which is then a function of k only, it is easy to see that

(69)

where we take advantage of the expression of convolutions in Fourier space and their resummations. The two-point correlation function is then given by two possible configurations: a neighbor can either be an ascendant or a descendant, so that the two-point correlation functions between positions r1 and r2 are given by

(70)

where n is the number density of points in the sample that can be associated with a typical length ℓn,

(71)

At large scale, this expression causes the power spectra to be power laws. They scale like kα, and the resulting two-point correlation function then takes the form in the large separation limit,

(72)

(73)

It is to be noted, however, that this expression does not take into account the boundary conditions, in particular if they are assumed to be periodic. This case is examined in some detail in the next paragraph. It is to be noted, however, that in this case, the function ξ(r) has a more complex form. It is in particular no more isotropic.

Higher-order correlation functions can also be computed in this model: n points are correlated when they are embedded in a chronological sequence that can be run in one direction or the other. Thus the three-point function is simply given by

(74)

with five other terms obtained by all combinations of the indices. Expressing the result in terms of the two-point function, we have

(75)

corresponding to a tree structure with ν2 = 1/2.

Higher-order correlation functions can be computed similarly. They follow a tree structure in the sense above, with ν2 = 1/2 and νp = 0 for p ≥ 3.

4.4. Periodic boundary conditions

We briefly explore here the case of periodic boundary conditions. Then the multipoint density field gPBC(ri) for periodic boundary conditions can be expressed in terms of the former density field g(ri) as sums of all copies of the sample, that is,

(76)

where ni are vectors whose components are integers, and the sums run over all integer values for all i; L is the size of the sample (assumed to be the same in all directions).

When it is applied in this context, we can construct the n-point density function out of the density function f computed previously. Thus the two-point density function is given by

(77)

where n12 = n2 − n1 and nPBC is the resulting one-point (and therefore homogeneous) density in the sample. This expression is therefore written in terms of the function

(78)

We can now compute its expression in terms of the power spectra, or more specifically, the function ψ(k) defined previously. We have

(79)

and the latter sum ensures that the contributing wave modes k are only those that are periodic in all three directions, that is, those whose components are multiples of 2π/L so that

(80)

with

(81)

and where the sum is all over possible integer triplets but n = (0, 0, 0). The two-point correlation function is now given by

(82)

A similar result can be obtained for the three-point correlation function with

(83)

As a consequence, the functional relation between the three-point correlation function and the two-point correlation function is left unchanged. This is also the case at all orders.

4.5. Covariance matrix of the minimal tree model

Remarkably, in case of the minimal tree model, the derivation of the CGF can also be made for multiple cells, and in particular, for two cells. Its expression is derived in the appendix. It takes the form

(84)

In this case, it is then possible to expand its expression in powers of ξ12 for distant cells or in powers of for close cells, and in both cases, closed forms can be obtained to any order. It leads to the possibility of computing the joint PDF for any configuration (see the appendix for details) and finally to evaluate the covariance matrix directly. This is even possible for any of the thre sets of variables we consider, {ρi}, , or .

We performed these computations for the minimal tree model with a power-law behavior ξ(r)∼r−1.5 (α = 0.5), a 2D survey with a size of 2002 pixels, and a top-hat smoothing radius of 4.25 pixels. The amplitude of the correlation function was fixed to give at the smoothing scale. It precisely corresponds to the setting of the numerical simulations of Rayleigh-Levy flights we also performed, as described in the next section. It allows us to compare the two approaches. These analytic results have two limitations: the results are based on the mean-field approximation for the computation of the two-variable GCF, and the covariance elements are computed ignoring the bin sizes (i.e., by evaluating the expression of the covariance for their central values). Although in most cases this should not be an issue, it still might have a non-negligible impact when the PDF varies rapidly with the density.

5. Simplified models of the covariance matrix

The purpose of this section is then to propose two levels of modeling of the covariance matrix based on the previous results and to compare these propositions with results of either the full analytic results presented before or with the results of numerical experiments based on Rayleigh-Levy flights.

5.1. Modeling the covariance matrix

More specifically, we considered two approximate forms for the full covariance. The first approximation is fully analytic. It makes use of the large-scale contributions and those from the short distance expression (38). It reads as the sum of the two contributions

(85)

In this expression, the only free parameter is rmax. This is indeed a crucial parameter as it determines to a large extent the amplitude of the short-distance effects. In the following, we take rmax = R, that is, the filtering scale. It is found to give a good result for the 2D case and for ns = −1/2, but this choice is likely to depend on the shape of the power spectrum. In general, this formula is intended to give a good account of the general properties of the covariance matrix, it cannot provide reliable quantitative results a priori.

The other form we propose is intended to be much more precise quantitatively. Is is given by the following expression:

(86)

where CovPBC(ρi, ρj) is the expression of the covariant matrix for periodic boundary conditions. It is obtained here simply by replacing by before integrating over rd so that the averaged joint correlations vanish identically. The rationale for this proposition is that CovPBC(ρi, ρj) could be more easily estimated from specific numerical experiments. In both cases, the short-distance contributions are the same for the three types of observables ρi, , and . These forms are then compared to numerical results.

5.2. Numerical experiments with the Rayleigh-Levy flight model

A series of experiments of 2D walks with a large number of samples were performed as described below. We restricted our analysis to α = 0.5 with l0 = 0.003 pixel size (the dependence on l0 was tested as illustrated on Fig. 2, where l0 = .006 was also used, but the analyses were made for a fixed value of l0). Figure 1 illustrates how points are distributed in these samples. The point distribution does not show the filamentary structure of realistic cosmological simulations. It exhibits the presence of concentrated halos surrounded by empty regions, however, which are reminiscent of the structure of the largest matter concentrations of the cosmic web.

thumbnail Fig. 1.

Example of a realization of a Rayleigh-Levy walk. Points mark the end point of each displacement. They are clearly correlated.

Two different setting were employed to explore different aspects of the results that were found:

  • Set 𝒜: 1600 samples extracted from a single numerical realization (with periodic boundary conditions) with a size of 8000 × 8000 pixels2 containing 64 × 106 points. Each sample then has 200 × 200 pixels2 containing an average of 2002 points each. For this set of samples, the average and covariance of the PDF were extracted following the three procedures mentioned before: either the density was taken with respect to the mean density of the realization, with respect to the density of each sample, or by subtracting the sample density. It therefore corresponds to an evaluation of the mean and covariance of the PDF of ρi, , and , respectively.

  • Set 𝓑: 1600 samples, each with periodic boundary conditions, with a size of 200 × 200 pixels2 containing 2002 points each. By construction, the average two-point function in the sample, , vanishes in this case, and covariance is entirely due to proximity effects.

In each case, the local density was obtained after a filtering procedure. The point positions were first pixelized, that is, each point was attributed to a pixel so that the mean number of points per pixel was one. The field was then filtered by a (quasi) circular top-hat functions. In practice, the number of pixels in the window function was 57. This makes the effective smoothing radius about 4.25 in pixel units. The resulting density was then measured at each pixel location. Their histograms were then computed after density binning. To avoid large undue discrete effects, the bin width was chosen to be a multiple of 1/57, and in order to ensure that the requirement (49) was met at the pixel distance, we chose a bin width of about 1/4, more precisely, of 14/57.

Figure 2 shows the resulting PDF as measured in the simulations and how it compares to the theoretical prediction, Eq. (59), for two different choices of l0. The expected scaling for is recovered. The measured PDF also follows theoretical predictions for a wide range of probabilities remarkably well. It gives us confidence in the whole procedure and in the approaches used to compute PDFs in this model. The detailed comparisons were made for l0 = 0.003, leading to and a sample density variance in sets 𝒜 given by .

thumbnail Fig. 2.

One-point density PDF obtained with top-hat filters compared with the theoretical predictions, Eq. (59). The values of are 0.8 and 1.09 for the blue and red curves, respectively, corresponding to two different values of l0. The bottom panel shows the residuals. Departure from theory might be due to binning and/or to the finite number of samples.

The measured variance of the density PDF is obtained from 1600 samples in each case. The resulting shapes are presented in Figs. 35 for the different cases, density in a supersample realization, and in samples with periodic boundary conditions. The results show the comparison between results obtained in the numerical simulations with yellow symbols, and results derived from the analytic prescriptions as blue dots, based on the mean-field approximation. The agreement between the two is very good. The overall shape of the variance and its dependence on the density is well reproduced. Discrepancies can be observed for densities above 4 or 5, however, where the theoretical predictions are seen to underestimate the results. The reasons for these discrepancies are not clear at this stage. A possible explanation might be the finite number of samples that is used to infer the variances9. The variance of the density PDF is also compared with the large-scale contributions (19), (29), and (34) for set 𝒜 depending on the cases (at this order, the covariance vanishes for set 𝒝). It shows that this formula captures some features of the variance (especially at low and moderate densities), but does not account for all. This is also illustrated in Fig. 4, which shows the reduced covariance. The fact that the covariance is determined to a large extent by its leading large-scale contribution leads to values of the reduced covariance close to 1 or −1, leading to these butterfly patterns. Proximity effects, not captured in these forms, also contribute to the covariances at a significant level, however. This is already apparent in Fig. 3.

thumbnail Fig. 3.

Measured variance of the density PDF, i.e., diagonal elements of the covariance matrix, in sets 𝒜 for α = 0.5 and different prescription of the measured density. From left to right, raw density ρi, scaled density , and scaled density . The blue dots and solid lines are from the mean-field analytical expressions, and the large gold symbols are from the numerical simulations. The dashed black lines are what is expected from the large-scale leading contribution. The variance at cell scale is about 1.09, and the variance at sample scale, , is about 0.09.

thumbnail Fig. 4.

Resulting reduced covariance matrix for the three types of observables for set 𝒜. The covariance matrix is dominated by its leading eigenvalue and direction, leading to this typical butterfly shape of the reduced covariance matrix.

thumbnail Fig. 5.

Measured variance of the density PDF obtained for set 𝒝). Symbols are the same as in Fig. 3.

5.3. Testing models of covariance matrices

Expressions (85) and (86) are precise propositions to show how the large-scale contributions can be completed to account for the full form of the covariance. The comparisons between the predicted form and those obtained from the numerical experiments are explored in detail at different levels and using the following criteria:

  • Amplitude of the PDF variance,

  • Density dependence of the first eigenvalue of the covariance matrix,

  • Amplitude of the eigenvalues of the covariance matrix, and

  • Resulting χ2 distribution of a set of data vectors drawn from the original covariance.

These comparisons are shown in Figs. 68. For model (86), the term CovPBC(ρi, ρj) is taken from the measured covariance of set 𝒝. Figures 6 and 7 show that these two prescriptions give a good account of the leading behavior of the covariance matrix. The conclusion is quite sensitive for the choice of rmax for prescription (85). On the other hand, there is no free parameter that can be adjusted for model (86). Interestingly, Fig. 7 shows that the PDF variance also departs significantly from the large-scale term. The first eigenvector reproduces the functional form of the large-scale density-bias functions very faithfully.

thumbnail Fig. 6.

Measured variance of the density PDF, i.e., diagonal elements of the covariance matrix, in sets 𝒜 and comparisons with proposed approximate forms. The yellow line and symbols are the results obtained in the numerical experiments. The dot-dashed line is the prediction derived from relation (86), and the dashed gray line shows the prediction from Eq. (85). The dot-dashed black lines correspond to the large-scale contributions.

thumbnail Fig. 7.

Behavior of the first eigenvector with the same color-coding as in Fig. 6. The dashed black lines are the large-scale prediction, b#(ρi)P(ρi) appropriately normalized. The size of the data vector is 30.

thumbnail Fig. 8.

Performances of the approximate forms of the covariance matrix in terms of rigenvalues and χ2-distributions. Top panel: eigenvalues of the covariance matrices (rebinned into six bins) compared to what can be obtained from the proposed approximate forms; same color-coding as for Fig. 6. The χ2 distributions are shown in the bottom panel. Model (86) reproduces the very same χ2 distributions. Model (85), in gray, is not as accurate and tends to slightly overestimate the χ2. This latter behavior is amplified when a larger number of bins is used.

The last two criteria are designed to verify that the reconstructed covariances also capture the subleading behavior of the matrix and can eventually be safely inverted and used as a model of likelihood. To avoid numerical uncertainties and make the comparison tractable, we chose to reduce the binning to six bins (through a rebinning of the histograms and densities ranging from 0.5 to 6.5). The resulting eigenvalues are shown in the top panel of Fig. 8. It shows that the eigenvalues decrease rapidly in amplitude, suggesting that the eigendirections are well sequenced and that the approximate form captures their values rather accurately. Form (86) in particular reproduces all six eigenvalues almost exactly.

Finally, χ2 distributions were computed from a set of random values drawn in each case from a Gaussian likelihood built from the measured covariance (with six bins). The values of were then computed for each data vector, and their histogram was computed from each of the proposed models (including the original model for reference),

(87)

where is the inverse of the covariance matrix, either computed from Eq. (85) or from Eq. (86). For the original model, the expected distribution of the χ2 values is then expected to be precisely that of a χ2 distribution with six degrees of freedom. This is indeed what is almost exactly obtained for model (86). Results obtained from prescription (85) are not quite as good. This is expected as the short-distance effects are estimated rather crudely in Eq. (85). The performance of this prescription deteriorates when the dimension of the data vector (i.e., the number of bins) increases.

6. Conclusions and lessons

We presented key relations that give the large-scale behavior of the joint PDF, and hence the leading behavior of the covariance matrix of the density PDF. These contributing terms do not give the expression of a covariance matrix that can be used to build a likelihood function, however, as it is not invertible. Further significant contributions are found to be due to small separation effects, and an approximate form is proposed in Eq. (38). The latter is found to encapsulate most of the proximity effects, that is, it informs on the fact that nearby regions are likely to be correlated. They also give an indication on the minimal grid size that can be used the maximum bin size that can be used without information loss for a given bin width.

We then used a toy model for which numerical experiments can easily be performed and for which the exact PDF and large-scale covariance can be derived. It allows us to evaluate the efficiency of approximate schemes precisely. The conclusions of these comparisons are listed below.

  • The theoretical forms Eqs. (19), (29), and (34) give the leading-order expression of the covariance elements when supersample effects are taken into account. It gives an accurate prediction of the leading eigenvalue and eigendirection of the covariance matrix.

  • Whether subdominant effects can be accounted for by subsequent terms depends on the behavior of the two-point function: if the rms of the two-point function is dominated by large separations, then next-to-leading-order effect need to be taken into account; otherwise, short-distance effects will be the dominant contributor.

  • In case short-distance effects dominate, the covariance matrix can be accessed from small simulations provided the relevant dominant large-scale contributions are added.

  • This suggests that in realistic situations, the supersample effects, that is, the effects due to modes whose wavelength is larger than the size of the survey, have limited impact on the structure of the covariance matrix and that they can be captured by the only leading large-scale contribution. This is supported by a further analysis of the behavior of ξ(rd) in realistic cosmological settings. For the standard model of cosmology (as derived from cosmic microwave background observations, Planck Collaboration VI 2020), the behavior of the matter correlation function can be derived. This is illustrated in Fig. 9, which illustrates the scales that are the main contributors to the first two moments of the two-point correlation function. Whether in 2D or in 3D, the first moment is dominated by large-scale contributions, whereas the second moment is dominated by small-scale contributions.

    thumbnail Fig. 9.

    Scale dependence of the matter correlation functions for a realistic cosmological model (cosmological parameters derived from Plank, Planck Collaboration VI 2020) for the 3D density and the projected density (for a uniformly sampled survey with a depth of about 800 h−1Mpc between z = 0.75 and z = 1.25). The top panel shows (solid blue line) and (dashed red line) for the 3D density field, and the bottom panel shows and for the projected density. In both cases, the average value of the first moment of the two-point correlation function is dominated by large-distance contributions, whereas short-distance contributions dominate the second moment, assuming survey sizes of about 100 h−1 Mpc or above.

  • In the context of this study, we assumed that the measured Pi were Gaussian distributed. Although it is difficult to assess the accuracy of this hypothesis, the structure uncovered in Sect. 3 can be used to make such an attempt. In tree models, higher-order expressions of the joint density PDFs are expected to preserve the tree structure; see Bernardeau & Schaeffer (1999). The connected part of the three-point joint density PDF is then expected to take the form

    (88)

    where b2(ρ) is the two-line bias function of amplitude similar to b2(ρ). This implies in particular that the third-order cumulant is about , much smaller than , making the distribution of the measured values of P(ρ) (quasi-) Gaussian distributed. There might be some combination of ρi and values of , however, for which a higher-order term could play a role in the expression of the likelihood function.

For the application of these formulae in practical cases, some limitations have to be noted. We list them below.

  • In the proposed form, the fact that in practice, PDFs are generally measured on a grid, that is, on a finite set of locations, is not taken into account. For instance, the exclusion of nonoverlapping cells is not considered. this is expected to introduce additional noise in the PDF estimates. The covariance matrix for these constructions cannot then be derived from general formulae (6) even when the integral in rd is restricted above a given threshold.

  • Relation (38) has been derived in a specific regime (using saddle point approximations) for tree hierarchical models. They are expected to capture the phenomenon at play for “typical” values of the densities, but they may not perform so well in the rare event tails (the exception being the minimum model, for which it is exact). Further checks of the validity of (38) should therefore certainly be done.

  • The general formulae (19), (29), and (34) are valid for any type of filtering schemes, even for a compensated filter. This is not the case for relation (38). The proximity effects for compensated filters ought to be considered specifically.

  • Prescription (86) is found to give a very precise account of the properties of the covariance matrix. It is based on the proposition that large-scale (supersample) effects can be added separately from the proximity effects and that the latter can be evaluated with small-scale mocks in which supersample effects are absent (with periodic boundary conditions). This is not an exact result, however,. It relies in particular on the fact that the r.m.s. of the ξs is dominated by scales much smaller than the sample size.

  • Prescription (85) is less solid. It can be used for a quick assessment of the different contributing terms, or to build fully invertible covariance matrices, but it is unlikely to give reliable predictions at the χ2 level.

In all cases, prescriptions (85) and (86) can be the starting point of a more precise evaluation of the covariance from specific numerical experiments that can complement its evaluation following the approach presented in Friedrich & Eifler (2018), for instance. The authors also showed that some strategies could be adopted to limit the number of realizations required to reach a specific precision. This point is not discussed here.


1

This is an ideal estimate in the sense that μ is evaluated in an infinite number of locations. We therefore neglect here the impact of measuring μ on a finite number of locations when evaluating Pi. Regarding this aspect, a specific derivation that takes a finite number of measurements into account can be found in Codis et al. (2016a).

2

Whether this is a correct assumption is difficult to assess in general. It probably depends on the detailed properties of the setting. The Conclusion section contains further comments on this aspect.

3

For a square survey of unit size (with nonperiodic boundary conditions), the distance distribution function Ps(rd) is given by

as can be obtained after integrating over three of four of the position coordinates.

4

This is not necessarily so, as exemplified in Carron (2011), Carron & Neyrinck (2012).

5

In the minimal tree model, it is possible to compute these terms in the so-called mean-field approximation (see appendix), but they do not lead to a positive definite covariance matrix and therefore cannot be the sole, or dominant, contribution to the covariances.

6

This limited form would induce a minimum contribution to Cov(ρi, ρi) given by Δρi/℘(ρi), where Δρi is the bin size in density.

7

Perturbation theory results do not exactly follow this construction as vertices are then dependent on the geometry of the incoming lines. However, in this case, Qp values are indeed obtained from a product of vertices.

8

it is minimal in the sense that it can be shown that ν2 cannot be smaller than 1/2 in the strongly nonlinear regime (Peebles 1980).

9

Although the number of samples is large, the number of haloes contained in each sample is finite leading to discretization errors in the estimate of the covariance. Estimate of the minimal number of realizations required to make such estimates is beyond the scope of this paper.

10

The formulae are written here for a top-hat profile, but can be extended to arbitrary profiles.

11

It takes the very same form as the stationary equation in the context of the large deviation principle, although the results here do not correspond to this regime.

12

A quicker approach is to view φ0(λ) as a function of the “leaf weight”, that is, the value of ν1; φ1(λ), φ2(λ) are then obtained by successive derivatives of φ0(λ) with respect to ν1.

13

it is possible to fully solve the consistency relations in case of two cells; the conclusions remain unchanged.

Acknowledgments

The author of this article is indebted to Cora Uhlemann, Alex Gough, Oliver Friedrich, Sandrine Codis, Aoife Boyle and Alexandre Barthelemy for many comments and careful examination of the preparatory notes of this manuscript.

References

  1. Abbott, T. M. C., Abdalla, F. B., Alarcon, A., et al. 2018, Phys. Rev. D, 98, 043526 [NASA ADS] [CrossRef] [Google Scholar]
  2. Amendola, L., Appleby, S., Avgoustidis, A., et al. 2018, Liv. Rev. Rel., 21, 2 [Google Scholar]
  3. Balian, R., & Schaeffer, R. 1989, A&A, 220, 1 [NASA ADS] [Google Scholar]
  4. Barthelemy, A., Codis, S., & Bernardeau, F. 2021, MNRAS, 503, 5204 [Google Scholar]
  5. Bernardeau, F. 1992, ApJ, 392, 1 [NASA ADS] [CrossRef] [Google Scholar]
  6. Bernardeau, F. 1996, A&A, 312, 11 [NASA ADS] [Google Scholar]
  7. Bernardeau, F. 2013, ArXiv e-prints [arXiv:1311.2724] [Google Scholar]
  8. Bernardeau, F., & Reimberg, P. 2016, Phys. Rev. D, 94, 063520 [NASA ADS] [CrossRef] [Google Scholar]
  9. Bernardeau, F., & Schaeffer, R. 1992, A&A, 255, 1 [NASA ADS] [Google Scholar]
  10. Bernardeau, F., & Schaeffer, R. 1999, A&A, 349, 697 [NASA ADS] [Google Scholar]
  11. Bernardeau, F., & Valageas, P. 2000, A&A, 364, 1 [NASA ADS] [Google Scholar]
  12. Bernardeau, F., Pichon, C., & Codis, S. 2014, Phys. Rev. D, 90, 103519 [NASA ADS] [CrossRef] [Google Scholar]
  13. Boyle, A., Uhlemann, C., Friedrich, O., et al. 2021, MNRAS, 505, 2886 [NASA ADS] [CrossRef] [Google Scholar]
  14. Brouwer, M. M., Demchenko, V., Harnois-Déraps, J., et al. 2018, MNRAS, 481, 5189 [NASA ADS] [CrossRef] [Google Scholar]
  15. Carron, J. 2011, ApJ, 738, 86 [NASA ADS] [CrossRef] [Google Scholar]
  16. Carron, J., & Neyrinck, M. C. 2012, ApJ, 750, 28 [NASA ADS] [CrossRef] [Google Scholar]
  17. Codis, S., Pichon, C., Bernardeau, F., Uhlemann, C., & Prunet, S. 2016a, MNRAS, 460, 1549 [NASA ADS] [CrossRef] [Google Scholar]
  18. Codis, S., Bernardeau, F., & Pichon, C. 2016b, MNRAS, 460, 1598 [NASA ADS] [CrossRef] [Google Scholar]
  19. Colombi, S., Bouchet, F. R., & Schaeffer, R. 1995, ApJS, 96, 401 [NASA ADS] [CrossRef] [Google Scholar]
  20. Friedrich, O., & Eifler, T. 2018, MNRAS, 473, 4150 [NASA ADS] [CrossRef] [Google Scholar]
  21. Friedrich, O., Gruen, D., DeRose, J., et al. 2018, Phys. Rev. D, 98, 023508 [Google Scholar]
  22. Friedrich, O., Uhlemann, C., Villaescusa-Navarro, F., et al. 2020, MNRAS, 498, 464 [NASA ADS] [CrossRef] [Google Scholar]
  23. Fry, J. N. 1984a, ApJ, 277, L5 [NASA ADS] [CrossRef] [Google Scholar]
  24. Fry, J. N. 1984b, ApJ, 279, 499 [NASA ADS] [CrossRef] [Google Scholar]
  25. Gruen, D., Friedrich, O., Krause, E., et al. 2018, Phys. Rev. D, 98 [Google Scholar]
  26. Heymans, C., Tröster, T., Asgari, M., et al. 2021, A&A, 646, A140 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  27. Hikage, C., Oguri, M., Hamana, T., et al. 2019, PASJ, 71, 43 [Google Scholar]
  28. Ivezić, Ž., Kahn, S. M., Tyson, J. A., et al. 2019, ApJ, 873, 111 [Google Scholar]
  29. Kaiser, N. 1998, ApJ, 498, 26 [NASA ADS] [CrossRef] [Google Scholar]
  30. Laureijs, R., Amiaux, J., Arduini, S., et al. 2011, ArXiv e-prints [arXiv:1110.3193] [Google Scholar]
  31. Mo, H. J., & White, S. D. M. 1996, MNRAS, 282, 347 [Google Scholar]
  32. Peebles, P. J. E. 1980, The Large-scale Structure of the Universe (Princeton University Press) [Google Scholar]
  33. Planck Collaboration VI. 2020, A&A, 641, A6 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  34. Repp, A., & Szapudi, I. 2020, MNRAS, 498, L125 [NASA ADS] [CrossRef] [Google Scholar]
  35. Repp, A., & Szapudi, I. 2021, MNRAS, 500, 3631 [Google Scholar]
  36. Schneider, P. 1996, MNRAS, 283, 837 [Google Scholar]
  37. Szapudi, I., & Colombi, S. 1996, ApJ, 470, 131 [NASA ADS] [CrossRef] [Google Scholar]
  38. Szapudi, I., Colombi, S., & Bernardeau, F. 1999, MNRAS, 310, 428 [NASA ADS] [CrossRef] [Google Scholar]
  39. Takada, M., & Hu, W. 2013, Phys. Rev. D, 87 [CrossRef] [Google Scholar]
  40. Touchette, H. 2011, ArXiv e-prints [arXiv:1106.4146] [Google Scholar]
  41. Uhlemann, C., Friedrich, O., Villaescusa-Navarro, F., Banerjee, A., & Codis, S. 2020, MNRAS, 495, 4006 [NASA ADS] [CrossRef] [Google Scholar]
  42. Valageas, P. 2002, A&A, 382, 412 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  43. White, S. D. M. 1979, MNRAS, 186, 145 [NASA ADS] [CrossRef] [Google Scholar]

Appendix A: Hierarchical tree models

In hierarchical tree models, the p-point matter correlation functions are assumed to follow tree structures in the sense described in the main text. They are thus entirely defined by the two-point functions ξ(r) and the vertex-generating function ζ(τ). The exact generating function of multiple cell correlation functions can be built through simple transforms. We therefore consider a set of n cells Vi. These cells can overlap. The joint cumulants we consider are those of the average densities in cells Vi that can be expressed in terms of spatial averages10 of correlation functions,

(A.1)

We then wish to build the cumulant-generating function,

(A.2)

This function represents the generating function of (averaged) tree diagrams where λi counts the number of points in each cells. As shown in Bernardeau & Schaeffer (1999), this is obtained with the help of the intermediate function τ(x) solution of the consistency equation11,

(A.3)

and then

(A.4)

This is an exact result based on pure combinatorics.

For cases of interest, it is possible to do a mean-field approximation that consists of assuming that τ(x) is constant within each cell. We then have the system of equations for τi,

(A.5)

where

(A.6)

and

(A.7)

Bernardeau & Schaeffer (1999) found this to be very accurate, and we extensively use this approximation in the following, in particular for the minimal tree model.

Appendix B: Joint PDF, density-bias function in PT, and hierarchical tree models

Here we consider the joint distribution of densities in two cells whose centers are at distance d. The calculation is based on the inverse Laplace transform of the joint cumulant-generating function φ(λ1, λ2),

(B.1)

where are the cumulants of the local density fields. They depend on the size and distance d between the cells. We assume in particular that the density correlation function between cells ξ(d) is small compared to unity and can serve as a small parameter.

B.1. Leading order in the mean-field approximation

Expanding with respect to ξ(d) then leads to the following form:

(B.2)

that is, a factorization of the linear term in ξ(d). This can explicitly be shown in case of tree models (as described in the main text). This is also the case in perturbation theory.

In case of the tree models, we have

(B.3)

where is the averaged correlation function within the cells. This is obtained assuming a mean-field approximation. We expect subleading corrections when d becomes comparable to the size of the cells.

B.2. Extending the previous case to the three variable case

In addition to the two variables ρ1 and ρ2, we introduce here the variable ρs, which is the density within the sample. We assume that the correlation functions are built with the same model. Here the small parameter is the correlation between two cells V1 and V2 (at positions x1 and x2) and the correlation function at sample size. It is natural in the context we are interested in to assume that these two quantities are on the same order.

We further assume we can use the mean-field approximation for the two cells V1 and V2. This is not a priori the case for the density in the whole sample, however. We therefore derive the results without this approximation. The general expression is then

(B.4)

with the consistency relations

(B.5)

(B.6)

(B.7)

We therefore derive the expression of φ(λs, λ1, λ2) up to linear order jointly in , ξ(xs, x1), ξ(xs, x1) and ξ12 ≡ ξ(x1, x2). At zeroth order, τ(xs) vanishes and

(B.8)

with a similar relation for . At linear order, we have

(B.9)

(B.10)

(B.11)

where

(B.12)

The resulting cumulant-generating function reads

(B.13)

where the first three terms are at zeroth order and the last three are at linear order. Using the previous expression, we obtain

(B.14)

This relation is used to derive the expression of the sample bias functions in the next section.

B.3. Second order in the mean-field approximation

Results of Sect. B.1 can naturally be extended to any order in the cross-cell correlation function in the context of the tree-hierarchical models (as illustrated on Fig B.1). Up to second order, it takes the form

(B.15)

thumbnail Fig. B.1.

Diagrammatic visualization of the bias functions. The function φ0(λ) is the generating function of all trees within one cell, φ1(λ) of all trees within one cell with one external lines, and φ2(λ) with two external lines. The resulting connected diagrams up to second order in ξ12 are thus those presented here. Two φ2 generating functions cannot be conntected as that would induce a loop contribution. This reflects the underlying tree structure.

where the function φ2(λ) takes the form

(B.16)

This last expression can be directly obtained through a perturbative expansion as presented in the previous subsection12.

In case of the minimal tree model, this perturbative expansion can be directly compared with exact mean-field results13. This is shown in Fig. B.2. It shows that for a large regime in λ, the relation (B.15) provides a very accurate description of the joint cumulant-generating function down to a distance corresponding to overlapping cells. For overlapping cells, relation (B.15) continues to be accurate except for high values of λ. In general,

thumbnail Fig. B.2.

Joint CGF as a function of d and for different values of λ1, λ2: from bottom to top, λ1 = λ2 = −0.1, λ1 = λ2 = 0.2, λ1 = λ2 = 0.5. The solid blue lines correspond to the two-cell mean-field expression, (C.13). The predictions given in Eqs. (B.2) and (B.15) are shown as dashed red lines and dotted black lines. The shaded area is the region of overlapping cells.

(B.17)

when d → 0 (more accurate results in case of the minimal tree model are given in Appendix D.).

Using (B.15), we then derive corrective perturbative terms to the joint density PDF. More specifically, we have

(B.18)

where

(B.19)

(B.20)

For a sample with periodic boundary conditions, the average of ξ12 vanishes, which a priori makes the other terms the leading contributors to the covariance elements. Equation (B.18) can be written as a sum of symmetric factorized terms,

(B.21)

showing the eigenstructure of the resulting matrix and showing that it defines three different eigendirections at most.

B.4. Relative density joint PDFs and bias functions

We wish to compute the joint PDF of the density when expressed in terms of the survey average density ρs. In order to do so, we consider the joint density P(ρs, ρi, ρj), where ρs is the density in the sample and ρi and ρj are the densities in two cells at distance d.

We wish to compute the joint probability distribution function of and , defined as

(B.22)

and the joint distribution functions of and , defined as

(B.23)

From these changes of variables, we have

(B.24)

and

(B.25)

Similarly, we also have

(B.26)

and

(B.27)

We continue the calculations by expressing the joint PDF with the help of inverse Laplace transforms,

(B.28)

(B.29)

As a result,

(B.30)

(B.31)

after integration over ρs. The latter expressions can be expressed as

(B.32)

(B.33)

In a similar manner, we can obtain the form of the joint PDF for ,

(B.34)

(B.35)

We then use the relation (B.14) to compute the form of these functions.

Noting that the expressions , ∫dx0ξ(x0, x1) take all the same averaged value when integrated over the sample, which we note , then at linear order in ,

(B.36)

(B.37)

At the same order, we then have

(B.38)

Combining both the expressions of and and expanding all terms at linear order in , we obtain

(B.39)

This leads to the definition of the first sample bias function,

(B.40)

which can be re-expressed in terms of the density-bias function defined in Eq. (18) and the derivative of with respect to

(B.41)

The second sample-bias function can be obtained in a similar manner. We indeed have

(B.42)

which eventually leads to

(B.43)

B.5. Response to a change in amplitude in ξ

A close notion related to the density-bias function is how the PDF is changed when the parameters of the simulations are changed. In particular for tree models, the statistical properties are entirely determined by the amplitude of the two-point function, for instance, at cell size. This dependence can be made explicit by writing Eq. (17) as

(B.44)

after the change of variable and function,

(B.45)

where then the expression does not depend on (only on the functional form of ζ). It follows that

(B.46)

This expression can be used to defined the function bξ(ρi) as

(B.47)

It appears that bξ(ρi) is very similar to bs1(ρi), but the two are not equal in general.

B.6. Close cell results

B.6.1. CGF for 2 close cells

Saddle point approximation and close cell results. In case of two cells, the general system in the mean-field approximation leads to

(B.48)

(B.49)

and

(B.50)

We are interested here in the behavior of φ(λ1, λ2) when the two cells are close, that is, when . When , τ1 and τ2 are also equal, making φ(λ1, λ2) a sole function of λ1 + λ2 and therefore forcing the joint PDF to be proportional to δDirac(ρi − ρj). To be more precise, in this regime, , δρ ≡ (ρi − ρj)/2 is expected to be distributed with a width of about . This suggests that in this limit, the difference λ1 − λ2 should be treated as a large quantity of about . The limit behavior of the joint CGF can then be explicitly computed in terms of

(B.51)

(B.52)

(B.53)

In this limit, we obtain

(B.54)

(B.55)

leading to

(B.56)

and

(B.57)

The joint PDF of ρm = (ρi + ρj)/2 and δρ then reads

(B.58)

for which there is in general no closed form. We then need to rely on approximation schemes to complete the calculations.

B.6.2. Saddle point approximation

One of the approximations that can be used to evaluate Eq. (B.58) is to use the saddle point approximation. It has been used in the literature to compute the PDF (see Balian & Schaeffer 1989; Bernardeau 1992; Valageas 2002; Codis et al. 2016a). It is a priori valid when is small (and not for too high values of the density). In this approach, the expression under the exponential is approximated by a quadratic form at its minimum. In practice, the latter is obtained from the resolution of the system

(B.59)

(B.60)

which leads to the implicit or explicit values of λ, μ, and τ at the saddle point position (we hereafter denote this with the subscript s),

(B.61)

It is to be noted here that the value of τs is independent of δρ. At the saddle point position, we then have

(B.62)

This then suggests that the joint PDF is given by the product of the one-point PDF of ρm and a Gaussian distribution of δρ of width . For usual models, as described above, ζ′(τs)2 scales like a power of ρm so that one suggested form for the joint PDF is the following:

(B.63)

Interestingly, the value of α can be related to the reduced skewness of the density field from the computation of ⟨δρ2ρmc, and in the context of tree hierarchical models, it leads to

(B.64)

The validity of this form clearly ought to be checked. Its simplicity nonetheless offers a good grasp of the contribution of close cells to the covariance matrix.

Appendix C: Minimal tree model

In the previous section, general formulae were given. The aim of this section is to account for more precise results obtained in the case of a specific hierarchical model, namely the minimal tree model, as described below. It can then serve as a toy model for the construction of the approximate form for the covariance matrix. We first recall that this model describes the Rayleigh Levy flights model.

C.1. One-point results in the mean-field approximation

The Rayleigh-Levy flight model makes it possible to build synthetic samples whose statistical properties follow the minimal model, that is, a hierarchical model with the following vertex-generating function:

(C.1)

In the one-cell mean-field approximation, we have the equation

(C.2)

which can be solved in

(C.3)

which leads to the following expression for the CGF:

(C.4)

The one-point PDF of the density can then be computed explicitly. It takes the form in the continuous limit of

(C.5)

For this particular model, the void probability distribution (VPF) is nonzero even in the continuous limit. We recall here that the general expression of the VPF is given by exp(φ(−)), which for the minimal model leads to

(C.6)

when  → ∞.

The density-bias function can also be computed explicitly. For the minimal model, we have φ1(λ)=φ(λ) so that

(C.7)

for ρ >  0. For this model, the fact that φ1(λ)=φ(λ) implies that

(C.8)

This means that in the case of the minimal model, the density-bias function can be extracted from the functional form of the one-point PDF as

(C.9)

This is a somewhat remarkable identity (which can be extended to higher orders, as shown below.)

In this case, the second-order expansion leads to the form φ2(λ) given by

(C.10)

and we note that φp(λ) all vanish for p ≥ 3.

C.2. Two-cell results in the mean-field approximation

These mean-field calculations can be extended to the two-cell case. In this case, we have the system

(C.11)

(C.12)

when the two cells are of the same size. This leads to the following expression for the joint CGF:

(C.13)

Remarkably, this expression can be written in terms of the one-cell CGF as

(C.14)

This opens the possibility of computing the joint PDF to any order of ξ12. The calculation of this expansion is made simple by the following observations: The corrective terms will make intervene functions of the forms

(C.15)

We further note that

(C.16)

with

(C.17)

We the define cn(ρ) as

(C.18)

We then have on one side

(C.19)

and on the other side

(C.20)

which derives from the fact that

(C.21)

after integration by parts. As a result, the expression of the join PDF to any order can be written as polynomials making intervene P(ρ1), b(ρ1) P(ρ2) and b(ρ1) only.

C.3. Perturbative expansion for close cells

Another interesting result is when the cell centers are close (so that cells overlap), as described above. In this case, the limit behavior of the joint CGF is given by

(C.22)

with

(C.23)

(C.24)

(C.25)

It is then remarkable to see the result can be expressed with the sole one-cell CGF,

(C.26)

In other words, the GFC of the variables ρm = (ρ1 + ρ2)/2 and δρ = (ρ1 − ρ2)/2 is given by Eq. (C.26). It is possible to compute the joint PDF,

(C.27)

with the change of variable

(C.28)

The integral in leads to the one-cell PDF of the density ρm, whereas the integral in μ can be done explicitly as it is a quadratic form in μ, leading to a Gaussian distribution in δρ. The final PDF is given by

(C.29)

This shows that the joint PDF peaks for ρ1 ∼ ρ2 with a width that depends on the distance between the cells through the difference . Moreover, this form has no overlapping regime with the previous expansions of the joint PDF. It captures different aspects of the covariance calculations as listed below.

  • The previous expression says that close cells contribute more specifically to the covariance when ρ1 and ρ2 are close. This suggests that Eq. (C.29) contributes mostly to the near diagonal terms, whereas off diagonal terms could still be well described by perturbative expansions, as described before.

  • As noted before, perturbative expansions are closely related to supersample effects. They encode the way in which the local densities are jointly correlated with long-wavelength modes. This is not the case in Eq. (C.29). It rather captures how a rare event, such as a peak, can contribute to the covariance elements: if there is a peak somewhere, nearby cells are likely to have a similar density up to distances for which remains small enough.

The above development can be pursued to any order in provided the following recipe is applied:

(C.30)

Then the joint density can be computed to any order in , making use of the very same expressions bn(ρ).

The next-to-leading order in is thus given by

(C.31)

and the expansion can be extended in any (even) order in . Fig. C.1 illustrates the convergence properties of these expansions. Depending on , either the expansion in ξ12 or that in gives a very accurate estimate of the joint PDF. It opens the way to computing the covariance matrix starting in the two-cell mean-field approximation (C.13).

thumbnail Fig. C.1.

Performances of the perturbative expansions of the joint PDF in the mean-field approximation either for the ξ12 expansion (open blue dots) or the expansion (red dots) up to 11th and 10th order, respectively. The comparisons are made for ρ1 = ρ2 = 2 and (top panel) and for ρ1 = 0.5,  ρ2 = 3.5 and (bottom panel) and for ξ12 equalling 0.1, 0.3, 0.5, 0.7, and 0.9.

C.4. Construction of the theoretical covariance matrix for the minimal tree model

The previous form can be used to compute the covariance matrix for the minimal tree model in simple implementations. It relies on analytic forms for both the two-point cell correlation functions, which can formally be written as

(C.32)

for a given power spectrum. We also make use of the form Ps(rd) given in footnote 2 to derive the PDF of cell distances. We then have all the required ingredients to compute the elements of the covariance matrix in the mean-field approximation,

(C.33)

In practice, is computed from the eighth-order expansion either in ξ12 when or in Δξ when . This is used to explore the detailed properties of the covariance matrix and the validity of approximate schemes.

C.5. Joined PDF for relative densities

The minimal model allows us also to pursue the computation of the joint PDF for the variables or in all regimes. The first step is to extend Eq. (B.14) to a regime in which ξ12 is not assumed to be small. We find that

(C.34)

where φc(λ1, λ2) is given by

(C.35)

and we can note that

(C.36)

At leading order in ξs, that is, when we assume that the density fluctuations at sample size are much smaller than at smoothing scale, this expression then reduces to

(C.37)

We can then exploit this relation to compute the and from Eqs. (25) and (31), respectively. We then have at leading order in ξs

(C.38)

and

(C.39)

To complete the formal calculation of these expressions, we introduce the function

(C.40)

We can first note that Eq. (C.9) can be extended to

(C.41)

This comes from the observation that

can also be written

(C.42)

and that

(C.43)

The final expression of the PDF of the relative densities can then be obtained by noting that applying a multiplicative factor λi to the moment-generating function is equivalent to the application of the operator ∂/∂ρi to the final expression this finally leads to the following forms:

(C.44)

(C.45)

for the and

(C.46)

(C.47)

for . These relations can then be applied to the expressions of the joint density such as found in the previous subsection.

All Figures

thumbnail Fig. 1.

Example of a realization of a Rayleigh-Levy walk. Points mark the end point of each displacement. They are clearly correlated.

In the text
thumbnail Fig. 2.

One-point density PDF obtained with top-hat filters compared with the theoretical predictions, Eq. (59). The values of are 0.8 and 1.09 for the blue and red curves, respectively, corresponding to two different values of l0. The bottom panel shows the residuals. Departure from theory might be due to binning and/or to the finite number of samples.

In the text
thumbnail Fig. 3.

Measured variance of the density PDF, i.e., diagonal elements of the covariance matrix, in sets 𝒜 for α = 0.5 and different prescription of the measured density. From left to right, raw density ρi, scaled density , and scaled density . The blue dots and solid lines are from the mean-field analytical expressions, and the large gold symbols are from the numerical simulations. The dashed black lines are what is expected from the large-scale leading contribution. The variance at cell scale is about 1.09, and the variance at sample scale, , is about 0.09.

In the text
thumbnail Fig. 4.

Resulting reduced covariance matrix for the three types of observables for set 𝒜. The covariance matrix is dominated by its leading eigenvalue and direction, leading to this typical butterfly shape of the reduced covariance matrix.

In the text
thumbnail Fig. 5.

Measured variance of the density PDF obtained for set 𝒝). Symbols are the same as in Fig. 3.

In the text
thumbnail Fig. 6.

Measured variance of the density PDF, i.e., diagonal elements of the covariance matrix, in sets 𝒜 and comparisons with proposed approximate forms. The yellow line and symbols are the results obtained in the numerical experiments. The dot-dashed line is the prediction derived from relation (86), and the dashed gray line shows the prediction from Eq. (85). The dot-dashed black lines correspond to the large-scale contributions.

In the text
thumbnail Fig. 7.

Behavior of the first eigenvector with the same color-coding as in Fig. 6. The dashed black lines are the large-scale prediction, b#(ρi)P(ρi) appropriately normalized. The size of the data vector is 30.

In the text
thumbnail Fig. 8.

Performances of the approximate forms of the covariance matrix in terms of rigenvalues and χ2-distributions. Top panel: eigenvalues of the covariance matrices (rebinned into six bins) compared to what can be obtained from the proposed approximate forms; same color-coding as for Fig. 6. The χ2 distributions are shown in the bottom panel. Model (86) reproduces the very same χ2 distributions. Model (85), in gray, is not as accurate and tends to slightly overestimate the χ2. This latter behavior is amplified when a larger number of bins is used.

In the text
thumbnail Fig. 9.

Scale dependence of the matter correlation functions for a realistic cosmological model (cosmological parameters derived from Plank, Planck Collaboration VI 2020) for the 3D density and the projected density (for a uniformly sampled survey with a depth of about 800 h−1Mpc between z = 0.75 and z = 1.25). The top panel shows (solid blue line) and (dashed red line) for the 3D density field, and the bottom panel shows and for the projected density. In both cases, the average value of the first moment of the two-point correlation function is dominated by large-distance contributions, whereas short-distance contributions dominate the second moment, assuming survey sizes of about 100 h−1 Mpc or above.

In the text
thumbnail Fig. B.1.

Diagrammatic visualization of the bias functions. The function φ0(λ) is the generating function of all trees within one cell, φ1(λ) of all trees within one cell with one external lines, and φ2(λ) with two external lines. The resulting connected diagrams up to second order in ξ12 are thus those presented here. Two φ2 generating functions cannot be conntected as that would induce a loop contribution. This reflects the underlying tree structure.

In the text
thumbnail Fig. B.2.

Joint CGF as a function of d and for different values of λ1, λ2: from bottom to top, λ1 = λ2 = −0.1, λ1 = λ2 = 0.2, λ1 = λ2 = 0.5. The solid blue lines correspond to the two-cell mean-field expression, (C.13). The predictions given in Eqs. (B.2) and (B.15) are shown as dashed red lines and dotted black lines. The shaded area is the region of overlapping cells.

In the text
thumbnail Fig. C.1.

Performances of the perturbative expansions of the joint PDF in the mean-field approximation either for the ξ12 expansion (open blue dots) or the expansion (red dots) up to 11th and 10th order, respectively. The comparisons are made for ρ1 = ρ2 = 2 and (top panel) and for ρ1 = 0.5,  ρ2 = 3.5 and (bottom panel) and for ξ12 equalling 0.1, 0.3, 0.5, 0.7, and 0.9.

In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.