Testing the cosmological Poisson equation in a model-independent way

We show how one can test the cosmological Poisson equation by requiring only the validity of three main assumptions: the energy-momentum conservation equations of matter, the equivalence principle, and the cosmological principle. We first point out that one can only measure the combination ${\mathcal M}\equiv \Omega_m^{(0)}\mu$, where $\mu$ quantifies the deviation of the Poisson equation from the standard one and $\Omega_m^{(0)}$ is the fraction of matter density at present. Then we employ a recent model-independent forecast for the growth rate $f(z)$ and the expansion rate $E(z)$ to obtain constraints on ${\mathcal M}$ for a survey that approximates a combination of the Dark Energy Spectroscopic Instrument (DESI) and Euclid. We conclude that a constant ${\mathcal M}$ can be measured with a relative error $\sigma_{\mathcal{M}}=4.5\%$, while if ${\mathcal M}$ is arbitrarily varying in redshift, it can be measured only to within $13.4\%$ (1 $\sigma$ c.l.) at redshift $z=0.9$, and 15-22\% up to $z=1.5$. We also project our constraints on some parametrizations of ${\mathcal M}$ proposed in literature, while still maintaining model-independence for the background expansion, the power spectrum shape, and the non-linear corrections. Generally speaking, as expected, we find much weaker model-independent constraints than found so far for such models. This means that the cosmological Poisson equation remains quite open to various alternative gravity and dark energy models.


