The VIMOS Public Extragalactic Redshift Survey (VIPERS) based on observations collected at the European Southern Observatory, Cerro Paranal, Chile, using the Very Large Telescope under programmes 182.A-0886 and partly 070.A-9007. Also based on observations obtained with MegaPrime/MegaCam, a joint project of CFHT and CEA/DAPNIA, at the Canada-France-Hawaii Telescope (CFHT), which is operated by the National Research Council (NRC) of Canada, the Institut National des Science de l’Univers of the Centre National de la Recherche Scientifique (CNRS) of France, and the University of Hawaii. This work is based in part on data products produced at TERAPIX and the Canadian Astronomy Data Centre as part of the Canada-France-Hawaii Telescope Legacy Survey, a collaborative project of NRC and CNRS. The VIPERS web site is http://www.vipers.inaf.it/.

# The VIMOS Public Extragalactic Redshift Survey (VIPERS) ††thanks: based on observations collected at the European Southern Observatory, Cerro Paranal, Chile, using the Very Large Telescope under programmes 182.A-0886 and partly 070.A-9007. Also based on observations obtained with MegaPrime/MegaCam, a joint project of CFHT and CEA/DAPNIA, at the Canada-France-Hawaii Telescope (CFHT), which is operated by the National Research Council (NRC) of Canada, the Institut National des Science de l’Univers of the Centre National de la Recherche Scientifique (CNRS) of France, and the University of Hawaii. This work is based in part on data products produced at TERAPIX and the Canadian Astronomy Data Centre as part of the Canada-France-Hawaii Telescope Legacy Survey, a collaborative project of NRC and CNRS. The VIPERS web site is http://www.vipers.inaf.it/.

On the correct recovery of the count-in-cell probability distribution function.
###### Key Words.:
Cosmology: cosmological parameters – cosmology: large scale structure of the Universe – Galaxies: high-redshift – Galaxies: statistics
offprints: Bel., J.

We compare three methods to measure the count-in-cell probability density function of galaxies in a spectroscopic redshift survey. From this comparison we found that when the sampling is low (the average number of object per cell is around unity) it is necessary to use a parametric method to model the galaxy distribution. We used a set of mock catalogues of VIPERS, in order to verify if we were able to reconstruct the cell-count probability distribution once the observational strategy is applied. We find that in the simulated catalogues, the probability distribution of galaxies is better represented by a Gamma expansion than a Skewed Log-Normal. Finally, we correct the cell-count probability distribution function from the angular selection effect of the VIMOS instrument and study the redshift and absolute magnitude dependency of the underlying galaxy density function in VIPERS from redshift to . We found very weak evolution of the probability density distribution function and that it is well approximated, independently from the chosen tracers, by a Gamma distribution.

## 1 Introduction

The galaxy clustering offers a formidable playground to try to understand how structures have been growing during the evolution of the universe. A number of statistical tools have been developed and used over the past thirty years (see Bernardeau et al. 2002, for a review). In general, these statistical methods use the fact that the clustering of galaxies is due to the gravitational pull of the underlying matter distribution. Hence, the study of the spatial distribution of galaxies in the universe allows us to get information about the statistical properties of its matter content. As a result, it is of paramount importance to be able to measure the statistical quantities describing the galaxy distribution from a redshift survey.

The development of multi-object spectrographs on 8-m class telescopes during the 1990s triggered a number of deep redshift surveys with measured distances beyond over areas of 1–2 deg (e.g. VVDS Le Fevre et al. 2005, DEEP2 Newman et al. 2012 and zCOSMOS Lilly et al. 2009). Even so, it was not until the wide extension of VVDS was produced (Garilli et al. 2008), that a survey existed with sufficient volume to attempt cosmologically meaningful computations at (Guzzo et al. 2008). In general, clustering measurements at from these samples remained dominated by cosmic variance, as dramatically shown by the discrepancy observed between the VVDS and zCOSMOS correlation functions at (de la Torre et al. 2010).

The VIMOS Public Extragalactic Redshift Survey (VIPERS) is part of this global attempt to take cosmological measurements at to a new level in terms of statistical significance. In contrast to the BOSS and WiggleZ surveys, which use large-field-of-view ( deg) fibre optic positioners to probe huge volumes at low sampling density, VIPERS exploits the features of VIMOS at the ESO VLT to yield a dense galaxy sampling over a moderately large field of view ( deg). It reaches a volume at comparable to that of the 2dFGRS (Colless et al. 2001) at , allowing the cosmological evolution to be tested with small statistical errors.

The VIPERS redshifts are being collected by tiling the selected sky areas with a uniform mosaic of VIMOS fields. The area covered is not contiguous, but presents regular gaps due to the specific footprint of the instrument field of view, in addition to intrinsic unobserved areas due to bright stars or defects in the original photometric catalogue. The VIMOS field of view has four rectangular regions of about square arcminutes each, separated by an unobserved cross (Guzzo et al. 2014; de la Torre et al. 2013). This creates a regular pattern of gaps in the angular distribution of the measured galaxies. Additionally, the Target Sampling Rate and the Survey Success Rate vary among the quadrants, and a few of the latter were lost because of mechanical problems within VIMOS (Garilli et al. 2014). Finally, the slit-positioning algorithm, SPOC (see Bottini et al. 2005), also introduces some small-scale angular selection effects, with different constraints along the dispersion and spatial directions of the spectra, as thoroughly discussed in de la Torre et al. (2013). Clearly, this combination of angular selection effects has to be taken properly into account when estimating any clustering statistics.

In this paper we measure the probability distribution function from the VIPERS Public Data Release 1 (PDR-1) redshift catalogue, including of the final number of redshifts expected at completion (see Guzzo et al. 2014; Garilli et al. 2014 for a detailed description of the survey data set). The paper is organized as follows. In §2 we introduce the VIPERS survey and the features of the PDR-1 sample. In §3 we review the basics of the three methods we compared. In §4 we present a null test of the three method on a synthetic galaxy catalogue. In §5 we use galaxy mock catalogues to assess performances of two of the methods. Magnitude and redshift dependance of the probability distribution function of VIPERS PDR-1 galaxies are presented in §6 and conclusions are drawn in §7.

