Testing one-loop galaxy bias: cosmological constraints from the power spectrum
Andrea Pezzotta, Martin Crocce, Alexander Eggemeier, Ariel G. Sánchez, Román Scoccimarro
TTesting one-loop galaxy bias: cosmological constraints from the power spectrum
Andrea Pezzotta,
1, 2, 3, ∗ Martin Crocce,
2, 3
Alexander Eggemeier, Ariel G. S´anchez, and Rom´an Scoccimarro Max-Planck-Institut f¨ur extraterrestrische Physik,Postfach 1312, Giessenbachstr., 85741 Garching, Germany Institute of Space Sciences (ICE, CSIC), Campus UAB,Carrer de Can Magrans, s/n, 08193 Barcelona, Spain Institut d’Estudis Espacials de Catalunya (IEEC), 08034 Barcelona, Spain Institute for Computational Cosmology, Department of Physics,Durham University, South Road, Durham DH1 3LE, United Kingdom Center for Cosmology and Particle Physics, Department of Physics,New York University, NY 10003, New York, USA (Dated: February 17, 2021)We investigate the impact of different assumptions in the modeling of one-loop galaxy bias onthe recovery of cosmological parameters, as a follow up of the analysis done in the first paper ofthe series at fixed cosmology. We use three different synthetic galaxy samples whose clusteringproperties match the ones of the CMASS and LOWZ catalogues of BOSS and the SDSS MainGalaxy Sample. We investigate the relevance of allowing for either short range non-locality or scale-dependent stochasticity by fitting the real-space galaxy auto power spectrum or the combinationof galaxy-galaxy and galaxy-matter power spectrum. From a comparison among the goodness-of-fit ( χ ), unbiasedness of cosmological parameters (FoB), and figure-of-merit (FoM), we find that afour-parameter model (linear, quadratic, cubic non-local bias, and constant shot-noise) with fixedquadratic tidal bias provides a robust modelling choice for the auto power spectrum of the threesamples, up to k max = 0 . h Mpc − and for an effective volume of 6 h − Gpc . Instead, a joint analy-sis of the two observables fails at larger scales, and a model extension with either higher derivativesor scale-dependent shot-noise is necessary to reach a similar k max , with the latter providing themost stable results. These findings are obtained with three, either hybrid or perturbative, pre-scriptions for the matter power spectrum, RESPRESSO , gRPT and EFT. In all cases, the inclusionof scale-dependent shot-noise increases the range of validity of the model in terms of FoB and χ .Interestingly, these model extensions with additional free parameters do not necessarily lead to anincrease in the maximally achievable FoM for the cosmological parameters (cid:0) h, Ω c h , A s (cid:1) , which aregenerally consistent to those of the simpler model at smaller k max . I. INTRODUCTION
Over the past decades galaxy redshift surveys have pro-vided a wealth of information on the large-scale distribu-tion of galaxies across the Universe. Clustering measure-ments of two-point statistics – the galaxy power spectrum P ( k ) and the two-point correlation function ξ ( s ) – fromlarge data samples can indeed provide precise measure-ments about the underlying cosmological model [1–4].The inference of cosmological parameters from large-scale structure is made intrinsically more difficult by therealisation that galaxies are a biased tracer of the to-tal matter density field [5–11]. This translates into thenecessity of having an accurate description of the rela-tionship between the galaxy and matter density fields,a phenomenon which is commonly referred to as galaxybias . Even though the latter can be partially understoodin a phenomenological way, e.g. using results from N-body simulations, the complexity of dealing with tracersfeaturing different morphological properties makes desir-able to develop an analytic formulation that is based ona more theoretical background. From this point of view, ∗ [email protected] perturbative approaches stand as a natural way of de-scribing galaxy bias in a physically motivated way.The main idea behind this formulation is that thegalaxy overdensity δ g can be expressed in terms of a se-ries of operators involving spatial derivatives of the grav-itational and velocity potentials. At leading order, thisrelation is captured by a single multiplicative factor, i.e. δ g = b δ , where δ is the matter density contrast and b isa multiplicative factor called linear bias [6]. Higher-ordercontributions become progressively more important onmildly non-linear scales, as expected from a spherically-symmetric gravitational collapse [12, 13], in a way thatthe expression for δ g can be expanded to include higherpowers of δ , i.e. δ g = (cid:80) n b n /n ! δ n [6, 14]. At the sametime it has been shown that anisotropies in the processof gravitational collapse are responsible for the genera-tion of non-negligible tidal effects, which also contributeto the local distribution of galaxies [15, 16]. This findingfollowed the realization that the local-in-matter-densitybias model had limitations in providing a proper descrip-tion of the clustering of dark matter halos [17, 18] andwas leading to incompatible constraints on the quadraticbias b from measurements of the power spectrum andbispectrum [19, 20].The most important aspect of the previously describedmodel is that galaxy bias is treated as a spatially local a r X i v : . [ a s t r o - ph . C O ] F e b quantity. However, it is well known that the formation ofhalos and galaxies is triggered by the gravitational col-lapse of matter from a spatially finite region, and there-fore the local assumption is bound to fail when approach-ing scales that roughly correspond to the Lagrangian ra-dius R of the host halos. In terms of δ g , this effect canbe taken into account by considering not only its depen-dency on density and tidal fields, but also on functionalsof δ [11]. This is equivalent to introducing higher deriva-tives of the matter density field, that at leading order pro-vide a contribution to δ g of the form R ∇ δ [21–23]. Afurther ingredient to the galaxy-matter relation is repre-sented by stochastic terms, that on large scales behave asan additional contribution to Poissonian shot-noise [24].Stochasticity is the direct result of small-scale perturba-tions, which are not correlated over long distances underthe assumption of Gaussian initial conditions [25–27]. Athigher wave modes the halo-exclusion effect [12, 28, 29]imprints a scale dependence on the stochastic contribu-tions, whose strength is controlled by the Lagrangianradius [30], similarly to higher derivatives. At next-to-leading order, this contribution scales as k , and mightbecome relevant even for clustering analyses based ontwo-point statistics, as shown in [31].Galaxy bias models based on the one-loop perturba-tive expansion have been used to extract cosmologicalconstraints from big data collaborations, such as BOSS,using clustering measurements both in configuration [32]and in Fourier space [33–35]. More recently, the samedata have been re-analysed in a number of works [36–38]with novel techniques that nevertheless assume the samebiasing scheme for the galaxy-matter relationship. Themajority of these analyses assume a different ansatz interms of the degrees of freedom of the model, by fixingsome of the free bias parameters to some physically mo-tivated relations. For example [32][35] and [38] fix thecubic non-local bias parameter to the Local Lagrangian (LL) relation, while [34] also fix the quadratic non-localparameter. On the other side, [36] set the cubic term tozero and leave the quadratic tidal bias completely free.This paper is part of a set of works [31][39] in whichwe explore the impact of different bias modeling choicesusing a set of three performance metrics, namely thegoodness-of-fit, the unbiasedness of sampling parameters,and the merit of the model. To do that, we make use ofa set of three different simulated galaxy samples with aneffective volume of 6 h − Gpc , whose clustering prop-erties and number densities reproduce the ones of threereal catalogues, the CMASS and LOWZ samples of BOSS[40–42], and the Main Galaxy Sample of SDSS [43]. Inorder to exclusively concentrate on the modeling of one-loop galaxy bias, we carry out this analysis in real-spaceremoving the impact of redshift-space distortions, whichwould require an additional modeling layer. In [31] weperformed a fixed cosmology analysis by using the mea-sured non-linear matter power spectrum. In that case, weused the linear bias parameter as a proxy for the goodnessof our model. Here, we additionally model the impact of matter non-linear evolution and sample over cosmologicalparameters to determine the impact of one-loop galaxybias on the recovery of such parameters. In the first partof the paper we use the hybrid perturbative-simulatedapproach, RESPRESSO , while comparing results for differ-ent PT-based dark matter models in a later section.Our paper is organised as follows. In Section II wesummarise the main ingredients of our galaxy bias modelincluding a description of all the terms contributing atone-loop in perturbation theory. This includes also ashort review of the three non-linear matter predictionswe employ in this work. In Section III we describe thesimulated galaxy samples we use to test one-loop galaxybias, along with fitting procedure, parameter priors andthe performance metrics we introduced above. Main re-sults from fit of the auto galaxy power spectrum and thecombination between the latter and the galaxy-mattercross power spectrum are given in Section IV. We finallydraw our conclusions in Section V.
II. ONE-LOOP PERTURBATION THEORY FORBIASED TRACERS
The theory describing the evolution of the clustering ofbiased tracers on mildly non-linear scales is a well estab-lished topic (for a review of galaxy bias see [11]) that canbe naturally described in the framework of perturbationtheory. This section stands as an overview of the mostimportant results obtained at one-loop in standard per-turbation theory (SPT) and derived approaches. Noticethat for sake of readability we omit all the dependenceson redshift z . A. Galaxy bias expansion
The general perturbative expansion of galaxy bias canbe interpreted as a sum of different operators that arefunctions of the gravitational potential Φ and velocitypotential Φ v . Focusing on one-loop contributions, therelationship between biased tracers and the underlyingdark matter density field can be described consideringterms up to third order in the matter perturbation δ .Following the notation of [44] we can write δ g ( x ) = b δ ( x ) + b δ ( x ) + γ G (Φ v | x )+ γ G ( ϕ , ϕ | x ) + β ∇ δ ( x ) + . . . , (1)where the first two terms on the right-hand side are partof the standard local expansion in powers of δ .In the previous equation G is a Galilean invariant op-erator, representing the tidal stress tensor generated bythe velocity potential Φ v , and it is given by G (Φ v ) = ( ∇ ij Φ v ) − ( ∇ Φ v ) . (2)In Fourier space this translates into G ( k ) = (cid:90) q (cid:20) q · ( k − q ) q | k − q | − (cid:21) θ ( q ) θ ( k − q ) , (3)where θ is the divergence of the matter velocity field, suchthat ∇ Φ v ≡ θ .Differently from the first two terms of Equation (1),terms involving G incorporate non-local (in matter den-sity δ ) contributions that spontaneously arise from thenon-linear evolution of the matter density field. The sec-ond of these terms is obtained by expressing the non-linear velocity potential up to second order (i.e. Φ v =Φ (1)v + Φ (2)v , ϕ = − Φ (1)v , ϕ = − Φ (2)v ), and keeping thenext-to-leading order correction [45], leading to G ( ϕ , ϕ ) = ∇ ij ϕ ∇ ij ϕ − ∇ ϕ ∇ ϕ , (4)where ∇ ϕ = − θ is the linear velocity divergence field,and ∇ ϕ = −G ( ϕ ) is the next-to-leading order. Thenet result of adding this higher-order correction is that δ g collects contributions up to third order in the matterperturbation δ .In addition to the standard expansion, any biasedtracer also comes with a physical scale which regulatesthe importance of higher-derivative operators [8, 21, 23,27, 46], whose leading order scales as ∇ δ . This scalequantifies the size of the region in which galaxy formationoccurs, and it is therefore influenced by any short-rangegravitational effect and baryonic corrections. For halos,this scale is close to the Lagrangian radius [21, 47], whileit differs for other kinds of tracers, such as galaxies andquasars, depending on their type.An important aspect of the previously described bias-ing scheme is the renormalization of the parameters onwhich it is based. As a matter of fact, the bare bias pa-rameters listed in Equation (1) are sensitive to the UVcutoff scale used to make one-loop integrals convergent.In particular, the current basis carries a dependence onthe variance of the matter density field σ = (cid:104) δ ( x ) (cid:105) ,and more generally to any higher-order correlator (cid:104) δ ( n ) (cid:105) .This dependence is completely non-physical and can bereabsorbed by means of an adequate renormalization ofthe bias parameters [48, 49]. The new basis (for con-vention this is denoted without the overscript hat) canbe identified using the peak-background split formalism[50], and it is constructed in a way that its componentsquantify the response of the cosmic mean abundance oftracers to a change in the background density, with nodependence on the one-loop cutoff scale.An alternative approach to renormalization can be ob-tained expanding δ g in terms of the galaxy multi-pointpropagators , a topic which was firstly described in [44]. We adopt the usual convention for the Fourier transform, δ ( x ) = (cid:82) k exp ( − i k · x ) δ ( k ), and use the short-hand notation for the 3dintegrals, (cid:82) k ,..., k n ≡ (cid:82) d k / (2 π ) · · · d k n / (2 π ) . Consistently with the naming convention adopted in thecontext of renormalized perturbation theory [51], thesequantities are defined as the ensemble averaged deriva-tives of δ g with respect to δ L , (cid:28) ∂δ g ( k ) ∂δ L ( k ) · · · ∂δ L ( k n ) (cid:29) ≡ (2 π ) Γ ( n )g ( k , . . . , k n ) × δ D ( k − k ··· n ) , (5)and, being observables themselves, they are already nor-malized by construction. For this reason, the multi-pointpropagators can be identified with the scale-dependentbias parameters in a way that δ g = Γ (1)g ⊗ H + Γ (2)g ⊗ H + . . . , (6)where the H n are the Wiener-Hermite functionals [44,52], and the operator ⊗ is defined as (cid:104) Γ ( n )g ⊗ H n (cid:105) ( k ) ≡ (2 π ) (cid:90) k ,..., k n δ D ( k − k ··· n ) × Γ ( n )g ( k , . . . , k n ) H n ( k , . . . , k n ) , (7)Using this approach, we can define our observables ina fully consistent framework. In this paper we are inter-ested in the galaxy auto power spectrum P gg and galaxy-matter cross power spectrum P gm , which are defined asthe auto and cross-correlation of the two density fields δ g and δ m , as (cid:104) δ g ( k ) δ g ( k (cid:48) ) (cid:105) ≡ (2 π ) P gg ( k ) δ D ( k + k (cid:48) ) , (8) (cid:104) δ g ( k ) δ m ( k (cid:48) ) (cid:105) ≡ (2 π ) P gm ( k ) δ D ( k + k (cid:48) ) . (9)Substituting Equation (1) in the previous set, we get tothe full expressions for the galaxy auto and cross powerspectra at one-loop, P gg ( k ) = b P mm ( k ) + b b P b b ( k ) + b γ P b γ ( k )+ b P b b ( k ) + b γ P b γ ( k ) + γ P γ γ ( k )+ b γ P b γ ( k ) − b β k P L ( k ) , (10) P gm ( k ) = b P mm ( k ) + 12 (cid:20) b P b b ( k ) + γ P b γ ( k )+ γ P b γ ( k ) (cid:21) − β k P L ( k ) , (11)where P L is the linear matter power spectrum, P mm isthe non-linear matter power spectrum, and the one-loopbias corrections read P b b ( k ) = 2 (cid:90) q F ( k − q , q ) P L ( | k − q | ) P L ( q ) , (12) P b γ ( k ) = P mc b γ ( k ) + P prop b γ ( k )= 4 (cid:90) q F ( k − q , q ) S ( k − q , q ) P L ( | k − q | ) P L ( q )+ 8 P L ( k ) (cid:90) q G ( k , q ) S ( k − q , q ) P L ( q ) , (13) P b b ( k ) = 12 (cid:90) q (cid:2) P L ( | k − q | ) P L ( q ) − P ( q ) (cid:3) , (14) P b γ ( k ) = 2 (cid:90) q S ( k − q , q ) P L ( | k − q | ) P L ( q ) , (15) P γ γ ( k ) = 2 (cid:90) q S ( k − q , q ) P L ( | k − q | ) P L ( q ) , (16) P b γ ( k ) = 4 P L ( k ) (cid:90) q S ( k − q , q ) S ( k , q ) P L ( q ) . (17)Here F and G are the symmetrised second-order mode-coupling kernels [52], and S ( k , k ) = (cid:16) ˆ k · ˆ k (cid:17) − G (Φ v ), asshown in Equation (3).The only two propagator-like contributions are per-fectly degenerate with each other, and follow the relation P b γ ( k ) = − P prop b γ ( k ) . (18)For this reason, in a real analysis, it is common prac-tice to either neglect one of the two tidal field relatedparameters (e.g. [36]) or to assume perfectly local-in-matter-density relations (see Section II D) to express oneof them in terms of lower order local bias parameters (e.g.[32, 53]).Since the P b b contribution does not asymptote to 0in the large-scale limit, we renormalise it as in Equation(14), and absorb the constant low- k amplitude as an ad-ditional contribution to the shot noise error, that will bediscussed in Section II C. B. Matter modeling
In this section we discuss the modeling options for thenon-linear matter power spectrum P mm in Equations (10)and (11). As a matter of fact an accurate modeling of P mm is essential, as any systematic effect in the descrip-tion of the matter density field in the range of scales weare considering might lead to invalid interpretations ofthe galaxy-matter bias relationship. Differently from [31], where one-loop bias was inves-tigated at fixed cosmology and adopting the measuredmatter power spectrum as reference, here the main goalis to assess the level of accuracy of our model in termsof cosmological parameters. For this reason, we have toexplicitly assume a model for P mm . In the rest of thissection we provide a description of the three models weare going to test. In particular we use 1) a refined RPT-derived model, based on the preservation of Galilean in-variance (dubbed gRPT), 2) an EFT-like approach basedon BAO damping and a non-trivial stress tensor, and3) a mixed approach, RESPRESSO , based on accurate N-body simulations and a perturbative expansion aroundthe fiducial cosmology.
1. Standard perturbation theory
The basic assumption of SPT is that dark matter be-haves as a perfect pressureless fluid on large enoughscales, where matter is not subject to shell crossing asin multi-streaming regions. Under this assumption, andafter having expanded the matter density contrast in aTaylor series, i.e. δ = δ (1) + δ (2) + δ (3) + . . . , we findsolutions at every order in perturbation as [52] δ ( n ) ( k ) = (cid:90) q ... q n δ D ( k − q − . . . − q n ) F n ( q , . . . , q n ) × δ L ( q ) . . . δ L ( q n ) , (19)where F n is the n -th order symmetrized kernel describ-ing the non-linear mode coupling between fluctuations atdifferent wavelengths.Moving to two-point statistics, the expansion for thematter power spectrum can be written as P mm ( k ) = P L ( k ) + P - loop ( k ) + . . . , (20)where at one-loop the only non-vanishing contributionsare P - loop ( k ) = P ( k ) + P ( k )= 2 (cid:90) q F ( k − q , q ) P L ( k − q ) P L ( q )+ 6 P L ( k ) (cid:90) q F ( q , − q , k ) P L ( q ) . (21)It is now well established that a SPT approach like theone described above lead to significant residuals if com-pared to the output of numerical simulations, even in-cluding higher-order corrections [54]. The source of thisinaccuracy can be mostly identified in two separate ef-fects, whose description is the subject of the next twosections.
2. BAO damping from large-scale “infrared” modes
One of the most acknowledged deviations between one-loop SPT predictions and the matter power spectrummeasured from numerical simulations is the shape of theBAO features. Since the characteristic scale of the BAOpeak is much larger than the scale at which non-linearcontributions become important, one may expect a stan-dard perturbative approach to provide accurate predic-tions on that scale. However, large-scale bulk motionsproduce a non-negligible effect on the amplitude of thepower spectrum at the BAO scale, the most significant ofwhich is a smearing of the BAO signal due to the large-scale relative displacement field [55–59].These corrections to the matter power spectrum werefirstly resummed in the context of RPT [55], and at lead-ing order the net effect is to apply a damping factor tothe BAO wiggles. Practically, we can express the matterpower spectrum as the sum of a smooth ( P nw ) and wiggly( P w ) term [60], P ( k ) = P nw ( k ) + P w ( k ) . (22)The smooth-wiggle split can be realised adopting sev-eral different recipes, and throughout this paper we fol-low the approach described in [61] and [62], where thesmooth component is defined as a rescaling of the fea-tureless spectrum firstly defined in [63] to account forbroadband difference with the linear power spectrum.At leading order, the damping factor can be calcu-lated assuming the Zel’dovich approximation, and subse-quently applied to the wiggly component, so that P LO ( k ) = P nw ( k ) + e − k Σ P w ( k ) , (23)whereΣ = (cid:90) k S P nw ( q ) (cid:20) − j (cid:18) qk BAO (cid:19) + 2 j (cid:18) qk BAO (cid:19)(cid:21) dq π , (24)is the relative displacement field two-point function at theBAO scale [56]. Here j n is the n -th order spherical Besselfunction, k BAO = π/l BAO is the wavemode correspond-ing to the reference BAO scale l BAO = 110 h − Mpc, and k S is the UV limit of integration. To properly accountfor the resummation of IR modes at any given scale k ,one should integrate over all modes q < k , in a way that k S = k S ( k ). However, pushing the integration to signif-icantly large values of k S would result in the breakingof the range of validity of the pertubative IR expansion.At the same time, it can be shown that the integrand ofEquation (24) gives significant contributions only up to k S ∼ . h Mpc − and therefore we restrict the integra-tion to the range [0 , .
2] when computing the value of Σ[36, 64].At next-to-leading order, the IR-resummed matterpower spectrum can be written as [64] P NLO ( k ) = P nw ( k ) + (cid:0) k Σ (cid:1) e − k Σ P w ( k )+ (cid:16) P - loopnw ( k ) + e − k Σ P - loopw ( k ) (cid:17) , (25) where P - loopnw is the one-loop matter correction definedin Equation (21) but evaluated using the smooth compo-nent P nw rather than the full linear power spectrum P L ,and P - loopw = P - loop − P - loopnw .
3. Small-scale corrections: non-trivial stress tensor
Along with the resummation of infrared modes, wealso have to consider the impact of small-scale physicson long wavelength fluctuations. This happens becausethe original assumption of a perfectly pressureless fluidis bound to fail on non-linear scales, where dark matterexperiences shell-crossing in multistreaming regions [52].Moreover, the effect of baryonic physics, such as galaxyformation, cooling and feedback, also contributes in gen-erating a baryonic pressure that impacts the clusteringof dark matter on larger scales.The net effect of UV scales on dark matter clusteringis to generate a non-zero stress tensor [65] whose leadingcontribution to the matter power spectrum is to add acounter-term of the form [65–67] P ctr ( k ) = − c k P L ( k ) . (26)Here c can be treated as an effective speed of sound, thatreflects the influence of short wavelength perturbations,and in particular of the complex physics beyond galaxyformation. Given the poor knowledge about these typesof processes, a standard assumption is to treat c as afree parameter (see e.g. [68, 69]) and marginalise over itto obtain the posterior distribution of the parameters ofinterest.By inspection of the individual terms contributing toone-loop formulas in Equations (10) and (11), we cannotice that the counter-term defined in Equation (26)is completely degenerate with the leading order higher-derivative contribution defined in Section II A, as theyboth scale as k P L ( k ). In principle we may remove thisdegeneracy when jointly fitting the auto and cross galaxypower spectra, as the non-vanishing stress tensor affectsonly the one-loop matter power spectrum P mm . Thisresults in different imprints on P gg and P gm , P gg ( k ) ⊃ − b ( b c + β ) k P L ( k ) ,P gm ( k ) ⊃ − (2 b c + β ) k P L ( k ) . (27)In practice, the sensitivity to UV modes might affect indifferent ways P gg and P gm , leading to inconsistent valuesof c between the two observables. For this reason inthe rest of the paper we will employ two independentfree parameters β P and β × P characterising the k P L ( k )contributions coming from P gg and P gm , respectively.
4. Modeling of the non-linear matter power spectrum
In this section we briefly summarise the three differentprescriptions we adopt to model the non-linear matterpower spectrum P mm throughout the rest of the paper.The first of such models is based on a perturbative-simulated mixed approach, which revolves around a fidu-cial high-resolution measurement of the non-linear mat-ter power spectrum from N-body simulations, and a two-loop perturbative expansion in the cosmological param-eter space for the response function [70, 71]. The latterquantifies the variation of the non-linear power spectrumat scale k induced by a variation of the linear power spec-trum at scale q , namely K ( k, q ) ≡ q ∂P mm ( k ) ∂P L ( q ) . (28)Under the assumption of having a reliable measurementof P mm for a fiducial cosmology θ fid , it is possible topredict the same observable at a generic position θ as P mm ( k | θ ) = P mm ( k | θ fid ) (cid:90) d (log q ) K ( k, q ) × [ P L ( q | θ ) − P L ( q | θ fid )] . (29)The range of validity of this mixed approach becomesprogressively less accurate for cosmologies that are farway from θ fid , but this issue can be overcome by employ-ing a multi-step reconstruction starting from the fiducialcosmology.The RESPRESSO public package [72] makes use of thisapproach, starting from a fiducial measurements of P mm from a set of high-resolution N-body simulations with thePlanck 2015 cosmology [73].The second model we consider is based on the EffectiveField Theory of Large Scale Structure [66, 74], and it isclose to what was recently used in the full shape analysisof the BOSS DR12 galaxy power spectrum [36]. At one-loop in real-space, this model is based on SPT results, butit also accounts for the effect of IR and UV modes on theevolution of the matter power spectrum, as described inthe previous two sections. Namely, we can write a simpleexpression for the one-loop matter power spectrum, thatreads P mm ( k ) = P NLO ( k ) + P ctr ( k ) , (30)where P NLO ( k ) and P ctr ( k ) are defined in Equations (25)and (26), respectively.The third model we consider is based on a particularflavour of Renormalised Perturbation Theory (RPT) [51,55]. In this kind of approach, the non-linear matter powerspectrum is typically separated into a component thatevolves the initial density contrast independently at eachwavelength, called propagator G ( k ), and a mode-couplingterm that accounts for the mixing of scales due to non-linear evolution, so that we can write P ( k ) = G ( k ) P L ( k ) + P MC ( k ) . (31)In a RPT-based approach (e.g. [75, 76]), the propagatoris resummed while keeping the mode-coupling term at afixed order, leading to a breaking of the Galilean invari-ance (GI) of equal-time correlators [77]. This translates into an unphysical damping of the broadband power,which becomes mostly significant in the UV regime. TheRPT flavour we consider here, known as gRPT, effec-tively attempts to resum the mode-coupling term in away that is consistent with the resummation of the prop-agator (see equation 17 of [31] for the explicit formula).At the same time, this approach naturally incorporatesIR resummation consistently with what is described inSection II B 3. We notice that, although even in this casethe impact of the non-zero stress tensor should be takeninto account with the addition of an effective speed ofsound, the broadband predicted by gRPT is slightly sup-pressed with respect to SPT predictions. This would inprinciple lead to an even smaller UV counter-term, andfor this reason we fix c = 0 when modeling P mm withgRPT.Being partially calibrated on numerical measurements, RESPRESSO provides much more accurate results on non-linear scales, as it intrinsically incorporates higher-ordercorrections with respect to the previous two models,which include only one-loop contributions to the mat-ter density field. For this reason, in the first part of thispaper we fix the description of dark matter non-linearevolution following this approach, and check the impactof using different matter models only in Section IV D.This choice is also motivated by the analysis performed in[31], that showed how
RESPRESSO is the model that mostclosely reproduces the performances of the true matterpower spectrum measured from simulations.
C. Stochasticity
One ingredient that is still missing from Equation (10)is the stochastic contribution to the galaxy power spec-trum. As a matter of fact, the previous relations arecompletely deterministic, and assume a one-to-one cor-respondence between the distribution of galaxies and thecombined effect of the matter density and tidal fields.However, galaxy formation is determined not only bylarge-scales perturbations, but also by short wavelengthmodes. Under the assumption of Gaussian initial con-ditions, these modes are completely uncorrelated fromlarge-scales fluctuations, and give birth to an additionalstochastic field ε g which is dependent on the local distri-bution of matter.In practice, we can define the stochastic contribution P ε g ε g to the power spectrum as [11] (cid:104) ε g ( k ) ε g ( k (cid:48) ) (cid:105) = (2 π ) P ε g ε g ( k ) δ D ( k + k (cid:48) ) , (32)and add this contribution to the one-loop galaxy powerspectrum in Equation (10). As for the modeling of P ε g ε g ,we notice that the relation between the distribution ofgalaxies and the high- k modes is not exactly local, asit depends on the distribution of matter within a smallfinite region, similarly to the higher-derivatives of thegalaxy density field. For this reason we can write P ε g ε g ( k ) = 1 n (cid:0) N + N k + . . . (cid:1) , (33)where n is the mean number density of galaxies in theconsidered volume. The constant contribution N repre-sents deviations from purely Poissonian shot-noise (1 /n ),while higher-order corrections, of which the leading termscales as k , are generated to account for the short-rangenon-locality described above. N is expected from the halo exclusion effect [28–30]for which two different dark matter halos cannot overlap(the same principle is what drives the effective matterpressure in multi-streaming regions). This implies a de-viation from Poissonian shot-noise that can be either pos-itive (super-Poisson) or negative (sub-Poisson), depend-ing on the considered tracer. Sub-Poissonian shot-noiseis more expected for central galaxies of massive halos,while super-Poissonian values are more typical of galaxypopulations with high satellite fractions [30, 78].Notice that assuming Gaussian initial conditions allcross-correlators of the form (cid:104) ε g δ (cid:105) are null by construc-tion, and therefore Equation (33) is the only stochas-tic contribution to the galaxy auto power spectrum P gg .In reality, non-linear gravitational evolution introduces adegree of correlation between long and short wavelengths,so that (cid:104) ε g δ (cid:105) (cid:54) = 0 at later times, but all of these contribu-tions are subdominant in the case of the power spectrum,and can thus be neglected.As anticipated in Section II A, the quadratic term P b b needs to be renormalized in order to provide a null con-tribution in the low- k limit. For this reason, we subtractfrom Equation (14) its large-scale asymptote [21, 48], de-fined by P noise b b = b (cid:90) q P ( q ) , (34)and reabsorb it into the constant shot-noise parameter N .Along with P gg it can be shown that also P gm re-quires an additional stochastic component. In this case,stochasticity is not sourced by ε g , as once again all corre-lators of the form (cid:104) ε g δ (cid:105) vanish, but rather by the matterdensity field itself via a new stochastic field ε m . Thereason of this is that dark matter ceases to behave as anideal pressureless fluid on small-scales, where the dynam-ics of gravitational collapse is subject to shell crossing.This translates into an effective pressure exerted by thematter density field, whose contribution to the galaxy-matter cross power spectrum scales as k in the low- k limit [79]. For this reason, it follows [11] (cid:104) ε g ( k ) ε m ( k (cid:48) ) (cid:105) = (2 π ) P ε g ε m ( k ) δ D ( k + k (cid:48) ) , (35)where P ε g ε m ( k ) = 1 n ( N × k + . . . ) . (36) D. Co-evolution relations
Although the previously described galaxy bias expan-sion is complete at one-loop, the large number of freebias parameters makes its applicability to real datasetsdifficult, particularly when using the information con-tained only in the two-point statistics (using additionalconstraints from e.g. the galaxy bispectrum can breaksome of the degeneracies between parameters). For thisreason it is common practice to adopt empirical relationsamong bias parameters in order to reduce the degrees offreedom of our model.The most natural expressions are the so-called localLagrangian relations [21, 80–82]. The latter are based onthe assumption that galaxies (or more generally any bi-ased tracer of the matter density field) are formed instan-taneously at an infinite past time, that galaxy formationis driven exclusively by the local matter density field, andthat the number of tracers is conserved after their forma-tion. Under these assumptions, it is possible to describethe subsequent evolution of the galaxy density field un-der the effect of gravity, which leads to the appearanceof higher-order non-local operators even in the presenceof a purely local relation at the time of formation.It is possible to relax the Lagrangian local-in-matter-density assumption, while still requiring the conservationof the total number of tracers. In this way, we can ex-press the Eulerian non-local parameters γ and γ as afunction of the corresponding Lagrangian counterparts, γ = −
27 ( b −
1) + γ , L , (37) γ = −
221 ( b −
1) + 67 γ + γ , L , (38)where the subscript L indicates Lagrangian quantities,and the remaining terms in the right-hand sides are theresult of gravitational evolution [15, 16, 44].Although local Lagrangian relations (i.e. γ , L = γ , L = 0) have been proven to be much more accuratethan simply neglecting non-local terms, recent measure-ments of non-local bias parameters from a wide range ofhalo samples showed a slight deviation with respect toobservations [83, 84]. An alternative approach for thequadratic tidal parameter is based on the excursion settheory [85], for which γ can be predicted as a functionof b using a quadratic fit, γ , ex ( b ) = 0 . − . b + 0 . b . (39)As shown in Figure 1 of [31], this relation is more accuratethan the local Lagrangian approximation for tracers with b (cid:38) .
3. Therefore this relation should apply much bet-ter to our datasets, which show a linear bias consistentlylarger than this value (see Table I). For this reason, wefix γ to Equation (39) throughout the rest of the paper. III. DATA AND METHODOLOGYA. Simulated galaxy samples
The robustness of our biasing scheme can be assessedonly by validating the model over a wide range of trac-ers, featuring different host halo masses, galaxy bias, andredshifts. For this purpose, we make use of a set of threedifferent synthetic galaxy samples, whose main proper-ties are summarised in Table I.The three catalogues were generated by populatingdark matter halos with galaxies using Halo OccupationDistribution (HOD) prescriptions. Since the HOD pa-rameters were calibrated to obtain number densities con-sistent with the ones of pre-existing real observations, forsake of easiness we label our simulated catalogues withthe name of the corresponding data samples. Neverthe-less, we remind the reader that here we make use of thefull comoving snapshot volume, without considering anyselection effect, as the goal is to investigate the range ofvalidity of 1-loop galaxy bias, leaving aside the impact ofobservational systematics.The CMASS sample is based on the
Minerva simula-tions [53], which consist in a set of 100 realizations of a(1500 h − Mpc) cubic box with periodic boundary con-ditions. The simulations were run using the public Gad-get code [86, 87], that regulated the motion of 1000 dark matter particles within the aforementioned volume.Initial conditions were set up using the linear power spec-trum obtained from CAMB [88] and displacing particlesaccording to second-order Lagrangian Perturbation The-ory (2LPT) [89]. The comoving snapshot has a redshiftof z = 0 .
57 and it is meant to reproduce the propertiesof the BOSS CMASS galaxy sample [90].The LOWZ ( z = 0 . z = 0 . LasDamas
N-body simulations [91, 92]. Theseare a set of 40 independent realizations with periodicboundary conditions, with a volume of (2400 h − Mpc) and (1000 h − Mpc) , and a mass resolution of 4 . × h − M (cid:12) and 4 . × h − M (cid:12) , respectively. Ini-tial conditions were also set up using 2LPT, but in thiscase the initial power spectrum is computed using CMB-Fast [93]. The HOD parameters of these synthetic cat-alogues were selected to reproduce the properties of theBOSS LOWZ and SDSS Main Galaxy Sample (MGS) at M r < − Min-erva and
LasDamas is summarised in Table II.
B. Measurements of power spectra and theircovariances
We make use of the estimator described in [94], basedon fourth order particle assignment scheme and interlac-ing optimization to reduce the effect of large modes alias-ing, to measure the galaxy auto power spectrum P gg and the galaxy-matter cross power spectrum P gm for all thetracers described in Table I. We adopt a linear k binningfrom k min = k f , where k f = 2 π/L is the fundamentalfrequency of the box of size L , to k max = k Nyq , where k Nyq = πN grid /L is the Nyquist frequency correspondingto a FFT grid of size N grid . We use a linear binning withstep ∆ k = k f for CMASS and LOWZ, and ∆ k = 2 k f forMGS. Since we model the stochastic contribution to thegalaxy auto power spectrum in terms of deviations fromPoisson shot-noise, we correct our raw measurements of P gg by subtracting the constant factor P noise = 1 /n .Our final data vector consists in the ensemble averageover the number of independent realizations N R (100 forCMASS, 40 for both LOWZ and MGS). First, we definethe averaged observable X i = 1 N R N R (cid:88) n =1 X ( n ) i (40)where X ∈ { P gg , P gm } , i is the index running throughthe k binning, and the superscript n refers to the n -threalization of the considered observable. From this defini-tion we estimate the auto- and cross- covariance matricesas C X × Y,ij = 1 N R N R (cid:88) n =1 (cid:16) X ( n ) i − X i (cid:17) (cid:16) Y ( n ) j − Y j (cid:17) , (41)where once again X, Y ∈ { P gg , P gm } .Similarly to what was done in [31] we retain only thediagonal entries of the previously defined covariance ma-trices. This approximation is justified given the signif-icant low number density of the synthetic samples weconsider (actual values are listed in Table I). These num-bers translate in a significant shot-noise contribution tothe power spectrum, whose main effect is to boost thediagonal entries of the covariance matrix, leading to sub-dominant off-diagonal terms. Therefore we approximateour covariance matrices with a block-diagonal shape, andconsider only the auto- and cross-correlation at the samewavelength (i.e. C ij = 0 for k i (cid:54) = k j ). Additionally, inorder to reduce noise due to the limited number of in-dependent realizations, we compare the raw covariancematrix to Gaussian predictions [53] for each k bin, andretain the maximum of the two values.The error budgets of our galaxy samples are subse-quently rescaled in order to match the same effective vol-ume [95, 96], defined as V eff ( k (cid:63) ) = (cid:20) n P gg ( k (cid:63) )1 + n P gg ( k (cid:63) ) (cid:21) V , (42)where n and V are the mean number density and volumeof the considered sample, respectively, and we choosethe reference scale to evaluate the effective volume as k (cid:63) = 0 . h Mpc − . In this way the constraining powerand the signal-to-noise ratio of the three galaxy cata-logues are artificially set to roughly match the same am-plitude. We choose to rescale the covariance of both TABLE I. Main properties of the three synthetic galaxy samples used in this work. The table shows the label we use to identifyeach sample, the N-body run on which it is based, the total number of independent realizations, redshift, galaxy number density,effective volume scaling ratio (see Section III B), linear galaxy bias and large-scale deviation from Poisson shot-noise (in unitsof 1 /n ). In order to obtain fiducial values for the last two columns, we adopt the same strategy described in [31] that makesuse of the large-scale limit of the quantities P gg , P gm and P mm .Identifier Simulation N R z n (cid:2) ( h/ Mpc) (cid:3) η b N MGS
LasDamas
Carmen 40 0 .
132 1 . × − .
04 1 . ± . − . ± . LasDamas
Oriana 40 0 .
342 9 . × − . ± . − . ± . Minerva
100 0 .
57 4 . × − .
47 2 . ± . − . ± . h − Mpc.Simulation Ω m Ω Λ Ω b h n s σ LasDamas
Minerva
CMASS and MGS to match the effective volume of theLOWZ sample (cid:0) ≈ /h ) (cid:1) . Table I shows the rescal-ing factor for the three samples, that is simply defined asthe sample-to-LOWZ ratio between the respective effec-tive volumes. C. Fitting procedure and prior choices
We make use of the large number of independent real-izations for each synthetic galaxy sample to define theoverall likelihood function as the product of the indi-vidual N R likelihoods [31, 97]. In this way, under theassumption that the individual likelihoods are well de-scribed by a multivariate Gaussian distribution, we candefine the overall likelihood function as − L tot = − N R N R (cid:88) n =1 log L ( n ) = 1 N R N R (cid:88) n =1 χ n ) = 1 N R N R (cid:88) n =1 N b (cid:88) i,j =1 (cid:16) X ( n ) i − µ i (cid:17) C − X,ij (cid:16) X ( n ) j − µ j (cid:17) , (43)where X ( n ) i is the measurement from the n -th realizationof the i -th k bin, µ i is the corresponding model predic-tion, and C ij is the rescaled covariance matrix as de-scribed in Sec. III B. With this definition, our likelihoodends up coinciding with the likelihood of the mean of thedata with covariance C ij . In practice, there is a constant factor between the two definitions, that depends on thenumber of independent realizations N R and the numberof data points N b . This factor is taken into account whenderiving the goodness of fit for each of the tested config-urations.The inference of model parameters is carried outthrough a least- χ analysis based on a standardMetropolis-Hastings MCMC algorithm. The likelihoodincorporates a complete recipe for one-loop galaxy clus-tering and an interface to CAMB . For each tested con-figuration of the MCMC, we first run some preliminarychains that are needed to obtain a robust estimate ofthe parameter covariance, which is essential for the goodconvergence of the Markov chain in a highly dimensionalparameter space. We iterate this process twice, each timespecifying the parameter covariance obtained at the pre-vious step, before running the final set of chains. Theseare terminated as soon as the chains satisfy a standardGelman-Rubin convergence criterium, i.e. as soon as thebetween-chain and within-chain variances (we run 8 inde-pendent chains for each case, changing the initial randomseed) are in agreement within
R < . getdist [99] to extract the parameter posteriors.The models described in Section II embed a large num-ber of free parameters. In this analysis we vary bothcosmological and nuisance terms (the latter include bias,shot-noise and matter counter-term). Since we are con-straining two-point statistics only via galaxy clustering,we fix both the baryon density parameter Ω b h and theprimordial spectral index n s to their fiducial values (listedin Table II), and vary the CDM density parameter Ω c h and the Hubble constant h . We add the primordial scalaramplitude A s , but only when fitting the combination of P gg and P gm , since the fit of the auto galaxy power spec-trum alone cannot break the strong degeneracy between A s and b . On the contrary, while both P gg and P gm have the same dependency on A s , they scale differentlywith the linear bias, i.e. P gg ∝ b and P gm ∝ b in thelinear limit. We choose a flat prior for the three cosmo-logical parameters that is large enough to fully containtheir posterior distribution up to 2 σ even for the leastconstraining configuration.0 TABLE III. Adopted priors for the full list of cosmological andnuisance parameters of the fits. We impose a uniform prior(U) for all the sampling parameters. The scalar amplitude A s is not sampled over when fitting only P gg but it is whenadding the additional constraint from P gm . In all cases thesecond-order tidal bias γ is fixed to the excursion set relationdefined in Equation (39).COSMOLOGY h U[0 . , c h U[0 . , . A s P gg : P gg + P gm :Fixed to fiducial U[0 . , . b σ U[0 . , b σ U[ − , γ Fixed to γ ,ex ( b ) γ σ U[ − , β P [ k − ] U[ − , β × P [ k − ] U[ − , N [ n − ] U[ − , . N [ k − n − ] U[ − , N × [ k − n − ] U[ − , c U[0 , The same choice of uniform priors is chosen for the nui-sance parameters. We sample the combination of the biasparameters with the corresponding power of σ in orderto avoid strong degeneracies between the one-loop biasexpansion and the intrinsic non-linearity of the matterpower spectrum P mm , that at first order is well capturedby σ . As anticipated in Section II D, we fix γ to Equa-tion (39), and vary all other bias parameters freely. Wequote higher-derivatives parameters with respect to anarbitrary fiducial scale k HD = 0 . h Mpc − , while all theshot-noise parameters can be naturally expressed in unitsof 1 /n . Since terms involving N and N × both carry a k -dependency, we also express the latter in units of k HD .The full list of priors for the model parameters is shownin Table III. D. Rescaling of the input linear power spectrum
As discussed in Section III,
Minerva and
LasDamas adopt different Boltzmann solvers to obtain the powerspectra for the initial particle displacements. In partic-ular, the
Minerva runs assume initial conditions basedon CAMB, while the
LasDamas simulations make use of
FIG. 1. Ratio between the CMBFast and CAMB linear pre-diction at the reference cosmology of the
LasDamas simula-tions. The two shaded grey bands mark the intrinsic 1-sigmaerror of the LOWZ and MGS galaxy power spectrum P gg . the power spectrum computed with CMBFast . Figure 1shows the fractional deviation between the output of thetwo codes at the reference cosmology of the
LasDamas runs and at the redshift of the MGS sample. This dif-ference is compared to LOWZ and MGS intrinsic signal-to-noise ratio, represented by the two shaded areas. Thetwo spectra deviate at the level of ∼ .
5% in the range ofscales that are relevant to this analysis, and this mightrepresent an issue when modeling the galaxy-galaxy andgalaxy-matter power spectra. For this reason, when fit-ting LOWZ and MGS, we employ a rescaling schemebased on the renormalization of the input power spec-trum of
CAMB by the ratio shown in Figure 1. Thisapproximation is valid as long as the sampled cosmologyis not far from the fiducial one, and this becomes a com-pletely legitimate assumption when hitting exactly thetrue cosmology. However, we claim that even for cos-mologies that are slightly off with respect to the fiducialone, this rescaling scheme is already better than directlyusing the
CAMB prediction with no further correction.In order to perform this rescaling, we first apply theunits transformation h Mpc − → Mpc − . This is meantto assure that the power spectrum ratio is obtained atthe same set of scales, avoiding the dependency on theHubble parameter h [100]. At each step of the Markovchain, the current CAMB prediction is firstly trans-formed into Mpc − units, rescaled using the approach de-scribed above, and finally transformed back into h Mpc − units.1 E. Performance metrics
The validation of one-loop galaxy bias as described inSection II is mostly based on its range of validity. In orderto select the maximum mode at which our model stopsproviding good performances, we run multiple MCMCchains varying k max in the interval (0 . , . h Mpc − with a step of 0 . h Mpc − , leading to 9 different fit-ting ranges for each configuration. We stop at k max =0 . h Mpc − as pushing the model to even smaller scaleswould inevitably require accounting for two-loop contri-butions in the description of galaxy bias and non-linearmatter evolution. Clearly, exploiting the small-scale in-formation contained in the non-linear regime naturallyprovides a better constraining power for the model pa-rameters, but at the cost of introducing theory system-atics in their posterior distribution. As described in [31](see [62] for a similar approach), we make use of a setof three different performance metrics, whose goal is toprovide a quantitative way to compare the posteriors ex-tracted from our Markov chains, and to determine when aparticular model no longer gives an accurate descriptionof our datasets.
1. Figure of bias
The first quantity we are interested in evaluating is theability of the model to provide unbiased measurementsof its free parameters. In particular, in this work wefocus mostly on the systematic error of the cosmologicalparameters: the Hubble constant h , the cold dark matterdensity Ω c h , and, when considering the joint fit between P gg and P gm , the primordial scalar amplitude A s .We define the Figure of Bias (FoB) of the consideredmodel for a given parameter set θ asFoB( θ ) ≡ (cid:104)(cid:0) θ − θ fid (cid:1) (cid:124) S − (cid:0) θ − θ fid (cid:1)(cid:105) / , (44)where θ represents the mean of the posterior distribution, θ fid represents the fiducial position, and S is the param-eter covariance expressed in matrix form. In this way,we are simply quantifying the relative separation of themeasured parameter from its true value in terms of thevariance of the posterior distribution. If the FoB is evalu-ated only for one parameter, then the standard 68% − −
2, respec-tively. On the contrary, when more than one parameteris considered to compute the FoB, the two percentilescan be calculated by direct integration of a multivariateGaussian distribution with the corresponding dimension-ality. For the case we consider n = 2 ( n = 3) it followsthat the 68% −
95% percentiles correspond to a FoB of1 . − .
49 (1 . − .
2. Goodness-of-fit
Along with finding unbiased constraints on the pa-rameters of interest, we also ask our model to providea good description of the observables we use in the fit.We quantify the goodness-of-fit in terms of the standard χ extracted from the Markov chains (Equation (43)),but after rescaling it to account for the additional factor η introduced in the covariance of the data (see SectionIII B). With this rescaling, we can compare the recovered χ to the 68% and 95% percentiles of a χ distributionwith N dof degrees of freedom, where N dof = N R × N b − N p . (45)Here N R is the number of independent realizations foreach galaxy sample, N b is the total number of k binsincluded in the fit, and N p is the number of free modelparameters.
3. Figure of merit
The final metric we employ is based on the merit of theconsidered model with respect to the cosmological pa-rameters that are varied in the fit. We define the
Figureof Merit (FoM) of a given model for a subset of parame-ters θ as the inverse of the determinant of the parametercovariance matrix corresponding to the subset θ , i.e.FoM( θ ) ≡ (cid:113) det (cid:2) S ( θ ) (cid:3) . (46)Here S ( θ ) ≡ S ( θ ) / θ is the block corresponding tothe parameter subset θ rescaled by their fiducial values.With this last normalization we are explicitly asking fora relative FoM rather than an absolute quantity, thatwould be more difficult to compare to the one of otherparameters.
IV. TESTING ONE-LOOP GALAXY BIAS
In this section we analyse parameter posteriors ex-tracted from several Markov chains, where we vary boththe maximum fitting scale k max and the final expressionwe adopt to model one-loop galaxy power spectra. Thecomplete list of model configurations is summarised inTable IV, along with the total number of free parame-ters in each case. Based on results from [31], we decideto employ RESPRESSO to model P mm in this section, andto leave the comparison between different matter modelsto Section IV D. In this way we can focus exclusively onassessing the validity of one-loop bias in real-space leav-ing aside both the impact of non-linear matter evolutionand redshift-space distortions (that will be the topic of afuture paper).2 TABLE IV. List of all the considered combinations betweenthe modeling of P mm and the extensions to one-loop galaxybias, as described in Section II. For each model, we denotewith a check mark all the free parameters of the P gg fit, whilewe use a double check mark for the free parameters that areexclusive to the joint P gg + P gm fit. P mm RESPRESSO gRPT EFTBias model STD k2N HD STD k2N HD STD k2N h (cid:88) (cid:88) (cid:88) (cid:88) (cid:88) (cid:88) (cid:88) (cid:88) Ω c h (cid:88) (cid:88) (cid:88) (cid:88) (cid:88) (cid:88) (cid:88) (cid:88) A s (cid:88)(cid:88) (cid:88)(cid:88) (cid:88)(cid:88) (cid:88)(cid:88) (cid:88)(cid:88) (cid:88)(cid:88) (cid:88)(cid:88) (cid:88)(cid:88) b σ (cid:88) (cid:88) (cid:88) (cid:88) (cid:88) (cid:88) (cid:88) (cid:88) b σ (cid:88) (cid:88) (cid:88) (cid:88) (cid:88) (cid:88) (cid:88) (cid:88) γ σ (cid:88) (cid:88) (cid:88) (cid:88) (cid:88) (cid:88) (cid:88) (cid:88) β P (cid:88) (cid:88) β × P (cid:88)(cid:88) (cid:88)(cid:88) N (cid:88) (cid:88) (cid:88) (cid:88) (cid:88) (cid:88) (cid:88) (cid:88) N (cid:88) (cid:88) (cid:88) N × (cid:88)(cid:88) (cid:88)(cid:88) (cid:88)(cid:88) c (cid:88) (cid:88) For each configuration, we compute the three perfor-mance metrics defined in Section III E and show theirtrends as a function of k max . At the same time, weperform sanity checks on the linear bias b σ and thelarge-scale shot-noise parameter N , for which we havefiducial values to compare with (see Table I). As antici-pated in Section II D, throughout the rest of the paper wefix the second-order non-local bias parameter γ to theexcursion-set relation (Equation (39)), in order to breakthe strong degeneracy between the former and γ . Thischoice is also based on results from [31], where it wasshown how this choice led to stable and overall accuratemeasurement of the linear bias.In the next two sections we analyse the posterior de-rived from fits of the auto galaxy power spectrum P gg alone, and from the combination of the auto and crossgalaxy power spectra, P gg + P gm . In the latter case, theadditional information contained in P gm allows to breakthe degeneracy between the linear bias and the amplitudeof primordial fluctuations, allowing us to also sample over A s .In both cases, we adopt a criterion for defining therange of validity of a given model, which is based on acombination of FoB and goodness-of-fit. Following [31],we define the model-breaking statistics as σ MB ( k ) ≡ FoB( k )FoB + χ ( k ) − χ ( k ) − , (47)where the subscript percentages correspond to the per-centiles of the corresponding distribution (FoB and χ ). We say that a given model breaks down at a scale k † when the model-breaking statistics assumes a criticalvalue σ crit at k † , i.e. σ MB ( k † ) = σ crit . (48)We arbitrary choose the critical threshold σ crit = 1 .
5. Inthis way we would accept models with a maximum FoBone and a half times larger than the value correspondingto the 1-sigma of the FoB distribution, but only under aperfect recovery of the shape of the input dataset. Prac-tically, this case is unrealistic at high k max , as the χ pro-gressively deviates from the number of degrees of freedombecause of the weakening of the model. For this reason σ MB receives individual contributions from the FoB andthe goodness-of-fit. A. Validity of one-loop galaxy bias for the autopower spectrum
In this section we analyse results from the fits tothe auto power spectrum P gg . We fix the descriptionof the non-linear matter power spectrum and test dif-ferent extensions of one-loop galaxy bias. The stan-dard (STD) model is based on two free cosmological pa-rameters ( h, Ω c h ) plus four free nuisance parameters( b σ , b σ , γ σ , N ). The k -dependent shot-noise(k2N) and higher-derivatives (HD) extensions have oneadditional free parameter each, i.e. N and β P , respec-tively.In Figure 2 we show the three performance metricsextracted from the three different P gg measurements(CMASS, LOWZ and MGS), plotted against the max-imum scale k max included in the fit. In this case weshow the FoB and FoM corresponding to the combina-tion h, Ω c h ), i.e. to all the cosmological parametersthat are varied in the model. The two grey shaded areasin the panels corresponding to FoB and goodness-of-fitmark the 68th and 95th percentiles of the correspondingquantity. As described in Section III E 1, the 68% − n = 2 number ofdimensions.For the three galaxy samples, we find that the STDmodel (blue line) alone is enough to provide a good de-scription of the dataset in terms of goodness-of-fit, as thereduced χ constantly stays within the 68% (dark greyshaded area) of the corresponding χ distribution. Atthe same time, the combined constraint on the ( h, Ω c h )pair is unbiased, showing a multivariate posterior distri-bution that is able to capture the fiducial position at typ-ically ∼ . σ . Adding either an additional k -dependentshot-noise (red line) or the leading higher derivatives cor-rection (green line) does not significantly alter the twometrics, while reducing the amplitude of the correspond-ing FoM ( ∼
10% when adding the stochastic contribution3
FIG. 2. FoM, FoB and goodness-of-fit extracted from chains corresponding to the fits of the galaxy auto power spectrum P gg .Performance metrics refer to the combinations of h and Ω c h , which are the only two free cosmological parameters of the model.The standard (STD), k -dependent shot-noise (k2N), and higher derivatives (HD) models are colour-coded as shown in thelegend. In all three cases, the non-linear matter power spectrum is modelled with RESPRESSO , and the second-order tidal bias γ is fixed to Equation (39). The two shaded grey regions mark the 68th and 95th percentiles of the FoB and χ distribution. N , ∼ −
30% when adding the higher derivatives pa-rameter β P ). Even for the STD model, for which weexpect the model to fail earlier than the other two ex-tensions) the model-breaking scale k † is higher than themaximum scale we probe ( k max = 0 . h Mpc − ). How-ever, we decide to stop at this scale without exploringsmaller scales, even if the model-breaking criterion weremost likely satisfied at higher k max . As a matter of factit would be hard to separate the effective goodness ofthe model from the impact of two-loop contributions togalaxy bias.We notice how, similarly to what was obtained in [31]when focusing on the linear galaxy bias, there is a ten-dency of the FoM to flatten above an approximate scale of k max = 0 . h Mpc − . This is an indication that pushingthe non-linear model to increasingly higher wave modesdoes not necessarily imply a more stringent measurementof the cosmological parameters, as most of the additionalconstraining power is absorbed by nuisance parameterssuch as higher-order galaxy bias contributions.In Figure 3 we show the dependence on k max of the lin-ear bias parameter b σ and the constant stochastic con- tribution N , for which we have fiducial values obtainedexploiting the large-scale limit of the measured galaxyand matter power spectra. Although the extension toeither the k2N or the HD model enlarges the size of theerrorbars, we notice how the former model is the only onecapable of providing a simultaneous unbiased measure-ments of the two parameters, for the three galaxy sam-ples we are considering. In particular the LOWZ sam-ple shows a clear detection of N , that, if not accountedfor, can lead to a significant ( > σ ) systematic effect on N above k max ∼ . h Mpc − , while partially recover-ing the true amplitude of P gg with an underestimationof b σ . The only exception is represented by the MGSsample, whose linear bias is constantly overestimated bya ∼
2% factor for all the three modeling assumptions,and that is consistent with the fiducial value only up to2- σ .4 FIG. 3. Marginalised constraints on the linear bias b σ (top) and the constant stochastic term N (bottom) as a functionof k max , for the fits of the galaxy auto power spectrum P gg . In all cases dark matter non-linear evolution is modelled using RESPRESSO . Fiducial values of the linear bias and constant shot-noise contribution are reproduced with black dashed lines.Uncertainties on the fiducial values are marked with shaded grey bands (marking 1-and 2- σ confidence intervals). Fiducialmeasurements are listed in Table I. B. Consistency between auto and cross powerspectra
In this section we focus on a much more stringent test,corresponding to the simultaneous fit of both the galaxyauto power spectrum P gg and the galaxy-matter crosspower spectrum P gm . As a matter of fact, a good accu-racy for this combined statistics might become crucial inanalyses that aim to exploit the whole information con-tained in galaxy clustering and galaxy-galaxy weak lens-ing, as most of the upcoming large observational projectsare going to do (3 × P gg and P gm on the linear bias b leads to the breaking of the strong b − A s degeneracy.Therefore, we additionally sample the scalar amplitude A s , effectively extending by one the number of degrees offreedom of our model. Also in this case, we concentrateon the one-loop biasing scheme, fixing the description ofthe matter clustering to the output of RESPRESSO . Wetest the standard (STD) model against the k -dependentshot-noise (k2N) and higher derivatives (HD) templates,for which we introduce two more degrees of freedom (onefree parameter for P gg and P gm , each). In summary, theSTD model is based on three free cosmological param-eters ( h , Ω c h , A s ) plus four free nuisance parameters( b σ , b σ , γ σ , N ). The k2N and HD models havetwo additional free parameters each, i.e. ( N , N × ) and( β P , β × P ), respectively. Note that that the second-ordertidal bias γ is not a free parameter, but it is fixed to the excursion-set relation defined in Equation (39) for all thecases we consider.In Figure 4 we show the three performance metrics ina similar way as we did for the P gg -only fits. In this caseboth FoM and FoB are computed from the combinationsof the three cosmological parameters ( h, Ω c h , A s ). Wefirst notice how the use of the STD model is no longersufficient to provide an accurate recovery of the cosmo-logical parameters on the overall range of scales we areconsidering. Indeed, the FoB of both CMASS and LOWZquickly increases above 2- σ at an approximate scale of k max ∼ . h Mpc − and keeps getting worse at higherwave modes. The same trend is observed in panels re-ferring to the goodness-of-fit, as the reduced χ starts tofall outside the 95% confidence region at approximatelythe same scale cut. According to the model-breaking cri-terion, the range of validity of the STD model is limitedto scales below k ∼ . h Mpc − (this is reflected in theplot by the solid-to-dashed line transition).Adding k -dependent stochasticity increases both theaccuracy in the recovery of the cosmological parametersand the overall likelihood between data vectors and best-fitting model. Indeed, the k2N model shows a FoB whichis constantly within 1- σ for the three galaxy samples.Moreover, the goodness-of-fit is significantly improvedwith respect to the STD model, with the χ measured atthe highest k max still being consistent within 1- σ . Thecombination of FoB and goodness-of-fit can be observedvia the model-breaking criterion on the FoM, for which5 FIG. 4. Same as Figure 2 but for the combinations of the galaxy auto and cross power spectrum, P gg and P gm . In thiscase, performance metrics (FoM and FoB) refer to the combinations of the three cosmological parameters ( h, Ω c h , A s ). Thestandard (STD), k -dependent shot-noise (k2N), and higher derivatives (HD) models are colour-coded as shown in the legend.The solid-to-dashed transition mark the model-breaking scale k † . In all three cases, the non-linear matter power spectrum ismodelled with RESPRESSO . The two shaded grey regions mark the 68th and 95th percentiles of the FoB and χ distribution. the k2N model is the only one managing to be acceptedup to k max = 0 . h Mpc − .The HD model also provides a better recovery (com-pared to STD) of the cosmological parameters, within1- σ from the fiducial position, but it fails earlier thanthe k2N model to provide the required goodness-of-fitfor both CMASS and LOWZ, hinting for a stronger ne-cessity of adding k -dependent stochasticity with respectto short-range non-localities.In terms of FoM amplitude, the k2N and HD mod-els behave in a mostly similar way, significantly reducingthe statistical constraints of the cosmological parametersat any given k max cut. However, we notice how accord-ing to the model-breaking criterion both models allowto push the analysis to higher k max with respect to theSTD case. In other words, a proper comparison of themerit of the three models should be done considering theindividual scale k † after which the model breaks down.Considering the CMASS sample, the STD model pro-vides a FoM consistent with the one obtained from theHD case, and ∼
30% larger than the one of the k2Nmodel. However, the latter two models have yet to reach the maximum scale k † corresponding to the breaking ofthe model. A simple exercise of extending the k2N modelto slightly larger modes values shows that the latter turnsout to be competitive against the STD model at a scale of k max ∼ . h Mpc − . As for LOWZ, the HD case is sig-nificantly less constraining than the STD model ( ∼ k max = 0 . h Mpc − but stops earlier. In con-trast, the k2N model outperforms the standard model bya ∼
10% factor, although also in this case the modelbreaking scale k † has yet to be reached.A particular case is represented by MGS, that differ-ently from the previous two samples shows a stable andaccurate recovery of cosmological parameters all the wayup to k max = 0 . h Mpc − even with the STD model.The reason is most likely related to our choice of em-ploying RESPRESSO to model the non-linear matter evo-lution. This ultimately provides an accurate descriptionof P mm even at the low redshift of MGS, leaving only thegalaxy-matter relation to be modelled. Since MGS fea-tures the lowest large-scale bias among the samples weconsider, it is likely that the one-loop standard expan-6 FIG. 5. Same as Figure 3 but for the combined fits of the galaxy auto power spectrum P gg and the galaxy-matter cross powerspectrum P gm . In all cases dark matter non-linear evolution is modelled using RESPRESSO . Fiducial values of the linear bias andconstant shot-noise contribution are reproduced with black dashed lines. Uncertainties on the fiducial values are marked withshaded grey bands (marking 1- and 2- σ confidence interval). sion is enough to provide both a good χ and an accu-rate parameter recovery on the full range of scales up to k max = 0 . h Mpc − .Differently from the fit of the galaxy auto-power spec-trum, in this case the FoM for the three different galaxysamples do not show a flattening on the range of scaleswe are considering. Checking the independent contribu-tions coming from the individual cosmological parame-ters, we notice how the FoM does actually flatten forboth h and A s , whereas it monotonically increases forΩ c h . This is a hint showing that the matter densityparameter still benefits of additional informations con-tained in the mildly non-linear regime, when combiningthe two galaxy power spectra.In Figure 5 we show the dependence of the linear bias b σ and the constant shot-noise contribution on k max ,similarly to what we showed in Figure 3 for P gg alone.The more stringent requirement of simultaneously fitting P gg and P gm makes the standard model fail in providingaccurate measurements of both parameters. In particu-lar, b σ starts to be biased at approximately the samescale for which the model also provides incorrect mea-surements of the cosmological parameters. Adding either k -dependent shot noise terms or higher derivatives pa-rameters help in recovering the true values of the bias,with a slight preference for the former. More generallythe k2N model is the only one that is capable of simul-taneously providing unbiased measurements of the cos-mological parameters, the linear bias and the large-scalestochastic contribution to the galaxy power spectrum. This is also in agreement with the results obtained in[31] at fixed cosmology. C. Constraints on stochasticity andhigher-derivative parameters
After having focused on the marginalised posteriorsof the cosmological parameters, linear bias and constantshot noise, here we try to put internal constraints on boththe k -dependent stochastic parameters and the higherderivatives parameters. In particular, we want to checkwhether these terms can be detected with statistical sig-nificance under the assumption of using a reference effec-tive volume like the one described in Section III.Figure 6 shows marginalised constraints for the twosets of parameters (extracted from the P gg + P gm chains)as a function of k max . For consistency check, in the pan-els referring to parameters that enter the modeling of thegalaxy auto power spectrum, N and β P , we also showless stringent constraints extracted from the chains cor-responding to the fit of P gg (shaded grey band).The consistency of the k -dependent shot-noise modelis additionally reinforced by the trend observed in thetop panels of Figure 6. As a matter of fact, the re-lation between N /N × and k max can be well describedby a straight line with zero slope, this fact intrinsicallyhinting at a correct behaviour for both parameters, astheir profiles remain stable when adding constraints athigher wave modes. Moreover, there is a clear detec-7 FIG. 6. Marginalised posterior of the k -dependent stochastic parameters N and N × (top), and the higher derivativesparameters β P and β × P (bottom), as a function of k max , for the combination of P gg + P gm . In all cases dark matter non-linearevolution is modelled using RESPRESSO . For consistency check we also show the 1- σ standard deviations of the parameters N and β P obtained from the fit of P gg only as a shaded grey area. tion of the k -dependent shot-noise parameters, with astatistical significance that further reinforce the impactof adding N and N × in terms of goodness-of-fit. In-deed N and N × are detected at k max = 0 . h Mpc − with a statistical significance of 3 . σ and 2 . σ , respec-tively, for the CMASS sample. With LOWZ the detec-tion becomes even larger, with numbers that are closeto 6- σ and 4- σ , respectively. The MGS sample is theone showing the least significant detection of both pa-rameters, with a possible modification of the recoveredvalue of N × at k max = 0 . h Mpc − . However, we re-mind that the STD model performs surprisingly well onthis dataset, and therefore the weaker detection of thetwo parameters is partially expected.The higher derivatives parameters also show stable re-sults as a function of k max , with the highest detectionrepresented by β × P for LOWZ. However on these scalesthe HD model is already broken, as shown in Figure 4.For the other cases, the typical significance of the detec-tion is set to ∼ σ and ∼ σ for β P and β × P , respectively.Differently from [31], we do not observe incompatible results between the marginalised posteriors of N and β P from fits of P gg and P gg + P gm . However in this case wesignificantly extended the dimensionality of the param-eter space leading to weaker constraints on the modelparameters, so that the deviation may be hidden by thelarger statistical noise. D. Results for alternative matter models
So far we have fixed the modeling of the non-linearmatter power spectrum using a hybrid approach such as
RESPRESSO . Another way of proceeding is to employ aperturbative description not only for one-loop bias butalso for P mm . In this section, we analyse the impact of adifferent modeling of P mm inside Equation (10) and (11),and for this purpose we make use of the two additionalmodels described in Section II B 4.We remind the reader that these models correspond to1) an EFT-like model based on resummation of IR modesand non-trivial stress tensor that at leading order in real-8 FIG. 7. FoM, FoB and goodness-of-fit extracted from fits of P gg + P gm . Performance metrics are separated into each individualcomponent, i.e. h , Ω c h and A s . Different dark matter models and different one-loop modeling assumptions are color- andstyle-coded as shown in the legend. The two shaded grey regions mark the 68th and 95th percentiles of the FoB and χ distribution. c , that is completely degenerate with lead-ing order higher derivatives terms. For this reason, in thiscase the STD model already includes a higher-derivativecorrection. The full list of configurations, together withfree model parameters, is given in Table IV.Figure 7 shows the performance metrics for the threedifferent dark matter models when fitting the combina-tion of P gg and P gm . In this case, we show the individualcontributions to FoM and FoB from the three cosmologi-cal parameters, h , Ω c h and A s . In addition, for the sakeof simplicity, we identify the model breaking scale k † withthe end of the corresponding line in the FoM panels.We first notice how the recovered goodness-of-fit isseverely worsened when adopting gRPT to describe thenon-linear matter power spectrum, particularly whenconsidering the STD model. We separately tested theimpact of extending this model to include higher deriva-tives, and while this extension slightly improves the rangeof validity of the gRPT-based model, it still fails in be-ing accepted up to k max = 0 . h Mpc − . The corre-sponding k2N model definitely improves the recoveredtrends but the net result is a worse performance respectto RESPRESSO and EFT. The STD model of EFT breaksat higher values of k max , and this is somewhat expectedgiven the presence of the additional free parameter c .Nevertheless in order to provide a χ consistent within95% confidence interval at all scales we still have to in-clude additional stochastic terms.The FoB is overall consistent among the three differentparameters, showing a strong bias on the marginalisedposterior ( > σ at k max > . h Mpc − ) when adopt-ing the STD for either gRPT or RESPRESSO . Similarlyto the goodness-of-fit, the STD case combined with EFTprovides better results (consistent within 1- σ ) since it al-ready incorporates a k P L ( k ) correction which is absentfrom the other matter models. The overall interpretationis that all three matter models strongly hint for the ne-cessity of adding k -dependent stochastic contributions,both to produce an accurate description of the joint datavector P gg + P gm and to provide unbiased measurementsof the cosmological parameters.As for the FoM, we notice how its dependency on k max is significantly different between the ( h, A s ) pairand Ω c h . For the former, the difference in performancebetween STD and extended models is tiny ( ∼ k -dependent stochasticity(or higher derivatives) results in the suppression of theFoM at a reference scale of k max = 0 . h Mpc − by a factor of ∼
50% and ∼ k max ∼ . h Mpc − , when using the STD model (slightlylarger with the HD extension). Substituting gRPT with RESPRESSO helps in extending the range of validity ofthe HD model, as shown in Section IV B, while leavingalmost untouched the model breaking scale of the STDcase. However, this last configuration, i.e.
RESPRESSO with the STD biasing scheme systematically results inone of the best performing models in terms of maxi-mally achievable FoM. Adding k -dependent stochasticterms to both matter models extends the range of va-lidity up to k max = 0 . h Mpc − , with a value of FoMwhich is typically consistent with the one of the STDcase at the corresponding k † . Finally, the range of va-lidity of the EFT model is hitting our maximum valueof 0 . h Mpc − already with the STD model for bothCMASS and MGS, while being limited to intermediatevalues ( k max ∼ . h Mpc − ) for LOWZ. The latter canbe also be adjusted by considering k -dependent shot-noise parameters. E. Impact of P gm on the Figure of Merit After having assessed the impact of different modelingassumptions we can now check the statistical significanceof adding the galaxy-matter cross power spectrum to afull shape measurement fit like the one we carried out inthis analysis. Indeed, as we already pointed out in Sec-tion IV B 3 × P gm .In Figure 8 we compare the marginalised posteriordistribution of the cosmological parameters and linearbias obtained when fitting only P gg and the combina-tion P gg + P gm . For this comparison we model P mm with RESPRESSO , include k -dependent stochastic contri-butions, and include all scales up to k max = 0 . h Mpc − .We remind that the comparison is not properly fair, as wefix the scalar amplitude A s in the fits of the auto powerspectrum alone, while this is treated as a free parameterin the fits of the combined statistics.The observed trend when including P gm in the fit is asignificant reduction of the statistical uncertainty of allparameters. In particular we measure a 1- σ standard de-viation for h and Ω c h that is 1 .
4, 2, and 1 . b σ by a factor 3 for all three samples.0 FIG. 8. Comparison between the posterior parameter distribution (68% and 95% confidence intervals) obtained by fitting thegalaxy auto power spectrum P gg (blue) and the combination of galaxy-galaxy and galaxy matter power spectrum P gg + P gm (red). We show constraints obtained at k max = 0 . h Mpc − for the cosmological parameters h , Ω c h , A s , and the linear bias b σ . In all cases we assume non-linear dark-matter evolution to be described by RESPRESSO . Dashed solid lines correspond tothe fiducial values of the corresponding parameter. Grey solid bands mark the 1- σ and 2- σ error on the fiducial linear biasparameter. V. CONCLUSIONS
In this work we performed a full-shape analysis of thegalaxy power spectrum meant to assess the robustnessof one-loop galaxy bias models, extending the investiga- tion carried out in [31] to include also the sampling overcosmological parameters. Since we deal with galaxy clus-tering in real space, we decided to fix both the baryondensity Ω b h and the spectral index n s while leaving asfree parameters the cold dark matter density Ω c h , theHubble constant h , and the scalar amplitude A s . The1latter is kept as a free parameter only when consideringjoint fits to the galaxy-galaxy and galaxy-matter powerspectra.We measured both observables from a set of three dif-ferent synthetic galaxy samples, whose clustering prop-erties are meant to reproduce the ones of three real datacatalogues, i.e. CMASS and LOWZ from BOSS, andMGS from SDSS. We rescaled the statistical uncertain-ties of each of these samples to match an effective vol-ume of 6 (cid:0) h − Gpc (cid:1) (i.e. the one of the LOWZ sample),which is representative of the volume of tomographic binsfrom next-generation galaxy surveys (e.g. Euclid [104]).The analytical recipes we adopted to model these ob-servables are based on a standard one-loop expansion ofthe galaxy density field on the matter density field, col-lecting terms related to both spherically-symmetric grav-itational collapse and tidal fields up to third-order in per-turbations. In the first part of the paper, we fixed the de-scription of the non-linear matter power spectrum using ahybrid perturbative-simulated approach, i.e. RESPRESSO ,which was already shown in [31] to provide the closestbehaviour to the true measured matter power spectrumon the overall range of scales and redshifts that we con-sider. In a later section we explicitly tested the impactof different matter modeling assumptions, by employingan EFT-like model and a RPT-derived model (gRPT).In all the considered cases we fixed the quadratic tidalbias γ using an excursion-set-derived relation, in orderto break the strong degeneracy between the latter andthe cubic non-local parameter γ .Our main interest was in understanding the range ofvalidity of the standard one-loop expansion, and if oursynthetic galaxy sample would hint for the necessity ofadding additional terms to the standard recipe. For this,we considered two extensions of the standard model thattake into account either higher-derivatives of the gravi-tational potential or scale-dependent stochasticity, as ex-pected from short-range non-locality due to galaxy for-mation and the halo-exclusion effect, respectively.In order to quantify the goodness of each model wemake use of a set of three different performance metrics,which are 1) the Figure of Bias, to assess the level of biasintroduced in the recovery of cosmological parameters,2) the goodness-of-fit, to quantify the likelihood betweenthe input data vector and the best fit model, and 3) theFigure of Merit, to compare the constraining power ofeach model at any given k max . All of these quantitiesare easily obtained by post-processing the MCMC chainsthat we run for each combination of matter modeling,one-loop bias extensions, and maximum mode k max .Our results can be summarised as follows:(i) When using the fiducial description of the non-linear matter power spectrum (i.e RESPRESSO ), wefind that a standard one-loop bias model (linearbias, quadratic bias, cubic non-local bias, con-stant shot-noise parameter) with fixed quadratictidal bias can provide a good description of thegalaxy power spectrum for all our samples up to k max = 0 . h Mpc − while also returning unbiasedvalue of the cosmological parameter set (cid:0) h, Ω c h (cid:1) .(ii) Similarly to the fixed cosmology analysis of [31]we notice that when we consider the additionalinformation from the galaxy-matter cross powerspectrum, the standard model is no longer capa-ble of providing a good performance on the overallrange of scales we consider for all the three sam-ples. In particular, this model breaks at a scaleof k max = 0 . h Mpc − and k max = 0 . h Mpc − for CMASS and LOWZ, respectively, while beingaccepted all the way up to k max = 0 . h Mpc − for MGS. This might be representative of the levelof non-linearities in the galaxy-matter relationship,which at first order can be captured by the param-eter b σ (that is ∼ .
52 for LOWZ, ∼ .
26 forCMASS and ∼ .
06 for MGS).(iii) Extending the standard model to account also forthe presence of either higher derivatives or scale-dependent stochasticity provides a better perfor-mance both in terms of Figure of Bias (in this casefor the parameter set ( h, Ω c h , A s )) and goodness-of-fit. Adding scale-dependent stochastic terms re-stores the range of validity of the model up to k max = 0 . h Mpc − for both CMASS and LOWZ,while the importance of higher-derivatives is lim-ited by a worse χ for k max (cid:38) . h Mpc − , andslightly worse recovery of the linear bias parameterand the constant shot-noise correction. As an ad-ditional check, we verify that the marginalised pos-terior distribution of both k -dependent terms andhigher derivatives parameters is stable as a functionof k max . However we notice how the maximallyachievable FoM of the extended models is consis-tent with the one of the standard configuration atits lower k max , leading to equivalent statistical con-straints on the cosmological parameters.(iv) Changing the description of the non-linear mat-ter power spectrum to either EFT or gRPT in-duces modifications to the range of validity of one-loop galaxy bias models. Overall, we find that,as highlighted in [31], RESPRESSO is the best per-forming model, followed by EFT and gRPT. TheEFT-based model is penalised in terms of Figureof Merit because of the presence of an additionaldegree of freedom (two in case of the combinedfits), but it compensates this with an extendedrange of validity. On the contrary, gRPT breakssooner than the other models, at a typical scale of k max ∼ . h Mpc − , which is consistent with therange of validity observed in direct fits of the mat-ter power spectrum. This indicates that this modelwould also benefit of the presence of an additionalfree parameter representing a non-zero stress ten-sor such as in the case of the EFT-based model, asexpected. Although the gRPT-based model breaks2earlier than the other two models, its maximallyachievable FoM is consistent with the one recov-ered using RESPRESSO .(v) The additional constraints coming from the galaxy-matter cross power spectrum results in an improvedstatistical precision on measurements of cosmolog-ical parameters. Although we fix the scalar ampli-tude A s when fitting only the galaxy auto powerspectrum, we still find that in the combined casewe achieve better constraints on both h and Ω c h by a factor 1 .
4, 2 and 1 . h − Gpc . Al-though the results we obtained are partially based on thischoice, as shown in [31], this volume can be well represen-tative of individual tomographic redshift bins that will beadopted by next-generation galaxy surveys. In additionwe explore the additional constraining power brought byadding the galaxy-mass cross power spectrum to galaxypower spectrum fits, which is a relevant study in the con-text of 3 × ACKNOWLEDGMENTS
AP and MC acknowledge support from the Span-ish Ministry of Science and Innovation throughgrants PGC2018-102021-B-100 and ESP2017-89838-C3-1-R, and EU grants LACEGAL 734374 and EWC 776247with ERDF funds. IEEC is funded by the CERCA pro-gram of the Generalitat de Catalunya. AE acknowledgessupport from the European Research Council (grantnum- ber ERC-StG-716532-PUNCA. AGS acknowledgessupport from the Excellence Cluster ORIGINS, which isfunded by the Deutsche Forschungsgemeinschaft (DFG,German Research Foundation) under Germany’s Excel-lence Strategy - EXC-2094 - 390783311 AP and AGSwould like to thank Daniel Farrow, Jiamin Hou, MarthaLippich and Agne Semenaite for their help and usefuldiscussions. The analysis presented here has been per-formed on the high-performance computing resources ofthe Max Planck Computing and Data Facility (MPCDF)in Garching. This research made use of matplotlib, aPython library for publication quality graphics [105]. [1] M. Davis and P. J. E. Peebles, ApJ , 465 (1983).[2] S. J. Maddox, G. Efstathiou, W. J. Sutherland, andJ. Loveday, MNRAS , 43 (1990).[3] M. Tegmark, M. R. Blanton, M. A. Strauss, F. Hoyle,D. Schlegel, R. Scoccimarro, M. S. Vogeley, D. H. Wein-berg, I. Zehavi, A. Berlind, and et al., The Astrophys-ical Journal , 702–740 (2004).[4] S. Cole, W. J. Percival, J. A. Peacock, P. Norberg,C. M. Baugh, C. S. Frenk, I. Baldry, J. Bland-Hawthorn,T. Bridges, R. Cannon, and et al., Monthly Notices ofthe Royal Astronomical Society , 505–534 (2005).[5] M. G. Hauser and P. J. E. Peebles, ApJ , 757 (1973).[6] N. Kaiser, ApJ , L9 (1984).[7] M. Davis, G. Efstathiou, C. S. Frenk, and S. D. M.White, ApJ , 371 (1985).[8] J. M. Bardeen, J. R. Bond, N. Kaiser, and A. S. Szalay,ApJ , 15 (1986).[9] M. Tegmark, M. Blanton, M. Strauss, F. Hoyle,D. Schlegel, R. Scoccimarro, M. Vogeley, D. Wein-berg, I. Zehavi, A. Berlind, T. Budavari, A. Connolly,D. Eisenstein, D. Finkbeiner, J. Frieman, J. Gunn,A. Hamilton, L. Hui, B. Jain, D. Johnston, S. Kent,H. Lin, R. Nakajima, R. Nichol, J. Ostriker, A. Pope,R. Scranton, U. Seljak, R. Sheth, A. Stebbins, A. Sza-lay, I. Szapudi, L. Verde, Y. Xu, J. Annis, N. Bah-call, J. Brinkmann, S. Burles, F. Castander, I. Csabai,J. Loveday, M. Doi, M. Fukugita, J. Gott, G. Hennessy, D. Hogg, ˇZ. Ivezi´c, G. Knapp, D. Lamb, B. Lee, R. Lup-ton, T. McKay, P. Kunszt, J. Munn, L. O’Connell,J. Peoples, J. Pier, M. Richmond, C. Rockosi, D. Schnei-der, C. Stoughton, D. Tucker, D. Vanden Berk,B. Yanny, and D. York, ApJ , 702 (2004).[10] I. Zehavi, Z. Zheng, D. H. Weinberg, J. A. Frieman,A. A. Berlind, M. R. Blanton, R. Scoccimarro, R. K.Sheth, M. A. Strauss, I. Kayo, Y. Suto, M. Fukugita,O. Nakamura, N. A. Bahcall, J. Brinkmann, J. E. Gunn,G. S. Hennessy, ˇZ. Ivezi´c, G. R. Knapp, J. Loveday,A. Meiksin, D. J. Schlegel, D. P. Schneider, I. Szapudi,M. Tegmark, M. S. Vogeley, and D. G. York, ApJ ,1 (2005).[11] V. Desjacques, D. Jeong, and F. Schmidt, Phys. Rep. , 1 (2018).[12] H. J. Mo and S. D. M. White, MNRAS , 347 (1996),arXiv:astro-ph/9512127.[13] H. J. Mo, Y. P. Jing, and S. D. M. White, MNRAS , 189 (1997).[14] J. Frieman and E. Gaztanaga, ApJ , 392 (1994).[15] K. C. Chan, R. Scoccimarro, and R. K. Sheth,Phys. Rev. D , 083509 (2012), arXiv:1201.3614 [astro-ph.CO].[16] T. Baldauf, U. Seljak, V. Desjacques, and P. McDonald,Phys. Rev. D , 083540 (2012), arXiv:1201.4827 [astro-ph.CO]. [17] M. Manera and E. Gazta˜naga, MNRAS , 383 (2011),arXiv:0912.0446 [astro-ph.CO].[18] N. Roth and C. Porciani, MNRAS , 829 (2011),arXiv:1101.1520 [astro-ph.CO].[19] J. E. Pollack, R. E. Smith, and C. Porciani, MNRAS , 3469 (2012), arXiv:1109.3458.[20] J. E. Pollack, R. E. Smith, and C. Porciani, MNRAS , 555 (2014), arXiv:1309.0504 [astro-ph.CO].[21] P. McDonald and A. Roy, J. Cosmology Astropart.Phys. , 020 (2009), arXiv:0902.0991 [astro-ph.CO].[22] V. Desjacques, Phys. Rev. D , 103503 (2008),arXiv:0806.0007 [astro-ph].[23] V. Desjacques, M. Crocce, R. Scoccimarro, and R. K.Sheth, Phys. Rev. D , 103529 (2010), arXiv:1009.3449[astro-ph.CO].[24] R. J. Scherrer and D. H. Weinberg, ApJ , 607 (1998),arXiv:astro-ph/9712192 [astro-ph].[25] A. Dekel and O. Lahav, ApJ , 24 (1999),arXiv:astro-ph/9806193 [astro-ph].[26] A. Taruya and J. Soda, ApJ , 46 (1999), arXiv:astro-ph/9809204 [astro-ph].[27] T. Matsubara, ApJ , 543 (1999), arXiv:astro-ph/9906029 [astro-ph].[28] R. Sheth and G. Lemson, MNRAS , 767 (1999).[29] R. E. Smith, R. Scoccimarro, and R. K. Sheth,Phys. Rev. D , 063512 (2007), arXiv:astro-ph/0609547.[30] T. Baldauf, U. Seljak, R. E. Smith, N. Hamaus,and V. Desjacques, Phys. Rev. D , 083507 (2013),arXiv:1305.2917 [astro-ph.CO].[31] A. Eggemeier, R. Scoccimarro, M. Crocce, A. Pezzotta,and A. G. S´anchez, Phys. Rev. D , 103530 (2020).[32] A. G. S´anchez, R. Scoccimarro, M. Crocce, J. N. Grieb,S. Salazar-Albornoz, C. Dalla Vecchia, M. Lippich,F. Beutler, J. R. Brownstein, C.-H. Chuang, D. J. Eisen-stein, F.-S. Kitaura, M. D. Olmstead, W. J. Percival,F. Prada, S. Rodr´ıguez-Torres, A. J. Ross, L. Samushia,H.-J. Seo, J. Tinker, R. Tojeiro, M. Vargas-Maga˜na,Y. Wang, and G.-B. Zhao, MNRAS , 1640 (2017),arXiv:1607.03147.[33] H. Gil-Mar´ın, W. J. Percival, L. Verde, J. R. Brown-stein, C.-H. Chuang, F.-S. Kitaura, S. A. Rodr´ıguez-Torres, and M. D. Olmstead, MNRAS , 1757(2017), arXiv:1606.00439 [astro-ph.CO].[34] F. Beutler, H.-J. Seo, S. Saito, C.-H. Chuang, A. J.Cuesta, D. J. Eisenstein, H. Gil-Mar´ın, J. N. Grieb,N. Hand, F.-S. Kitaura, C. Modi, R. C. Nichol, M. D.Olmstead, W. J. Percival, F. Prada, A. G. S´anchez,S. Rodriguez-Torres, A. J. Ross, N. P. Ross, D. P.Schneider, J. Tinker, R. Tojeiro, and M. Vargas-Maga˜na, MNRAS , 2242 (2017), arXiv:1607.03150[astro-ph.CO].[35] J. N. Grieb, A. G. S´anchez, S. Salazar-Albornoz,R. Scoccimarro, M. Crocce, C. Dalla Vecchia, F. Monte-sano, H. Gil-Mar´ın, A. J. Ross, F. Beutler, S. Rodr´ıguez-Torres, C.-H. Chuang, F. Prada, F.-S. Kitaura, A. J.Cuesta, D. J. Eisenstein, W. J. Percival, M. Vargas-Maga˜na, J. L. Tinker, R. Tojeiro, J. R. Brown-stein, C. Maraston, R. C. Nichol, M. D. Olmstead,L. Samushia, H.-J. Seo, A. Streblyanska, and G.-b.Zhao, MNRAS , 2085 (2017), arXiv:1607.03143.[36] M. M. Ivanov, M. Simonovi´c, and M. Zaldarriaga, arXive-prints , arXiv:1909.05277 (2019), arXiv:1909.05277[astro-ph.CO]. [37] G. D’Amico, J. Gleyzes, N. Kokron, D. Markovic,L. Senatore, P. Zhang, F. Beutler, and H. Gil-Mar´ın, arXiv e-prints , arXiv:1909.05271 (2019),arXiv:1909.05271 [astro-ph.CO].[38] T. Tr¨oster, A. G. S´anchez, M. Asgari, C. Blake,M. Crocce, C. Heymans, H. Hildebrandt, B. Joachimi,S. Joudaki, A. Kannawadi, C.-A. Lin, and A. Wright,A&A , L10 (2020), arXiv:1909.11006 [astro-ph.CO].[39] A. Eggemeier, R. Scoccimarro, R. E. Smith, M. Crocce,A. Pezzotta, and A. G. S´anchez, “Testing one-loopgalaxy bias: joint analysis of power spectrum and bis-pectrum,” (2021), arXiv:2102.06902 [astro-ph.CO].[40] D. J. Eisenstein, D. H. Weinberg, E. Agol, H. Ai-hara, C. Allende Prieto, S. F. Anderson, J. A. Arns,E. Aubourg, S. Bailey, E. Balbinot, and et al., TheAstronomical Journal , 72 (2011).[41] K. S. Dawson, D. J. Schlegel, and et al., AJ , 10(2013), arXiv:1208.0022 [astro-ph.CO].[42] B. Reid, S. Ho, N. Padmanabhan, W. J. Percival, J. Tin-ker, R. Tojeiro, M. White, D. J. Eisenstein, C. Maras-ton, and et al., (2015), arXiv:1509.06529 [astro-ph.CO].[43] M. A. Strauss, D. H. Weinberg, R. H. Lupton,V. K. Narayanan, and et al., AJ , 1810 (2002),arXiv:astro-ph/0206225 [astro-ph].[44] A. Eggemeier, R. Scoccimarro, and R. E. Smith,Phys. Rev. D , 123514 (2019), arXiv:1812.03208[astro-ph.CO].[45] K. C. Chan, R. K. Sheth, and R. Scoccimarro,Phys. Rev. D , 103543 (2017), arXiv:1511.01909.[46] M. Musso and R. K. Sheth, MNRAS , L102 (2012),arXiv:1201.3876 [astro-ph.CO].[47] T. Lazeyras and F. Schmidt, J. Cosmology Astropart.Phys. , 041 (2019), arXiv:1904.11294 [astro-ph.CO].[48] P. McDonald, Phys. Rev. D , 103512 (2006), astro-ph/0609413.[49] V. Assassi, D. Baumann, D. Green, and M. Zaldar-riaga, J. Cosmology Astropart. Phys. , 056 (2014),arXiv:1402.5916.[50] F. Schmidt, D. Jeong, and V. Desjacques, Phys. Rev. D , 023515 (2013), arXiv:1212.0868 [astro-ph.CO].[51] M. Crocce and R. Scoccimarro, Phys. Rev. D , 063519(2006), arXiv:astro-ph/0509418.[52] F. Bernardeau, S. Colombi, E. Gazta˜naga, and R. Scoc-cimarro, Phys. Rep. , 1 (2002), astro-ph/0112551.[53] J. N. Grieb, A. G. S´anchez, S. Salazar-Albornoz,and C. Dalla Vecchia, MNRAS , 1577 (2016),arXiv:1509.04293 [astro-ph.CO].[54] D. Blas, M. Garny, and T. Konstandin, ArXiv e-prints(2013), arXiv:1309.3308 [astro-ph.CO].[55] M. Crocce and R. Scoccimarro, Phys. Rev. D , 063520(2006), arXiv:astro-ph/0509419.[56] D. J. Eisenstein, H.-J. Seo, and M. White, ApJ ,660 (2007), arXiv:astro-ph/0604361.[57] D. J. Eisenstein, H.-J. Seo, E. Sirko, and D. N. Spergel,ApJ , 675 (2007), arXiv:astro-ph/0604362 [astro-ph].[58] S. Matarrese and M. Pietroni, Modern Physics LettersA , 25 (2008), arXiv:astro-ph/0702653.[59] M. Crocce and R. Scoccimarro, Phys. Rev. D , 023533(2008), arXiv:0704.2783 [astro-ph].[60] H.-J. Seo, E. R. Siegel, D. J. Eisenstein, and M. White,ApJ , 13 (2008), arXiv:0805.0117. [61] Z. Vlah, U. Seljak, M. Yat Chu, and Y. Feng,J. Cosmology Astropart. Phys. , 057 (2016),arXiv:1509.02120 [astro-ph.CO].[62] K. Osato, T. Nishimichi, F. Bernardeau, andA. Taruya, Phys. Rev. D , 063530 (2019),arXiv:1810.10104 [astro-ph.CO].[63] D. J. Eisenstein and W. Hu, The Astrophysical Journal , 5 (1999).[64] D. Blas, M. Garny, M. M. Ivanov, and S. Sibiryakov,J. Cosmology Astropart. Phys. , 028 (2016),arXiv:1605.02149 [astro-ph.CO].[65] S. Pueblas and R. Scoccimarro, Phys. Rev. D , 043504(2009), arXiv:0809.4606.[66] J. J. M. Carrasco, M. P. Hertzberg, and L. Senatore,ArXiv e-prints (2012), arXiv:1206.2926 [astro-ph.CO].[67] D. Baumann, A. Nicolis, L. Senatore, and M. Zal-darriaga, J. Cosmology Astropart. Phys. , 051 (2012),arXiv:1004.2488 [astro-ph.CO].[68] J. J. M. Carrasco, S. Foreman, D. Green, and L. Sen-atore, ArXiv e-prints (2013), arXiv:1310.0464 [astro-ph.CO].[69] T. Baldauf, L. Mercolli, and M. Zaldarriaga,Phys. Rev. D , 123007 (2015), arXiv:1507.02256[astro-ph.CO].[70] F. Bernardeau, A. Taruya, and T. Nishimichi,Phys. Rev. D , 023502 (2014).[71] T. Nishimichi, F. Bernardeau, and A. Taruya, PhysicsLetters B , 247 (2016), arXiv:1411.2970 [astro-ph.CO].[72] T. Nishimichi, F. Bernardeau, and A. Taruya,Phys. Rev. D , 123515 (2017), arXiv:1708.08946[astro-ph.CO].[73] Planck Collaboration, A&A , A13 (2016),arXiv:1502.01589 [astro-ph.CO].[74] R. Angulo, M. Fasiello, L. Senatore, and Z. Vlah, Jour-nal of Cosmology and Astro-Particle Physics , 029(2015), arXiv:1503.08826 [astro-ph.CO].[75] F. Bernardeau, M. Crocce, and R. Scoccimarro,Phys. Rev. D , 103521 (2008), arXiv:0806.2334.[76] A. Taruya, T. Nishimichi, and F. Bernardeau,Phys. Rev. D , 083509 (2013), arXiv:1301.3624 [astro-ph.CO].[77] R. Scoccimarro and J. Frieman, ApJS , 37 (1996),arXiv:astro-ph/9509047.[78] R. Casas-Miranda, H. J. Mo, R. K. Sheth, andG. Boerner, MNRAS , 730 (2002), arXiv:astro-ph/0105008 [astro-ph].[79] A. Akbar Abolhasani, M. Mirbabayi, and E. Pa-jer, J. Cosmology Astropart. Phys. , 063 (2016),arXiv:1509.07886 [hep-th].[80] P. Catelan, F. Lucchin, S. Matarrese, and C. Porciani,MNRAS , 692 (1998), arXiv:astro-ph/9708067.[81] P. Catelan, C. Porciani, and M. Kamionkowski, MN-RAS , L39 (2000), astro-ph/0005544.[82] T. Matsubara, Phys. Rev. D , 083518 (2011),arXiv:1102.4619 [astro-ph.CO].[83] T. Lazeyras and F. Schmidt, Journal of Cosmol-ogy and Astro-Particle Physics , 008 (2018),arXiv:1712.07531 [astro-ph.CO].[84] M. M. Abidi and T. Baldauf, J. Cosmology Astropart.Phys. , 029 (2018), arXiv:1802.07622.[85] R. K. Sheth, K. C. Chan, and R. Scoccimarro,Phys. Rev. D , 083002 (2013), arXiv:1207.7117 [astro-ph.CO]. [86] V. Springel, N. Yoshida, and S. D. White, New Astron-omy , 79–117 (2001).[87] V. Springel, Monthly Notices of the Royal AstronomicalSociety , 1105–1134 (2005).[88] A. Lewis and S. Bridle, Phys. Rev. D , 103511 (2002),arXiv:astro-ph/0205436 [astro-ph].[89] M. Crocce, S. Pueblas, and R. Scoccimarro, MNRAS , 369 (2006), arXiv:astro-ph/0606505.[90] S. Alam, M. Ata, S. Bailey, F. Beutler, D. Bizyaev,J. A. Blazek, A. S. Bolton, J. R. Brownstein, A. Burden,C.-H. Chuang, J. Comparat, A. J. Cuesta, K. S. Daw-son, D. J. Eisenstein, S. Escoffier, H. Gil-Mar´ın, J. N.Grieb, N. Hand, S. Ho, K. Kinemuchi, D. Kirkby, F. Ki-taura, E. Malanushenko, V. Malanushenko, C. Maras-ton, C. K. McBride, R. C. Nichol, M. D. Olm-stead, D. Oravetz, N. Padmanabhan, N. Palanque-Delabrouille, K. Pan, M. Pellejero-Ibanez, W. J. Per-cival, P. Petitjean, F. Prada, A. M. Price-Whelan,B. A. Reid, S. A. Rodr´ıguez-Torres, N. A. Roe, A. J.Ross, N. P. Ross, G. Rossi, J. A. Rubi˜no-Mart´ın,S. Saito, S. Salazar-Albornoz, L. Samushia, A. G.S´anchez, S. Satpathy, D. J. Schlegel, D. P. Schneider,C. G. Sc´occola, H.-J. Seo, E. S. Sheldon, A. Simmons,A. Slosar, M. A. Strauss, M. E. C. Swanson, D. Thomas,J. L. Tinker, R. Tojeiro, M. V. Maga˜na, J. A. Vazquez,L. Verde, D. A. Wake, Y. Wang, D. H. Weinberg,M. White, W. M. Wood-Vasey, C. Y`eche, I. Zehavi,Z. Zhai, and G.-B. Zhao, MNRAS , 2617 (2017),arXiv:1607.03155 [astro-ph.CO].[91] C. McBride, A. Berlind, R. Scoccimarro, R. Wech-sler, M. Busha, J. Gardner, and F. van den Bosch,in American Astronomical Society Meeting Abstracts , American Astronomical Society Meeting Ab-stracts, Vol. 213 (2009) p. 425.06.[92] M. Sinha, A. A. Berlind, C. K. McBride, R. Scocci-marro, J. A. Piscionere, and B. D. Wibking, MNRAS , 1042 (2018), arXiv:1708.04892 [astro-ph.CO].[93] U. Seljak and M. Zaldarriaga, ApJ , 437 (1996).[94] E. Sefusatti, M. Crocce, R. Scoccimarro, andH. M. P. Couchman, MNRAS , 3624 (2016),arXiv:1512.07295 [astro-ph.CO].[95] H. A. Feldman, N. Kaiser, and J. A. Peacock, ApJ ,23 (1994), arXiv:astro-ph/9304022.[96] M. Tegmark, Phys. Rev. Lett. , 3806 (1997),arXiv:astro-ph/9706198 [astro-ph].[97] A. Oddo, E. Sefusatti, C. Porciani, P. Monaco, andA. G. S´anchez, J. Cosmology Astropart. Phys. ,056 (2020), arXiv:1908.01774 [astro-ph.CO].[98] A. Gelman and D. B. Rubin, Statistical Science , 457(1992).[99] A. Lewis, (2019), arXiv:1910.13970 [astro-ph.IM].[100] A. G. S´anchez, Phys. Rev. D , 123511 (2020).[101] T. Abbott, F. Abdalla, A. Alarcon, J. Aleksi´c, S. Al-lam, S. Allen, A. Amara, J. Annis, J. Asorey, S. Avila,and et al., Physical Review D (2018), 10.1103/phys-revd.98.043526.[102] C. Heymans, T. Tr¨oster, M. Asgari, C. Blake, H. Hilde-brandt, B. Joachimi, K. Kuijken, C.-A. Lin, A. G.S´anchez, J. L. van den Busch, A. H. Wright, A. Amon,M. Bilicki, J. de Jong, M. Crocce, A. Dvornik, T. Erben,M. C. Fortuna, F. Getman, B. Giblin, K. Glazebrook,H. Hoekstra, S. Joudaki, A. Kannawadi, F. K¨ohlinger,C. Lidman, L. Miller, N. R. Napolitano, D. Parkinson,P. Schneider, H. Shan, E. Valentijn, G. V. Kleijn, and C. Wolf, (2020), arXiv:2007.15632 [astro-ph.CO].[103] S. Pandey, E. Krause, B. Jain, N. MacCrann, J. Blazek,M. Crocce, J. DeRose, X. Fang, I. Ferrero, O. Friedrich,M. Aguena, S. Allam, J. Annis, S. Avila, G. M. Bern-stein, D. Brooks, D. L. Burke, A. Carnero Rosell,M. Carrasco Kind, J. Carretero, M. Costanzi, L. N.da Costa, J. De Vicente, S. Desai, J. Elvin-Poole,S. Everett, P. Fosalba, J. Frieman, J. Garc´ıa-Bellido,D. Gruen, R. A. Gruendl, J. Gschwend, G. Gutierrez,K. Honscheid, K. Kuehn, N. Kuropatkin, M. A. G.Maia, J. L. Marshall, F. Menanteau, R. Miquel, A. Palmese, F. Paz-Chinch´on, A. A. Plazas, A. Rood-man, E. Sanchez, V. Scarpine, M. Schubnell, S. Ser-rano, I. Sevilla-Noarbe, M. Smith, M. Soares-Santos,E. Suchyta, M. E. C. Swanson, G. Tarle, and J. Weller(DES Collaboration), Phys. Rev. D , 123522 (2020).[104] A. Blanchard, S. Camera, C. Carbone, V. F. Cardone,S. Casas, S. Clesse, S. Ili´c, M. Kilbinger, T. Kitching,and et al., A&A , A191 (2020).[105] J. D. Hunter, Computing in Science & Engineering9