Introduction
The Poisson equation is one of the fundamental theoretical construct on which most of cosmology and astrophysics is based upon.Although its validity has been confirmed to strong precision in the solar system and, indirectly, on small astrophysical scales, no model-independent test of it has been so far possible on cosmological scales.In this paper we present a way to achieve this.
When cosmological data are analysed to extract cosmological information, for instance about the cosmic expansion, the perturbation growth rate, or deviations from standard gravity, one usually assumes a set of specific models, for instance extensions or modifications of the ΛCDM standard cosmological model.Unavoidably, therefore, the conclusions one derives are, to some extent, valid only for that particular set, and should be obtained anew when a different theoretical set is chosen.One way to alleviate, but not overcome, this dependence is to encapsulate extensions to ΛCDM in more and more general parameterisations for, e.g. the equation of state (EoS) evolution for dark energy [1], the growth factor as function of the growth index γ [2] or the evolution of the Newtonian and Weyl potential as a function of the matter density contrast [3,4].However these parameterisations as well need to assume a model dependent evolution with time [5,6,7].In recent work ([8, 9, 10]), a way to free the data analysis from this builtin model dependence has been investigated and shown to be able to produce relatively strong constraints of various quantities, including the expansion rate, the distance, the perturbation growth rate, and the ratio of the two linear gravitational metric potentials, the so-called anisotropic stress η.The latter can be constructed by a suitable combination of quantities extracted from the galaxy clustering, the shear lensing, and the supernovae Ia Hubble diagram, without assumptions concerning the power spectrum shape or evolution, or the bias function.
There is another parameter that is often employed to express the deviation from standard gravity, sometimes denoted as G eff or µ = G eff /G N (where G eff = G N denotes a standard Poisson equation in Planck units).This parameter quantifies how much the Poisson equation differs from the standard case, and is the focus of this paper.If dark energy does not cluster on the observed scales -as expected in most models of dark energy (see [11] for a discussion on this assumption) -then µ is indeed a measure of deviation from Einstein's gravity.In a scalar-tensor theory, for instance, µ can be expressed in terms of an extra interaction between the scalar field and matter [12].Several papers tried to constrain µ with various degrees of model-independence.In [13] the authors used principal component analysis methods to constrain µ and the cosmological stress tensor using multiple cosmological probes.However, they also used f σ 8 (z) measurements obtained from galaxy correlations in a model dependent way assuming ΛCDM model for the geometrical distances, a fixed shape for the power spectrum entering the σ 8 and different parameterisations to model the bias for each measurements.
In [14], the authors followed similar methods to forecast bounds on µ and on the dark energy EoS parameter w from future SKA galaxy clustering data.However, they assumed a model for the galaxy bias that is a linear function of the scale factor.Recently, Refs.[15] and [16] constrained µ with current data for the growth f σ 8 (z) from redshift space distortions (RSD) and H(z) from cosmic chronometers.These papers used measurements of the growth f σ 8 (z) obtained assuming a model-dependent shape for the power spectrum; moreover, Ref. [15] included a constraint on the matter density that is also model dependent (we will comment more on this below).Finally, Ref. [7] relaxed σ (0) 8 from its model dependency in an attempt to constrain the growth index γ, a parameter that could be remapped into µ (see [17,18]), but assumed ΛCDM model for the background evolution.In other works, µ is derived from specific models [19,20,21,22] or is parameterized in a phenomenological way [5,23,24,6,18,25].This brief review of some of the previous work shows that so far all tests of the Poisson equation have relied in significant measure on assumptions on the power spectrum shape, on the background evolution, on the growth factor, or on the bias parametrization.In many of these works, tight constraints on the Poisson equation are produced, but they are valid only within such restrictive assumptions.In this paper we present a possible way to overcome these limitations and show that the cosmological constraints on the Poisson equation are much weaker than so far presumed.
Before continuing, it is necessary to be more precise about our definition of model-independence.We aim at measuring (or forecasting, for now) quantities in a way that is independent of the background evolution, of the shape of the power spectrum and of the bias function, and of all the other higher-order bias functions that enter the non-linear perturbation theory (for more details on this, see [10]).The three main underlying assumptions that are needed for this approach are the cosmological principle, the energy-momentum conservation, and the equivalence principle.The cosmological principle ensures statistical homogeneity and isotropy and therefore it allows us to employ the Alcock-Paczyński effect to measure the cosmic expansion rate.The energy-momentum conservation of matter allows us to use redshift-space distortions to measure the growth rate.Finally, the principle of equivalence allows us to use the generalized non-linear kernel derived in [26] and to assume that µ is a universal quantity, that is, the same for dark matter and baryons.We are still nonetheless making an important simplification, namely that the growth rate f is k-independent.This is of course not correct in some modified gravity models or in presence of massive neutrinos.In these cases, our results for f should be considered as an average over a range of scales.However, we remark that in most models the k-dependence is very small: for instance, Ref. [27] has shown that f varies by less than 2.5% in the range k = 0.01 − 0.1h/Mpc.In any case, this limitation can in principle be overcome by employing the same methodology at the price, of course, of a weakening of statistical power.
The main problem of constraining µ, as we discuss in detail below, is that it is degenerate with the matter fraction Ω m (z), the expansion rate E ≡ H(z)/H 0 where H 0 is the Hubble constant at present, and the linear growth rate where δ m is the linear density contrast of matter, a(t) is the cosmological scale factor, and we use primes to denote derivatives with respect to log a.So in order to measure µ one needs to measure Ω m (z), E(z) and f (z) in a modelindependent way.A method to achieve this for E and f has been recently proposed in [10], in which a Fisher matrix forecast was obtained from forthcoming power spectrum and bispectrum data.Since at the moment, as far as we know, there are no proposals to measure Ω m in such a model-independent way (see however [28,29] for an investigation of this aspect), we cannot disentangle µ from it, and therefore we constrain the combination Ω m µ in various redshift bins.Similar degeneracy has also been discussed in [16].If we add the standard hypothesis that matter is pressureless, and therefore its density ρ m scales as a −3 , then we can measure where Ω (0) m is the present value of the matter density fraction.Additionally, given M at two different redshift bins, we can as well take their ratio and obtain µ(z 2 )/µ(z 1 ).If this ratio is different from unity, the standard Poisson equation is violated (of course, µ could be redshift-independent; in this case, a constant ratio would be inconclusive).In this sense, and within the limits of our assumptions, we can state, as advertised in the title of this paper, that our method is indeed a model-independent test of the Poisson equation.
In this paper we present the general idea and produce a forecast for surveys that approximates a combined DESI [30] and Euclid [31] survey.In future work we will apply the method to real data.