Throughout, the Hubble constant is parameterized via , all magnitudes in this paper are in the AB system (Oke & Gunn 1983) and we will not give an explicit AB suffix. In order to convert redshifts into comoving distances we assume that the matter density parameter is and that the universe is spatially flat with a CDM cosmology without radiations.

## 2 Data

The VIMOS Public Extragalactic Redshift Survey (VIPERS) is a spectroscopic redshift survey being built using the VIMOS spectrograph at the ESO VLT. The survey target sample is selected from the Canada-France-Hawaii Telescope Legacy Survey Wide (CFHTLS-Wide) optical photometric catalogues (Mellier et al. 2009). The final VIPERS will cover deg on the sky, divided over two areas within the W1 and W4 CFHTLS fields. Galaxies are selected to a limit of , further applying a simple and robust colour pre-selection, as to effectively remove galaxies at . Coupled to an aggressive observing strategy (Scodeggio et al. 2009), this allows us to double the galaxy sampling rate in the redshift range of interest, with respect to a pure magnitude-limited sample (). At the same time, the area and depth of the survey result in a fairly large volume, h Mpc, analogous to that of the 2dFGRS at . Such combination of sampling and depth is quite unique over current redshift surveys at . The VIPERS spectra are collected with the VIMOS multi-object spectrograph (Le Fevre et al. 2003) at moderate resolution (), using the LR Red grism, providing a wavelength coverage of 5500-9500 and a typical redshift error of km sec . The full VIPERS area is covered through a mosaic of 288 VIMOS pointings (192 in the W1 area, and 96 in the W4 area). A discussion of the survey data reduction and management infrastructure is presented in Garilli et al. (2012). An early subset of the spectra used here is analyzed and classified through a Principal Component Analysis (PCA) in Marchetti et al. (2012).

A quality flag is assigned to each measured redshift, based on the quality of the corresponding spectrum. Here and in all parallel VIPERS science analyses we use only galaxies with flags 2 to 9 inclusive, corresponding to a global redshift confidence level of 98%. The redshift confirmation rate and redshift accuracy have been estimated using repeated spectroscopic observations in the VIPERS fields. A more complete description of the survey construction, from the definition of the target sample to the actual spectra and redshift measurements, is given in the parallel survey description paper (Guzzo et al. 2014).

The data set used in this paper and the other papers of this early science release is the VIPERS Public Data Release 1 (PDR-1) catalogue, which have been made publicly available in September 2013. This includes objects, spread over a global area of deg and deg respectively in W1 and W4. It corresponds to the data frozen in the VIPERS database at the end of the 2011/2012 observing campaign, i.e. 64% of the final expected survey. For the specific analysis presented here, the sample has been further limited to its higher-redshift part, selecting only galaxies with . The reason for this selection is related to minimizing the shot noise and maximizing the volume. This reduces the usable sample to and galaxies in W1 and W4 respectively (always with quality flags between and ). The corresponding effective volume of the two samples are and h Mpc. At redshift they spann respectively the angular comoving distances and . We divide the W1 and W4 fields in three redshift bins and we build magnitude limited sub-samples in each of them. For convenience, we use the magnitude limits listed in Table (1) of di Porto et al. (2014), which we recall in Tab. (1).

The VIMOS footprint has an important impact on the observed probability of finding galaxies in a randomly placed spherical cell in the survey volume. As a matter of fact, a direct appreciation of the masked area can be shown on the first moment of the probability distribution, i.e. the expectation value of the number count . On one hand, we can predict the mean number of objects per cells from the knowledge of the number density in each considered redshift bins and on the other hand we can estimate it by placing a regular grid of spherical cells of radius into the volume surveyed by VIPERS. In fact, given the solid angle of W1 and W4 and the corresponding number of galaxies and contained in a redshift bin extracted from each field, one can estimate the total number density as

 ¯ρ=N1+N4Ω1+Ω41Vk, (1)

where is defined as the volume corresponding to a sector of a spherical shell with solid angle equal to unity. In the case of VIPERS PDR-1 the effective solide angles corresponding to W1 and W4 are respectively and (in square radians). One can therefore predict the corresponding expected number of objects in each cell by multiplying the averaged number density by the volume of a cell. It reads

 ¯NR=43πR3¯ρ, (2)

in the case of the spherical cells of radius considered in this work. The expectation value with respect to the radius of the cells corresponding to each luminosity sub-sample extracted from VIPERS-PDR1 is represented by lines in Fig. (1). On the same figure we display the measured mean number of object in each redshift bins. Note that to perform this measurement we place a grid of equally separated (Mpc) spheres of radius Mpc and we reject spheres with more than % of their volume outside the observed region (see Bel et al. 2014). We quantify the effect of the mask using the quantity

 α≡¯N¯NR, (3)

in fact the botton panels of Fig. (1) shows that for all sub samples and at all redshifts the neat effect of the masks is to under-sample the galaxy field by roughly %. It also shows that the correction factor depends on the considered redshift, on the luminosity and on the cell-size. The scale dependency can be explained by the fact that the correction parameter depends on how the cells overlap with the masked regions. The left panel of Fig. (1) suggests that at low redshift the mask effect behaves in the same way for all the luminosity samples while the middle panel shows a clear dependency with respect to luminosity. The correction factor depends on the redshift distribution, as a result the apparent dependency with respect to the luminosity is due to the dependence of the number density with respect to the luminosity of the considered objects.

The mask not only modifies the mean number of object but it also modifies the higher order moments of the distribution, such that the measured will be systematically altered. In the present paper we show that this systematic effect can be taken into account by measuring the underlying probability density function of the galaxy density contrast . It has been shown (see Fig. (8) of Bel et al. 2014) that after rejecting spheres with more than % of their volume outside the survey, the local poisson process approximation holds. In particular, it allows to use the “wrong” probability distribution function in order to get reliable information on the underlying probability density function . Then applying the Poisson sampling one can recover the unaltered using that . For the sake of completeness we provide the reader with the measured probability function obtained after rejecting the cells with more than % of their volume outside the survey (see Fig. 8).

In particular, let and , respectively, be the observed and the true Counting Probability Distribution Function (CPDF). Assuming that from the knowledge of there exists a process to get the underlying probability density function of the stochastic field , which is associated to the random variable , one can compute the true CPDF applying

 PN=∫∞0P[N|Λ]p(Λ)dΛ, (4)

where is called the sampling conditional probability; it determines the sampling process from which the discrete cell-count arises. In the following we assume that this sampling conditional probability follows a Poisson law (Layzer 1956), as a result in Eq. (4) we substitute

 P[N|Λ]=K[N,Λ]≡ΛNN!e−Λ. (5)

It is also convenient to express Eq. (4) in terms of the density contrast of the stochastic field , , it follows that

 PN=∫∞−1K[N|¯N(1+δ)]p(δ)dδ, (6)

where we used that , which is a property of the Poisson sampling.

Continuing along this direction that we propose to compare three methods which aim at extracting the underlying probability density function (PDF) in order to correct the observed CPDF from the angular selection effects of VIPERS.

## 3 Methods

In this section we review the PDF estimators that we use and compare with each others in this paper. The purpose is to select the method which will be more adapted to the VIPERS characteristics.

### 3.1 The Richardson-Lucy deconvolution

This is an iterative method which aims at inverting Eq. (6) without parametrizing the underlying PDF, it has been investigated by Szapudi & Pan (2004). This method starts with an initial guess for the probability density function which is used to compute the corresponding expected observed via

 PN,0=∫∞−1^K[N,¯N(1+δ)]p0(δ)dδ,

where . The probability density function used at the next step is obtained using

 ^pi+1(δ)=^pi(δ)Nmax∑N=0PNPN,i^K[N,¯N(1+δ)],

where . For each step the agreement between the expected observed probability distribution and the true one is quantified by

 χ2i≡Nmax∑N=0(PNPN,i−1)2.

It is therefore possible to know the evolution of the cost function with respect to the steps .

In fact it has been shown by Szapudi & Pan (2004) that it converges toward a constant value which corresponds to the best evaluation of the probability density function given the observed probability distribution . Since these authors have shown that this convergence occurs after around iterations. We did our own convergence tests which have shown that adopting a value of iterations is enough. However, it happens that the evolution of the is not always monotonic. In practice, we store the result of each step and we look for the step for which the is minimum, i.e. . As an initial guess we set that the discret CPDF is equal to the continuous one ().

### 3.2 The Skewed Log-Normal

This is a parametric method where the shape of the probability density depends on a given number of parameters, in this case the probability density function is assumed to be well described by a Skewed Log-Normal (Colombi 1994) distribution. It is derived from the Log-Normal distribution (Coles & Jones 1991) but it is more flexible. It is indeed built upon an Edgeworth expansion; be the stochastic field , following a Normal distribution then the density contrast follows instead a Log-Normal distribution. In the case of the Skewed Log-Normal (SLN) density function, the field follows an Edgeworth expanded Normal distribution

 PΦ(Φ)≡{1+⟨ν3⟩c6H3(ν)+⟨ν4⟩c24H4(ν)+572⟨ν3⟩2cH6(ν)}G(ν)σΦ, (7)

where , is the central reduced Normal distribution and denotes the cumulant expectation value of . As a result, the SLN is parameterized by the four parameters , , , and which are related, respectively to the mean, the dispersion, the skewness and the kurtosis of the stochastic variable . They can all be expressed in terms of cumulants of order of the weakly non-Gaussian field . In Szapudi & Pan (2004) they use a best fit approach and determine these parameters by minimizing the difference between the measured counting probability and the one obtained from

 PthN= ∫∞−1K[N,¯N(1+δ)]PΦ[ln(1+δ),μΦ,σ2Φ,⟨Φ3⟩c,⟨Φ4⟩c] (8) ×dln(1+δ).

However, this requires us to perform the integral (Eq. 8) in a four dimensional parameter space which is numerically expensive.

In the present paper we use an alternative implementation which is computationally more efficient. Instead of trying to maximize the likelihood of the model given the observations, we rather use the observations to predict the parameters of the SLN. To do so we use the property of the local Poisson sampling (Bel & Marinoni 2012); the factorial moments of the discrete counts are equal to the moments of the underlying continuous distribution . Since the transformation between the density contrast and the Edgeworth expanded field is local and deterministic, it is possible to find a relation between the moments and the cumulants .

By definition, the moments of the positive continuous field are given by

 ⟨Λn⟩≡∫∞0ΛnP(Λ)dΛ,

then for a local deterministic transformation the conservation of probability imposes , it follows that the moments of can be recast in terms of

 ⟨Λn⟩=¯Λn∫∞0enΦPΦ(Φ)dΦ.

In the right hand side one can recognize the definition of the moment generating function we therefore obtain that

 MΦ(t=n)=⟨Λn⟩¯Λn≡An. (9)

This equation allows us to link the moment of to the cumulants of via the moment generating function .

Moreover, since the probability density is the product of a sum of Hermite polynomials with a Gaussian function it is straightforward to compute the explicit expression of the moment generating function we obtain

 MΦ(t)={1+⟨Φ3⟩ct36+⟨Φ4⟩ct424+⟨Φ3⟩2c572t6}etμΦ+t2σ2Φ2. (10)

As a matter of fact, Eq. (10) and Eq. (9) together allow to set up a system of four equations, for it reads

 Yn2XnBn=An, (11)

where , and . In the system of equations (Eq. 11) the right hand side is given by observations and the left hand side depends on the cumulants , , and parameterized in terms of , , and . In appendix (A), we detail the procedure to solve this non-linear system of equations.

We therefore, get the values of the four parameters of the SLN by simply measuring the moments of the counting variable up to the fourth order.

### 3.3 The Gamma expansion

The Gamma expansion method follows the same idea as described in §3.2 but it uses a Gamma distribution instead of a Gaussian one. It uses the orthogonality properties of the Laguerre polynomials in order to modify the moments of the Gamma PDF. Such an expansion has been investigated in Gaztañaga, Fosalba & Elizalde (2000) where they compared it to the Edgeworth expansion in order to model the one-point PDF of the matter density field. Then it has been further extended, in a more general context, to multi-point distributions by Mustapha & Dimitrakopoulos (2010).