Model-independent cosmological observables
We consider the line element for a linearly perturbed Friedmann-Lemaître-Robertson-Walker metric in the longitudinal gauge, with the two scalar degrees of freedom characterized by the gravitational potentials Ψ and Φ, We choose units such that c = G N = 1.Here we do not consider vector and tensor perturbations.Assuming a spatially flat Universe and a pressureless perfect fluid for the matter content, one can obtain the two Poisson equations for Ψ, Φ in Fourier space, from the perturbed Einstein equations [8,32].The gravitational potential Ψ, assuming it is slowly varying, can be mapped out through the equation of linear matter perturbations, where k is the comoving wavenumber, ρ m and δ m the average of the background matter density and the root-meansquare matter density contrast, respectively.The dimensionless effective gravitational constant µ, depending in principle both on time and scale, is the function quantifying modified gravity.In standard General Relativity, at deep subhorizon scales, one has µ = 1.Here we assume universal gravity, i.e. the equivalence principle is preserved.For other consideration see e.g. the work of [33,34].We also assume pressureless, conserved matter, described as a perfect fluid.Then the Euler and continuity equations at deep subhorizon scales (k/(aH) ≫1), in terms of the perturbed variables, are given by, where θ m = ik i v i /aH is the velocity divergence of matter.Combining these two equations we obtain the second order growth equation: Plugging Eq. (4) into Eq.( 8), and replacing ρ m with 3H 2 Ω m /8π, and assuming ρ m scales as (1 + z) 3 , we are now able to extract which shows that M depends only on f and E and their derivatives.From now on we only consider time dependence for M, f and E. Once we have model-independent constraints on f and E, it is straightforward to put constraints on M.Even though Eq. ( 9) has been obtained under the assumption of a flat Universe, it is in general a good approximation also for the non-flat case at sub-horizon scales 1 .Notice that we will not take E(z) from standard candle/ruler measurements, because they measure distances, and E(z) can be inferred from distances only in flat space.
We will first estimate the constraints on M(z) on individual redshift bins for a survey that approximates the specifics of DESI from z < 0.6 and of Euclid in the range z ∈ (0.6 − 2).We call this survey the DE combined survey.The constraints turn out to be relatively weak.Then, since several parametrizations for µ(z) have been proposed in the past, we move to put constraints on them.We consider explicitly three such models.The first, and perhaps the simplest one, arises when assuming a non-minimal constant coupling β between dark energy and matter (see e.g.[8]).In this case in fact The second model is by design sensitive to cosmic acceleration [5,23,24], and we refer to it as "parametrized dark energy" or PDE.Compared to ΛCDM, it has only one extra parameter µ 0 , quantifying the departure from the standard Poisson equation in which µ 0 = 0.In this model µ follows, The third model is the Jordan-Brans-Dicke (JBD) theory of gravity [36].It belongs to a class of general scalartensor theories but it is particularly simple since it depends only on one parameter, w BD .In the GR limit one has w BD → ∞.From a well-known solution of Brans-Dicke gravity in [37,38,39], µ follows a power law: Finally, we also briefly consider the normal-branch of Dvali-Gabadadze-Porrati (nDGP) braneworld gravity ( [40]).In this model matter is confined to a four-dimensional brane embedded in five-dimensional Minkowski spacetime.The cross-over scale r c = G 5 /(2G N ), with G 5 being the Newton's constant in five-dimensional spacetime, is the only extra parameter with respect to GR and characterises the behavior of gravity.On scales above r c gravity is five-dimensional.
In the GR limit one has r c → ∞.It is useful to define the cross-over energy density fraction, In this model, the Poisson parameter µ follows (see e.g.[21]), However, we find later that due to the strong degeneracy between Ω (0) m and Ω rc , the nDGP parameter Ω rc is essentially unconstrained. 1 For a non-flat Universe Eq. ( 4) is replaced by (see e.g.[35]), tot is the curvature fractional density today.Combining this modified Poisson equation with Eq. ( 8), Eq. ( 9) becomes: 3  .
So in non-flat case the quantity one can constrain is However, since at sub-horizon scales k ≫ H 0 , the term 3K/k 2 would be negligible also due to the small Ω (0) k constrained from observations in standard models.

Methodology
In [10] a method, called FreePower, has been presented to estimate some cosmological functions regardless of the cosmological model.The method makes use of the galaxy power spectrum with one-loop correction and of the tree-level bispectrum.We briefly review the main results of Ref. [10] since we are going to use them here.Some more information is provided in Appendix B.
As shown in [26], the one-loop power spectrum and the tree-level bispectrum depend on the linear power spectrum shape, on its growth function, and on two time-dependent, but space-independent, bias functions and five so-called bootstrap parameters, also time-dependent.The bootstrap parameters define general kernel functions for the oneloop correction that depend only on the assumption of the equivalence principle, the non-relativistic limit of the GR diffeomorphism (the so-called extended Galilean invariance [26]), and on the conservation of the matter energymomentum tensor.They do not depend therefore on a specific cosmological model and can be used therefore to parametrize a vast class of cosmologies, including forms of modified gravity.Moreover, we do not rely on a particular shape of the power spectrum (e.g. based on inflationary initial conditions), because in our Fisher matrix approach we take the power spectrum P(k) free to vary in many wavebands.We are also independent of the background cosmology because we leave the expansion rate E(z) and the dimensionless angular-diameter distance H 0 D A (z) free to vary at every redshift.Finally, as already mentioned, we leave the bias and bootstrap parameters free to vary at every redshift.This approach frees ourselves from any specific cosmological model, at least within a large class of models.
For what follows, it is in particular important to stress that we leave the growth rate f (z), the dimensionless expansion rate E(z), and the angular diameter distance D A (z) free.The Alcock-Paczyński (AP) effect, based upon the assumption of the cosmological principle, allows one then to obtain constraints on E(z) and H 0 D A (z).It is important to remark that at linear level one can only get the combination f σ 8 (z) (up to a multiplicative constant), while to constrain f alone the extension to non-linear power spectra and/or bispectra is necessary.In [10] the forecasts for a survey that approximates the Euclid specification [41] assuming as fiducial a standard ΛCDM (the same we assume here) have been obtained.Here we follow with some changes the specifications in [42] that join a DESI-like survey at low redshift to a Euclid-like one at higher redshift, according to Table 1.The DESI-like survey reproduces the specifications for the DESI Bright Galaxy Survey for z ≤ 0.6 based on [30], and the DESI Emission Line Galaxies (ELG) survey for 0.6 ≤ z ≤ 0.9 based on [43], while for 0.9 ≤ z ≤ 2.0 we assume an Euclid-like survey based on [10] and [44].We call this the DE combined survey.
The fiducials for all parameters are based on ΛCDM and we refer to [10] and [42] for more detailed information.Since here we need only the constraints on E and f , we selected from the resulting full parameter covariance matrix only the corresponding entries for n b = 10 redshift bins, with size ∆z = 0.2 and central redshifts z = {0.1,0.3, 0.5, 0.7, 0.9, 1.1, 1.3, 1.5, 1.7, 1.9}.We adopt the "aggressive" case [10], which implies k max = 0.25h/Mpc for the one-loop power spectrum and k max = 0.1h/Mpc for the bispectrum.
Taking advantage of the Fisher matrix containing f i and E i with i = 0, 1, • • • , 9 corresponding to the redshift bin value in the array defined above, our goal is to project this (2×10) 2 = 20 2 matrix onto M to obtain its constraints.From Eq. ( 9), M depends on the derivatives f ′ and E ′ , which can be numerically approximated as and therefore, the discretized values M i are only defined in n b − 2 bins.Accordingly, we can now discretize (9), obtaining M at each bin, where i = 1, 2, • • • , n b − 2. We assume now that M i is distributed as a Gaussian variable.The expansion of M i around its fiducial reads in which for brevity we define ⃗ }, and we use the subscript (F ) to denote the values at the fiducial.Finally, the elements of the covariance matrix are given by (see Appendix A), We choose ΛCDM as our fiducial model, with Ω (0) m = 0.32.In our redshift range, f i and E i for the fiducial cosmology can be approximated with Since the relation between M i and f, E is non-linear, we also checked all the results numerically.That is, we generated 5,000 samples of f, E distributed as multivariate Gaussian variables according to their covariance matrix (inverse of the Fisher matrix) and their mean (the ΛCDM fiducial values), and then obtained the variance of M i from Eq. 15.The numerical results confirm that the Gaussian approximation for M i and the analytical calculation of the next section are very accurate.

Model-independent forecasting on M(z)
Following the method detailed above, we present the forecasted 1σ errors for M (green shades and lines) in Fig. 1 with details on the uncertainty values in each bin for M, E and f in Table 2.We observe that the error for M increases on both sides of the z axis, as a consequence of the same behavior for the errors on f and E of Table 2. To understand this, one has to consider that the Fisher matrix constraints on f and E at linear order are proportional to σ ∼ [(n g P + 1)/Vn g P] −1/2 , where P is the power spectrum.The density n g is small at high redshifts, so n g P ≪ 1, and the product Vn g P is smaller than at around z = 1.Therefore σ gets larger.At lower redshifts, instead, n g P ≫ 1, so σ ∼ V −1/2 , and since the volumes are smaller, σ increases again.So we can expect the best constraints somewhere at intermediate distances.As we see in Table 2, although the errors on E and f are at the percent level, the ones on M are one order of magnitude higher, reaching more than 30% for the farthest bins.We include in our analysis the correlations among the various f 's,E's, and between the f 's and the E's, both within the same bins and among different bins.It is interesting to highlight that these correlations actually help reducing the errors for M by more than six percent.In fact, performing the same analysis but removing the correlations, we find relative errors for M equal to 15.0% and 32.2% in the most tightly constrained and the last bin, respectively.In Fig. 2 we use GetDist [45]  to compare posterior distributions for these two cases.Also, as expected from Eq.15, we find (negative) correlations between M i and M i−2 .
In Fig. 1 we also report our results with priors.Since E is already strongly constrained around the 1% level, we explore how the errors on M would change if we choose strong priors σ (P)  f on f .These could be obtained for example if one manages to improve the higher-order correlators and reach a higher k max in the modeling of large scale structure.We show in the same plot in Fig. 1 the new errors for σ (P)  f equal to 1, 2 or 3% in blue, orange or black lines, respectively.We observe that the errors improve proportionally to the prior.In the best case, namely with a prior σ f = 1%, our constraints on M tighten by ≈ 37% and the relative error in the bin z = 0.9 is reduced to 8.38%. 1 We discuss now the ratios R i j ≡ M i /M j = µ i /µ j .As mentioned in Sec. 1, its deviation from unity indicates directly a violation of the standard Poisson equation, without assumptions on Ω (0) m .Errors on σ i j = σ(R i j ) can be easily estimated using error propagation while also taking the correlations between the nine M i 's into consideration.We visualise their values in Fig. 3.In the best case, we find the constraint on R i j around 20%.Models with timevarying µ would be excluded if their predictions are outside the regions in Table 2 and/or Fig. 3.If one assumes a power law G eff = µ ∼ a α−1 , then one can relate R i j to the exponent α where d i j = ln[(1 + z j )/(1 + z i )], and a dot denotes derivative with respect to the cosmic time.Then the relative error on α is σ α = σ i j /d i j .Therefore upper limits on R i j imply upper limits to the classic test of relativity | Ġeff /G eff |.We will discuss such a power-law G eff in the context of Jordan-Brans-Dicke gravity below.

Modified gravity models
We discuss now some parameterizations of modified gravity that have been proposed in the literature and described at the end of Sect.2.In all cases, we project our nine constraints on M i onto the specific parameters X j by a    transformation Jacobian J i j ≡ ∂M i /∂X j .It is important to remark that although we now choose specific parametrizations for M, we are still general as far as the cosmological model is concerned.Therefore we expect our results will be less strong, but more general, with respect to similar analyses that are carried out within specific models.
In Table 3 and Fig. 4 we present our constraints on the aforementioned models.Except for the constant M case (in which Ω (0) m and µ are fully degenerate, so we only consider the M combination), for each other model, we present our result both varying all the parameters and fixing Ω (0) m to 0.32.In all cases, we find strong correlations between Ω (0) m and the other modified gravity parameters.Notice that we refrain from adopting a Planck prior, or any other similar prior, because they have been obtained for specific models.In the first model, one can assume that M is constant along all the bins.As already mentioned, such a constant M is in fact expected in the simplest models with a coupling constant β between matter and dark energy.A constant M will clearly lead to quite stronger constraints than when M depends freely on redshift.The "new" Fisher matrix F ′ (with only one entry for M) is related to the "old" 8 × 8 Fisher matrix for M i through the transformation, with . We find a relative error σ M = 4.55%.Since M = Ω (0) m (1 + 2β 2 ), the relative error for M corresponds to the absolute error on 2β 2 if we fix Ω (0) m to its Planck value.This corresponds to an absolute constraint σ β 2 ≈ 0.0073 on the coupling β 2 .This constraint is weaker than those from CMB, for instance, Ref. [20] constrained β = 0.0158 +0.0067 −0.0120 , which corresponds to σ β 2 ≈ 0.0003, using CMB data alone and without any prior included (for other constraints on β see e.g.[19]).Neglecting the impact of fiducial values, our error on β 2 is one order of magnitude larger, but complementary, since it is obtained in a completely different redshift range.

Conclusions
In this work we have constrained the deviation from the standard Poisson equation in a model-independent way for a future survey that approximates DESI and Euclid in the range z < 2 and, for the first time, we performed forecasting on the parameter combination M = Ω (0) m µ.Assuming the equivalence principle and the energy-momentum conservation of matter, M depends only on the growth rate f and on the dimensionless expansion rate E, which can be both obtained model-independently from combined power spectrum and bispectrum data.We consider both the cases when M is a free function of redshift, and when it is parametrized in some forms suggested in the literature.
The main message of this paper is that the Poisson equation can be constrained with near-coming cosmological data only up to 13-30%, depending on the redshift, so that there is still considerable space for alternative cosmologies, i.e. models of modified gravity or clustering dark energy.Only by further assumptions, e.g. that M is constant, can better bounds be achieved, although always much broader than usually obtained for specific models.This shows that model-independent analyses are a useful tool to discern and quantify the extent to which the constraints on fundamental cosmological quantities depend on specific assumptions.
On the other hand, one could further tighten the constraints on M by combining measurements from other surveys like SKAO 2 and VRO [53,54,55], as they involve additional survey volume.Moreover, modeling the non-linear galaxy spectra to higher k max will result in more stringent constraints on E and, most importantly, on f .In this way, the uncertainty on M can be significantly lowered.Finally, one could constrain a more general µ by leaving f free in k and z bins and obtain µ(k, z).These developments are left for future work.

Figure 1 :
Figure 1: Reconstructed M = Ω (0) m µ with the 1σ estimated error as a function of redshift, with the blue dashed line the fiducial value of M = 0.32 at each redshift bin, and the purple area assuming a constant M along all the bins.Error bars indicate constraints on M obtained when considering various priors on the growth rate.See the discussion in Sec. 4.

Figure 2 :
Figure 2: Comparison of the one-dimensional posterior distributions of M i , and the corresponding contour plots at 1σ and 2σ confidence levels.We observe a negative correlation for each pair of adjacent M's. z

Table 1 :
Bins' central redshift, survey volume, and galaxy density at each bin of the DE combined survey.