As mentioned above the Gamma expansion requires the use of the Gamma distribution defined as

 ϕG(u)≡uk−1θΓ(k)e−u, (12)

where is the Gamma function (for an integer , , and are two parameters which are related to the two first moments of the PDF. If the galaxy probability density function is well described by a Gamma expansion at order then it can be formally written as

 P(Λ)=ϕG(u)f(k−1)n(u), (13)

where by definition , , . The function represents the expansion aiming at tuning the moments of the Gamma distribution; note that the exponent is not the derivative of order . Since this expansion is built upon the orthogonal properties of products of Laguerre polynomials with the Gamma distribution, the function is given by the sum

 f(k−1)n(x)≡n∑i=0ciL(k−1)i(x), (14)

where are the generalized Laguerre polynomials of order and the coefficients represent the coefficients of the Gamma expansion and therefore depend on the moments of the galaxy field

 cn≡n∑i=0(ni)Γ(k)Γ(k+i)(−1)i⟨Λi⟩θi. (15)

The main interrest of the Gamma expansion with respect to the SLN is that the coefficients of the expansion are directly related to the moments of the distribution we want to model, i.e. it is not necessary to solve a complicated non-linear system of equations nor to perform a Likelihood estimation of the coefficients. Moreover, it can be easily performed at higher order to describe as best as possible the underlying probability density function of galaxies.

Another advantage of describing the galaxy field by a Gamma expansion probability density function is that the corresponding observed can be expressed analytically, which is not the case for the SLN which must be integrated numerically.

In Appendix (B) we demonstrate the previous statement, it follows that the CPDF can be calculated from

 PN=(−θ)NN!n∑i=0ciΓ(i+k)Γ(k)h(N)i(θ), (16)

where and in this case we use the notation . The successive derivatives of can be obtained from the recursive relation

 h(N)i(θ)=(i)NfθNhi(θ)−N∑m=1(Nm)(i+k)mf(1+θ)mh(N−m)i(θ).

In addition to the fact that having the possibility of computing the corresponding observed without requiring an infinite integral for each number is computationally more efficient, it is also practical to have the analytical calculation for some peculiar values of the parameter of the distribution. In fact, when is lower than which occurs on small scales (Mpc), the probability density function goes to infinity when goes to (although the distribution is still well defined). In particular, this numerical divergence would induce large numerical uncertainties in the computation of the void probability . In addition, one can see that for the void probability we have the simple relation

 P0=n∑i=0ciΓ(k+i)Γ(k)hi(θ), (17)

which can be used to recover the true void probability in VIPERS.

## 4 Application of the methods on a synthetic galaxy distribution

In this section we analyse a suite of synthetic galaxy distributions generated from realizations of a Gaussian stochastic field. The full process involved in generating these bench-mark catalogues is detailed in Appendix C. Each comoving volume has a cubical geometry of size Mpc. We generate the galaxies by discretizing the density field according to the sampling conditional probability which we assume to be a Poisson distribution with mean . In this way we know the true underlying galaxy density contrast . We can therefore perform a fair comparison between the methods introduced in §3.

In order to avoid the effect of the grid (Mpc) we smooth both the density field and the discrete field using a spherical Top-Hat filter of radius Mpc. We apply the three methods mentioned in §3 and compare the reconstructed probability density function to the expected one obtained directly from the density field .

The discrete distribution of points contains an average number of object per cell which is the one expected according to our sampling process. The corresponding is given by the black histogram in the lower panel of Fig. (2), from this measurement we apply the three methods R-L, SLN and and obtain an estimation of the probability density function corresponding to each method. In the upper panel of Fig. (2) we compare the performance of the three methods in recovering the true probability density function (black histogram referred to as reference in the inset). Note that, for this test case, we use a Gamma expansion at order in order to be coherent with the order of the expansion of the Skewed Log-Normal. We have also represented the probability density function estimated when neglecting the shot noise (red dotted line), which is used as the initial guess in the case of the R-L method.

From the top panel of Fig. (2) we can conclude that the three methods perform reasonably well. It seems that the method reproduces better the density distribution of under-dense regions () but this is expected in the sense that the distribution used to generate the synthetic catalogues is a Gamma distribution (see Appendix C). Although, it is not obvious because the scale on which the density field has been set up is one order of magnitude smaller than the scale of the reconstruction Mpc.

The performance of the three methods is also represented in the bottom panel of Fig. (2), in which we compare the expected observed from each method to the true one. One can see that they all agree at the % level, hence it is not possible to conclude that one is better than an other. This was actually expected, from the comparison on the underlying density field (Fig. 2). On the contrary if one of the methods would not agree with the PDF then we would expect also a disagreement on the observed CPDF (see §6).

In the following part we investigate the sensitivity of the three methods with respect to the shot noise. In fact, as shown in Fig. (1), in most of the sub-samples of VIPERS PDR-1 we will work with a high shot noise level (). We therefore randomly under-sample the fake galaxy distribution by keeping only % of the total number of object contained in each comoving volume. This process gives an average number per cell of , which is more representative in the context of the application of the reconstruction method. We perform the same comparison as in the ideal case () and found that the R-L method appears to be highly sensitive to the shot noise. In fact if the mean number of object per cell is too few then the output of the method depends too much on the initial guess. It follows that, if it is too far from the true PDF the process does not converge (see top panel of Fig. 3) and the corresponding does not match the observed one (see bottom panel of Fig. 3). Note that we explicitly checked this effect by increasing the number of iterations from to . While in the case of both, the SLN and the Gamma expansion, one can see in Fig. (3), the output probability density function is in agreement (with a larger scatter) with the one obtained in the case. This means that the sensitivity regarding to the shot noise is much smaller when considering parametric methods.

Considered the sensitivity of the R-L method to the initial guess, knowing that the average number of galaxies per cell can be lower than unity and finally taking into account computational time, we shall continue our analysis only using the two parametric methods SLN and . In the following, we will compare them using more realistic mock catalogues but for which we don’t know apriori the true underlying PDF.

## 5 Performances in realistic conditions

In this section we discuss how observational effects have been accounted for in our analysis and test the robustness of the reconstruction methods SLN and Gamma expansion. For this purpose we use a suite of mock catalogues created from the Millenium simulation, they are also used in the analysis performed by di Porto et al. (2014).

We shall compare the reconstruction methods between two catalogues, namely REFERENCE and MOCK. The reference is a galaxy catalogue obtained from semi-analytical models. We simulate the redshift errors of VIPERS PDR-1 by perturbing the redshift (including distortions due to peculiar motions) with a Normally distributed error with rms . Each MOCK catalogue is built from the corresponding REFERENCE catalogue by applying the same observational strategy (de la Torre et al. 2013) which is applied on VIPERS PDR-1; spectroscopic targets are selected from the REFERENCE catalogue by applying the slit-positioning algorithm (SPOC, Bottini et al. 2005) with the same setting as for the PDR-1. This allows us to reproduce the VIPERS footprint on the sky, the small-scale angular incompleteness due to spectra collisions and the variation of target sampling rate across the fields. Finally, we deplete each quadrant to reproduce the effect of the survey success rate (SSR, see de la Torre et al. 2013). In this way, we end up with 50 realistic mock catalogues (named MOCK hereafter), which simulate the detailed survey completeness function and observational biases of VIPERS in the W1 and W4 fields.

In order to perform a similar analysis as the one we aim at doing for VIPERS PDR-1, we construct sub-samples of galaxies selected according to their absolute magnitude in B-band; we take all objects brighter than a given luminosity. We list those samples in Tab. (2), we have in total galaxy samples. The highest luminosity cut (M) allows us to follow a single population of galaxies at three cosmic epocs.

In Fig. 4, 5 and 6 we show the reconstruction performances for the SLN and the method. We consider the same population () but in three redshift bins, , and . In order to test the stability of the methods we perform the reconstruction at three smoothing scales, , and Mpc. The comparison is done as follows, on one hand we estimate the true from the REFERENCE catalogue (before applying the observational selection) and we perform the reconstruction on it, in this way we can test the intrinsic biases due to the assumed parametric method (SLN or ). On the other hand, we estimate the observed in the MOCK catalogues, from which we perform the reconstruction to verify if we recover the expected from the REFERENCE catalogue.

Inspecting Fig. (4) we can first see that the intrinsic error due to the specific modeling of the methods is much larger for the SLN (cyan diamonds compared to the black histogram) than for the (magenta diamonds compared to the black histogram). From the top panels we see that the SLN does not reproduce the tail of the CPDF and from the bottom panel we see that even for low counts it is showing deviations as large as %. This intrinsic limitation is propagating when performing the reconstruction on the MOCK catalogue (orange triangles compared to the black histogram) while for the we see that the agreement is better than % (magenta triangles compared to the black histogram) in the low count regime and the tail is fairly well reproduced. In the second place, comparing the performed on the REFERENCE and the MOCK catalogues (blue diamonds with respect to magenta triangles) one can see the loss of information due to the observational strategy has at most an impact of % on the reconstructed CPDF which reduces when considering larger cells (less shot noise).

In general, examination of Fig. (5) and (6) confirms that for the considered galaxy population the same results hold at lower redshifts. However, in particular the reconstruction at Mpc can exhibit deviations larger than %, this is at odds with the fact that the shot noise contribution is expected to be the same for the three redshift bins (magnitude limited). We attribute this larger instability to the fact that not only the shot noise contribution is higher for Mpc but also the volume probed is smaller when decreasing the redshift.

The performances of the reconstruction for the last three galaxy samples are shown in Fig.(7) where each row corresponds to a galaxy sample (we only show the residual with respect to the REFERENCE). This last comparison allows to say that the reconstruction instability at Mpc was indeed due to the high level of shot noise. We can conclude that in the HOD galaxy mock catalogues, the galaxy distribution is more likely to be modelled by a instead of an SLN. Finally, for a chosen reconstruction method, the information contained in the MOCK catalogues is enough to be able to reconstruct the CPDF of the REFERENCE catalogue to better than 10%.

## 6 VIPERS PDR-1 data

In this section we apply the reconstruction method to the VIPERS PDR-1. We saw in the previous sections that the SLN and methods are sensitive to the assumption we make about the underlying PDF. In fact, we saw in §4 that if the underlying PDF is close to the chosen model then the reconstruction works. We found in §5 that the galaxy distribution arising from semi-analytic models is better described by a than an SLN distribution. However, in the following we will not take for granted that the same property holds for the galaxies in the PDR-1.

We want to choose which one of the two distributions (Log-Normal or Gamma) best describes the observed galaxy distribution in VIPERS PDR-1, when no expansion is applied. Thus, we compare the observed PDF to the one expected from the Poisson sampling of the Log-Normal probability density function (PS-LN) and to the one expected from the Poisson sampling of the Gamma distribution (the so-called Negative Binomial). Error bars are obtained by performing a Jack-knife resampling of subregions in each fields W1 and W4.

The SP-LN distribution does not have an analytic expression and must be obtained by numerically integrating Eq. (6) while the Poisson sampling of the Gamma distribution leads to the Negative Binomial distribution defined as

 PN=θNN!r(r+1)...(r+N−1)(1+θ)N+r, (18)

where and to ensure that the first two moments of the Negative Binomial match those of the observed distribution. We show in Fig. (8) the outcome of this comparison, it follows that the Negative Binomial is much closer to the observed PDF than the PS-LN. As a result, the underlying galaxy distribution is more likely to be described by a Gamma distribution than by a Log-Normal. Hence, we only use the Gamma expansion to model the galaxy distribution of VIPERS PDR-1.

Moreover, the use of the Gamma expansion instead of the SLN simplifies substantially the analysis. In Fig. (9) we provide the reconstructed probability distribution function of VIPERS PDR-1 together with the corresponding underlying probability density function for each redshift bin and luminosity cut. Each panel of Fig. (9) shows how the choice of a particular class of tracers (selected according to their absolute magnitude in B-band) influence the PDF of galaxies. When measuring specific properties of the intrinsic galaxy distribution for each luminosity cut, it is enough to look at the CPDF however, when comparing the distributions with each other it is necessary to take care about the averaged number of objects per cell which varies from sample to sample. As a result it appears more useful to compare the properties of the different galaxy samples using their underlying probability density function which, assuming Poisson sampling, is free from sampling rate variation between different type of tracers.

For the two first redshift bins, we can see that the probability density function is broadening when selecting more luminous galaxies, this goes in the direction of increasing the linear bias with respect to the matter distribution. However, despite a less significant trend, for the highest redshift bin it seems that it goes in the oposite direction. This trend might be an artifact; indeed by analyzing Fig. (1) we see that for all these samples the averaged number of object per cell is between - which shows that theses samples could be highly affected by shot noise effects. As a result, specific care should be taken when interpreting those three high redshift samples.

In the following we focus on the evolution of the underlying PDF for a particular class of objects on the wide redshift range probed by VIPERS PDR-1. The Fig. (10) displays the outcome of this study, it shows how the PDF, for three populations (the three highest magnitude cuts), evolves regarding to the redshift at which it is measured. The three populations (top, middle and bottom panels) exhibit non-monotonic evolution with respect to the redshift. In particular, the more luminous population is showing that the PDF at appears to be systematically different than in the two lower redshift bins. However, we see also that some instabilities are appearing in the reconstruction (see wiggles at high ). This might be due to the fact that we have fewer galaxies in this sample giving rise to a large shot-noise contribution (). We indeed verified that for the high mass bin and the two other galaxy populations we vary the order of the expansion from to the resulting PDF changes by less than - while for the most luminous population, truncating the expansion at order only removes the instability without changing significantly the overall behavior of the PDF. This consistency test shows that the radical change in the measured PDF for the highest redshift bin appears to be the true feature. Probably only the final VIPERS data set will be able to give a robust conclusion.

Finally, in Tab. (3), we list the relevant coefficients of the Gamma expansion which we measured from the VIPERS PDR-1 at the scale Mpc. They can be used in order to model both the CPDF (Eq. 16) and the PDF (Eq. 13).

## 7 Summary

The main goal of the present paper is to measure the probability of finding galaxies falling into a spherical cell randomly placed inside a sparse sampled (i.e. with masked areas or with low sampling rate) spectroscopic survey. Our general approach to this problem has been to use the underlying probability density distribution of the density contrast of galaxies in order to recover the counting probability corrected from sparseness effects. We therefore compared three ways (R-L, SLN and ) of measuring the probability density of galaxies classified in two categories; direct and parametric. We found that when the sampling is high () the direct method (Rychardson-Lucy deconvolution) performs well and avoids putting any prior on the shape of the distribution. On the other hand, we saw that when the sampling is low () the direct method fails to converge to the true underlying distribution. We thus concluded that, in such cases, the only alternative is to use a parametric method.

We presented two parametric forms aiming at describing the galaxy density distribution, the SLN which is often used in the literature to model the matter distribution and the . Despite the fact that the two distributions used in this paper have been already investigated in previous works, the approach we propose to estimate their parameters is completely new. Previously, fitting procedures were used in order to estimate them. Here we propose to measure directly the parameters of the distributions from the observations. The method can be applied to both distributions SLN and and decreases considerably the computational time of the process.

Relying on simulated galaxy catalogues of VIPERS PDR1, we tested the reconstruction scheme of the counting probability () under realistic conditions in the case of the SLN and expansions. We found, that the reconstruction depends on the choice of the model for the galaxy distribution. However, we have also shown that it is possible to test which distribution better describes the observations.

Using VIPERS PDR1, on the relevant scales investigated in this paper (Mpc), we found that the distribution gives a better description of the observed than the one provided by the Log-Normal (see Fig. 8). We therefore adopted the parametric form in order reconstruct the probability density functions of galaxies. From these reconstruction we studied how their PDF changes according to their absolute luminosity in B-band and we also studied their redshift evolution. We found that little evolution has been detected in the two first redshift bins while it seems that the density distribution of the galaxy field is strongly evolving in the last redshift bin.

We finally used, the measured pdf in order to reconstruct the counting probability (CPDF) one would observe if VIPERS was not masked by gaps between the VIMOS quadrants.

###### Acknowledgements.
JB acknowledges useful discussions with E. Gaztañaga. We acknowledge the crucial contribution of the ESO staff for the management of service observations. In particular, we are deeply grateful to M. Hilker for his constant help and support of this programme. Italian participation to VIPERS has been funded by INAF through PRIN 2008 and 2010 programmes. JB, LG and BJG acknowledge support of the European Research Council through the Darklight ERC Advanced Research Grant (# 291521). OLF acknowledges support of the European Research Council through the EARLY ERC Advanced Research Grant (# 268107). AP, KM, and JK have been supported by the National Science Centre (grants UMO-2012/07/B/ST9/04425 and UMO-2013/09/D/ST9/04030), the Polish-Swiss Astro Project (co-financed by a grant from Switzerland, through the Swiss Contribution to the enlarged European Union), the European Associated Laboratory Astrophysics Poland-France HECOLS and a Japan Society for the Promotion of Science (JSPS) Postdoctoral Fellowship for Foreign Researchers (P11802). GDL acknowledges financial support from the European Research Council under the European Community’s Seventh Framework Programme (FP7/2007-2013)/ERC grant agreement n. 202781. WJP and RT acknowledge financial support from the European Research Council under the European Community’s Seventh Framework Programme (FP7/2007-2013)/ERC grant agreement n. 202686. WJP is also grateful for support from the UK Science and Technology Facilities Council through the grant ST/I001204/1. EB, FM and LM acknowledge the support from grants ASI-INAF I/023/12/0 and PRIN MIUR 2010-2011. CM is grateful for support from specific project funding of the Institut Universitaire de France and the LABEX OCEVU.

## References

• Bel & Marinoni (2012) Bel, J. & Marinoni, C. 2012, MNRAS, 424, 971
• Bel et al. (2014) Bel, J., et al. (the VIPERS Team) 2014, A&A, 563, A37
• Bernardeau et al. (2002) Bernardeau, F., Colombi, S., Gaztaãga, E. & Scoccimarro, R. 2002, PR, 367, 1
• Bottini et al. (2005) Bottini D., Garilli, B., Maccagni, D., et al. 2005, PASP, 117, 996
• Coles & Jones (1991) Coles, P., Jones, B. 1991, MNRAS, 248, 1
• Colless et al. (2001) Colless, M., Dalton, G., Maddox, S., et al. 2001, MNRAS, 328, 1039
• Colombi (1994) Colombi, S. 1994, ApJ, 435, 536
• de la Torre et al. (2010) de la Torre, S., Guzzo, L., Kovac, K., et al. (the ZCOSMOS collaboration) 2010, MNRAS, 409, 867
• de la Torre et al. (2013) de la Torre, S., Guzzo, L., Peacock, J.A., et al. (VIPERS team) 2013, A&A, 557A, 54
• di Porto et al. (2014) di Porto, C., et al. (VIPERS team) 2014, submitted, arXiv:1406.6692D
• Eisenstein & Hu (1998) Eisenstein, D. J. & Hu, W. 1998, ApJ, 496, 605
• Garilli et al. (2008) Garilli, B., Le Fèvre, O., Guzzo, L., et al. (the VVDS collaboration) 2008, A&A, 486, 683
• Garilli et al. (2012) Garilli, B., Paioro, L., Scodeggio, M. et al. 2012, PASP, 124, 1232
• Garilli et al. (2014) Garilli, B., Guzzo, L., Scodeggio, M., et al. (the VIPERS team) 2014, A&A, 562, 23
• Gaztañaga, Fosalba & Elizalde (2000) Gaztañaga, E., Fosalba, P. & Elizalde, E. 2000, ApJ, 539, 522
• Greiner & Enlin (2015) Greiner, M., Enlin, T. A. 2015, A&A, 574, 86
• Guzzo et al. (2008) Guzzo, L., Pierleoni, M., Meneux, B., et al. (the VVDS team) 2008, Nature, 451, 541
• Guzzo et al. (2014) Guzzo, L., Scodeggio, M., Garilli, B., et al. (the VIPERS team) 2014, A&A, 566, 108
• Layzer (1956) Layzer, D. 1956, AJ, 61, 383
• Le Fevre et al. (2003) Le Fèvre, O., Saisse, M., Mancini, D., et al. 2003, Proc. SPIE, 4841, 1670
• Le Fevre et al. (2005) Le Fèvre, O. Vettolani, G., Garilli, B., et al 2005, A&A, 439, 845
• Lilly et al. (2009) Lilly, S. J., Le Brun, V., Maier, C., et al. (the ZCOSMOS collaboration) 2009, ApJS, 184, 218
• Marchetti et al. (2012) Marchetti, A., Granett, B.R., Guzzo, L., et al. (the VIPERS team) 2012, MNRAS, in press, arXiv:1207.4374
• Mellier et al. (2009) Mellier, Y., Bertin, E., Hudelot, P., et al. 2008, The CFHTLS T0005 Release, http://terapix.iap.fr/cplt/oldSite/Descart/ CFHTLS-T0005-Release.pdf
• Mustapha & Dimitrakopoulos (2010) Mustapha , H. & Dimitrakopoulos, R. 2010, C&M, 60, 2178
• Newman et al. (2012) Newman, J. A., Cooper, M.C., Davis, M., et al. (the DEEP2 collaboration) 2012, arXiv:1203.3192
• Oke & Gunn (1983) Oke, J. B. & Gunn, J. E. 1983, ApJ, 266, 713
• Scodeggio et al. (2009) Scodeggio M., Franzetti P., Garilli B., et al. 2009, Msngr, 135, 13
• Szapudi & Pan (2004) Szapudi, I. & Pan, J. 2004, ApJ, 602, 26

## Appendix A Non-linear system

The problem of this system of equations is that it is non-linear, it is therefore difficult to solve however it can be reduced to a one dimensional equation which can be solved numerically.

The two first equations ( and ) can be used to express the two first cumulants with respect to the third and fourth order ones

 σ2Φ = ln(A2)+ln(B21B2) (19) μΦ = −12[ln(A2)+ln(B41B2)] (20)

where and are both functions of and . Then using other combinations of equation one can express a system of two equations for and alone

 B23 = a1B21B4 (21) B3B31 = a2B32, (22)

where and . In order to solve properly the system we prefer to express it in term of one parameter , moreover one can see that polynomials to are not independent, as a result

 B4=d+aB1+bB2+cB3,

where and which can be substituted in Eq. (21). Combining Eq. (21) and Eq. (22) one obtains a parametric equation for

 (a+bη)B31+(d+cf(η))B21−g(η)=0, (23)

which can be solved for each value of the parameter and an independent parametric equation for

 B3=f(η).

As a result we can find a couple , for each value of the parameter , it follows that one can express and with respect to and given the definition of the possible solution and must satisfy the condition

 B2[x(η),y(η)]−ηB1[x(η),y(η)]=0,

which gives the possible values of from which one can recover and . Finally, from Eq. (19) and Eq. (20) we can compute the values of and corresponding to each couple (, ) of solutions. This allows us to select the solution which provides a value of closer to the observed one.

Once the values of the cumulants , , and are known from the process detailed above, we know that the moments of the corresponding will match those of the observed on up to order . At the end, one can check whether the SLN distribution provides a good match to data by integrating numerically the probability density function convolved with the Poisson kernel (see Eq. 5).

## Appendix B Generating function

We show that the CPDF associated to a Gamma expanded PDF can be calculated analytically from an expression which depends explicitly on the coefficients of the Gamma expansion.

Be the generating function associated to the probability distribution , it is defined as

 GN(λ)≡∞∑i=0λNPN. (24)

In case of the Poisson sampling of a Gamma distribution, after some algebra, one can show that it can be expressed with respect to the coefficients of the Gamma expansion as

 GN(λ)=1Γ(k)n∑i=0ciFi(γ), (25)

where and

 Fi(γ)≡∫∞0xk−1e−xL(k−1)i(x)e−γxdx.

Nevertheless, this integral can be computed using the Laguerre expansion of the exponential

 e−γx=∞∑i=0γi(1+γ)i+α+1L(α)i(x),

 Fi(γ)=γi(1+γ)i+kΓ(i+k)i!. (26)

The formal expression of the generating function is therefore given by

 GN(λ)=(1+γ)−kΓ(k)n∑i=0ciΓ(i+k)i!(γ1+γ), (27)

where we still use . From the explicit expression of the moment generating function (Eq. 27) one can get the probability distribution by iteratively deriving the generating function with respect to

 PN≡1N!dNGN(λ)dλN∣∣ ∣∣λ=0=(−θ)NN!dNGN(γ)dγN∣∣ ∣∣γ=θ.

These derivatives can be calculated explicitly.

## Appendix C Synthetic galaxy catalogues

In this Appendix we describe how we generate synthetic galaxy catalogues from Gaussian realizations. The first requirement of these catalogues is that they must be characterized by a known power spectrum and 1-point probability distribution function. The second requirement is that the probability distribution function must be measurable.

The basic idea is simple, we generate a Gaussian random field in Fourier space (assuming a power spectrum), we inverse Fourier transform it to get its analog in configuration space. We further apply a local transform in order to map the Gaussian field into a stochastic field characterized by the target PDF. The two crucial step of this process are the choice of the input power spectrum and the choice of the local transform.

Be a stochastic field following a centered () reduced () Gaussian distribution. From a realization of this field, one can generate a non-Gaussian density field by applying a local mapping between the two, hence

 δ=L(ν). (28)

The local transform must be chosen in order to match some target PDF for the density contrast . Assuming that the local transform is a monotonic function which maps the ensemble into then, due to the probability conservation , the local transform must verify the following matching

 Cδ[δ]=Cν[ν], (29)

where stands for the cumulative probability distribution function. Be the definition assemble of the variable then its cumulative probability distribution function is defined as , where is the PDF of . By definition a probability density function is positive, it follows that its cumulative is a monotonic function and therefore Eq. (29) can always be inverted, it reads

 δ=C−1δ[Cν(ν)],

where the exponent stands for the reciprocal function such that . For example, by definition the local mapping which allows transform a Normal distribution into a Log-Normal distribution is . Note that depending on the PDF that must be matched this inversion can require a numerical evaluation which can be tabulated.

Once a local transform is chosen, we need to adress the question of finding the appropriate power spectrum of the Gaussian field which, once locally mapped into the density field , will match the expected power spectrum. Following Greiner & Enlin (2015), who considered a log-transform we generalized their result to a generic local transformation. This mapping is not direct in Fourier space while it is in configuration. Writing the two point moment of order two of the density field and assuming the probability conservation leads to

 ξδ≡⟨δ1δ2⟩=∫L(ν1)L(ν2)B(ν1,ν2,ξν)dν1dν2, (30)

where is a bivariate Gaussian defined as

 B(ν1,ν2,ξν)≡12π|Cν|1/2exp{−12νTC−1νν}. (31)

One can notice that in our case (central reduced Gaussian) the covariance matrix takes the simple form . Once integrated over the definition domain of and , Eq. (30) provides a mapping between the 2-point correlation function of the Gaussian field and the 2-point correlation function of the density field . However, we prefer to rotate the coordinate system before performing the integral (30) because in case of high correlation () then the gaussian will be comparable with a straight line; most of the sampling of this function will be useless. That’s why we look for the rotation allowing to diagonalize the matrix and therefore convert into a new variable . It follows that

 Cx=[1−ξν001+ξν]

and the integral becomes

 ξδ=12π√1−ξ2ν∫L(x2−x12)L(x2+x12)e−12(x21σ21+x22σ22)dx1dx2, (32)

where and we can therefore integrate over a bounded domain corresponding to the , along the axis and , along the axis. An other possibility to perform the integral 30 is to use the Mehler’s formula, doing so, one can show that the 2-point correlation of the density field can be expressed as a Taylor expansion on the 2-point correlation function of the field. It reads,

 ξδ=λ(ξν)≡∞∑n=0n!c2nξnν, (33)

where the are the coefficients of the Hermit transform of the local mapping and they can be calculated using the orthogonal properties of Hermit polynomials

 (34)

The latter approach considerably speed up the numerical evaluation of Eq. (32), it allows to compute the -D integral as a finite sum of -D integrals. It also allows to verify that when the 2-point function of the field is positive then the derivative of with respect to is positive. Moreover, from Eq.(30) one can see that implies . This means that the function which transforms into is invertible as long as is positive. On the other hand we know that the zero-crossing of the 2-point correlation function occurs at very large scales at which one can safely assume that thus by continuity one can truncate the Eq.(33) at order one providing a linear relation between and . As a result, one can take the reciprocal of the function such that .

Once the local transform and the 2-point correlation mapping are known, then the input power spectrum of the Gaussian field can be obtained as follow. We choose a power spectrum , in the present case Eisenstein & Hu (1998), for the density field , we calculate its corresponding -point correlation function

 ξδ=∫P(k)ei\@veck⋅\@vecrd3\@veck. (35)

At each scale , one can deduce the 2-point correlation function of the Gaussian field and finally using a Fourier transform we obtain the input power spectrum

 Pin(k)=1(2π)3∫ξν(r)e−i\@veck⋅\@vecrd3\@vecr. (36)

Finally, in order to make sure that the PDF target will be reproduced, it is necessary to verify that, once the input power spectrum have been set up on regular k-space grid which will be used to generate the Gaussian field, its integral is indeed equal to the expected variance on the size of the mesh. should be equal to . In general, and are not equal, thus we renormalize the target power spectrum by the quantity , .

We generate a Gaussian field (with a flat power spectrum), on a regular mesh of Mpc and a comoving box of Mpc. We then Fourier transform with an FFT and keep only the phases of the field . We generate at each position the value of the modulus of , where and is a random number with a uniform probability distribution between and . We then inverse Fourier transform the field to get a centered reduced Gaussian field. In Fig. (11) we show the input power spectrum of the Gaussian field compared to the one measured using a FFT, and to the one expected from the local transformation applied to the field in order to generate the density field .

You are adding the first comment!
How to quickly get a good reply:
• Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
• Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
• Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters