© ESO 2018

1 Introduction

The standard paradigm of cosmology describes the large-scale distribution of matter and galaxies in an expanding Universe (Dodelson 2003, and references therein). Strongly supported by observations, this model assumes a statistically homogeneous and isotropic Universe with cold dark matter (CDM) as the dominating form of matter. Matter in total has the mean density Ω_m ≈ 0.3 of which ordinary baryonic matter is just Ω_b ≈ 0.05; as usual densities are in units of the critical density (or its energy equivalent). The largest fraction Ω_Λ ≈ 0.7 in the cosmological energy density is given by a cosmological constant Λ or so-called dark energy, resulting in a flat or approximately flat background geometry with curvature parameter K = 0 (Einstein 1917; Planck Collaboration XIII 2016, and references therein). The exact physical nature of dark matter is unknown but its presence is consistently inferred through visible tracers from galactic to cosmological scales at different epochs in the cosmic history (Bertone et al. 2005, for a review). In particular the coherent shear of distant galaxy images (background sources) by the tidal gravitational field of intervening matter gives direct evidence for the (projected) density field of dark matter (Clowe et al. 2004). The basic physics of galaxy formation inside dark-matter halos and the galaxy evolution seems to be identified and reasonably well matched by observations, although various processes, such as star formation and galaxy-gas feedback, are still not well understood or worked out in detail (Mo et al. 2010). Ultimately, the ability of the ΛCDM model to quantitatively describe the observed richness of galaxy properties from initial conditions will be a crucial validation test.

One galaxy property of importance is spatial distribution. Galaxies are known to be differently distributed than matter in general; they are so-called biased tracers of the matter density field (Kaiser 1984). The details of the biasing mechanism are related to galaxy physics (Springel et al. 2018, 2005; Weinberg et al. 2004; Somerville et al. 2001; Benson et al. 2000; Peacock & Smith 2000). An observed galaxy bias for different galaxy types and redshifts consequently provides input and tests for galaxy models. Additionally, the measurement of galaxy bias is practical for studies that rely on fiducial values for the biasing of a particular galaxy sample or on the observational support for a high galaxy-matter correlation on particular spatial scales (e.g. van Uitert et al. 2018; Hildebrandt et al. 2013; Simon 2013; Mehta et al. 2011; Reyes et al. 2010; Baldauf et al. 2010). In this context, we investigate the prospects of weak gravitational lensing to measure the galaxy bias (e.g. Kilbinger 2015; Munshi et al. 2008; Schneider et al. 2006, for a review).

There are clearly various ways to express the statistical relationship between the galaxy and matter distribution, which both can be seen as realisations of statistically homogeneous and isotropic random fields (Desjacques et al. 2016). With focus on second-order statistics we use the common parameterisation in Tegmark & Bromley (1999). This defines galaxy bias in terms of auto- and cross-correlation power spectra of the random fields for a given wave number k: (i) a bias factor b(k) for the relative strength between galaxy and matter clustering; and (ii) a factor r(k) for the galaxy-matter correlation. The second-order biasing functions can be constrained by combining galaxy clustering with cosmic-shear information in lensing surveys (Foreman et al. 2016; Cacciato et al. 2012; Simon 2012; Pen et al. 2003). In applications of these techniques, galaxy biasing is then known to depend on galaxy type, physical scale, and redshift, thus reflecting interesting galaxy physics (Chang et al. 2016; Buddendiek et al. 2016; Pujol et al. 2016; Prat et al. 2018; Comparat et al. 2013; Simon et al. 2013, 2007; Jullo et al. 2012; Pen et al. 2003; Hoekstra et al. 2002).

Our interest here is the quality of lensing measurements of galaxy bias. For this purpose, we focus on the method by van Waerbeke (1998) and Schneider (1998), first applied in Hoekstra et al. (2001) and Hoekstra et al. (2002), where one defines relative aperture measures of the galaxy number density and the lensing mass to observe b(k) and r(k) as projections on the sky, averaged in bands of radial and transverse direction. The advantage of this method is its model independence apart from a cosmology-dependent normalisation. As an improvement, we define a new procedure to de-project the lensing measurements of the projected biasing functions, giving direct estimates of b(k) and r(k) for a selected galaxy population. In addition, we account for the intrinsic alignment of source galaxies that are utilised in the lensing analysis (Kirk et al. 2015). To eventually assess the accuracy and precision of our de-projection technique, we compare the results to the true biasing functions for various galaxy samples in a simulated survey of about 1000 deg² angular area, constructed with a semi-analytic galaxy model by Henriques et al. (2015, H15 hereafter), and data from the Millennium Simulation (Springel et al. 2005). To this end, a large part of this paper deals with the construction of flexible template models of b(k) and r(k) that we forward-fit to the relative aperture measures. These templates are based on a flexible halo-model prescription, which additionally allows us a physical interpretation of the biasing functions (Cooray & Sheth 2002). Some time is therefore also spent on a discussion of the scale dependence of galaxy bias, which will be eminent in future applications of our technique.

The structure of this paper is as follows. In Sect. 2, we describe the construction of data for a mock lensing survey to which we apply our de-projection technique. With respect to number densities of lens and source galaxies on the sky, the mock data are similar to realistic galaxy samples in the Canada-France-Hawaii Telescope Lensing Survey (CFHTLenS; Heymans et al. 2012). We increase the simulated survey area, however, to ~ 1000deg² in order to assess the quality of our methodology for current (ground-based) surveys in future applications. In Sect. 3, we revise the relation of the spatial biasing functions to their projected counterparts, which are observable through the aperture statistics. This section also adds to the technique of Hoekstra et al. (2002) new, potentially relevant higher-order corrections in the lensing formalism. It also incorporates a treatment of the intrinsic alignment of sources into the aperture statistics. Section 4 derives our template models of the spatial biasing functions, applied for a de-projection. Section 5 summarises the template parameters and explores their impact on the scale dependence of galaxy bias. The methodological details for the statistical inference of b(k) and r(k) from noisy measurements are presented in Sect. 6. We apply this inference technique to the mock data in the results described in Sect. 7 and assess its accuracy, precision, and robustness. As a first demonstration, we apply our technique to previous measurements in Simon et al. (2007). We finally discuss our results in Sect. 8.

2 Data

This section details our mock data, that is lens and source catalogues, to which we apply our de-projection technique in the following sections. A reader more interested in the method details for the recovery of galaxy bias with lensing data could proceed to the next sections.

2.1 Samples of lens galaxies

Our galaxy samples use a semi-analytic model (SAM) according to H15, which is implemented on the Millennium Simulation (Springel et al. 2006). These SAMs are the H15 mocks that are also used in Saghiha et al. (2017). The Millennium Simulation (MS) is an N-body simulation of the CDM density field inside a comoving cubic volume of 500 h⁻¹ Mpc side length, and it has a spatial resolution of 5 h⁻¹ kpc sampled by 10¹⁰ mass particles. The fiducial cosmology of the MS has the density parameters Ω_m = 0.25 = 1 −Ω_Λ and Ω_b = 0.045, σ₈ = 0.9 for the normalisation of the linear matter power spectrum, a Hubble parameter H₀ = 100 h km s⁻¹ Mpc⁻¹ with h = 0.73, and a spectral index for the primordial matter power spectrum of n_spec = 1.0. All density parameters are in units of the critical density ${\overline{\rho }}_{\text{crit}}=3H_0^2/8\pi G_{\text{N}}$, where G_N denotes Newton’s constant of gravity. The galaxy mocks are constructed by populating dark matter halos in the simulation based on the merger history of halos and in accordance with the SAM details. We project the positions of the SAMs inside 64 independent light cones onto a 4 × 4 deg² piece of sky. The resulting total survey area is hence 1024 deg².

We then select galaxies from the mocks to emulate the selection in redshift and stellar mass in Simon et al. (2013), SES13 henceforth. Details on the emulation process can be found in Saghiha et al. (2017). We give only a brief summary here. The mock galaxy and source samples are constructed to be compatible with those in recent lensing studies, dealing with data from the CFHTLenS (Saghiha et al. 2017; Velander et al. 2014; Erben et al. 2013; Heymans et al. 2012). Our selection proceeds in two steps. First, we split the galaxy catalogues into stellar mass, including emulated measurement errors, and i′ -band brightness to produce the stellar-mass samples SM1 to SM6; the photometry uses the AB-magnitude system. Second, we randomly discard galaxies in each stellar-mass sample to obtain a redshift distribution that is comparable to a given target distribution. As targets, we employ the photometric redshift bins “low-z” and “high-z” in SES13, which are the redshift distributions in CFHTLenS after a cut in photometric redshift z_p. The low-z bin applies 0.2 ≤ z_p < 0.44, and the high-z bin applies 0.44 ≤ z_p < 0.6. Figure 5 in SES13 gives the different target distributions. Our selection criteria for SM1 to SM6 are listed in Table 1. We note here that randomly removing galaxies at redshift z adds shot noise but does not change the matter-galaxy correlations and the (shot-noise corrected) galaxy clustering.

In addition to SM1-6, we define two more samples, RED and BLUE, based on the characteristic bimodal distribution of u − r colours (Table 1). Both samples initially consist of all galaxies in SM1 to SM6 but are then split depending on the u − r colours of galaxies: the division is at (u − r)(z) = 1.93 z + 1.85, which varies with z to account for the reddening with redshift. We crudely found (u − r)(z) by identifying by eye the mid-points ${(u-r)}_i$ between the red and blue mode in u − r histograms of CFHTLenS¹ SM1-6 galaxies in four photometric-redshift bins with means {z_i} = {0.25, 0.35, 0.45, 0.55} and width Δz = 0.1 (Hildebrandt et al. 2012). Then we fit a straight line to the four empirical data points $\{(z_i,{(u-r)}_i)\}$ and obtain the above (u − r)(z) as best fit. To split the mocks, we identify the precise redshifts z in H15 with the photometric redshifts z_p in CFHTLenS which, for the scope of this work, is a sufficient approximation. Similar to the previous stellar-mass samples, we combine theredshift posteriors of all CFHTLenS-galaxies RED or BLUE to define the target distributions for our corresponding mock samples.

For the following galaxy-bias analysis, we estimate the probability density function (PDF) p_d (z) of each galaxy sample from the mock catalogues in the foregoing step. Simply using histograms of the sample redshifts may seem like a good idea but is, in fact, problematic because the histograms depend on the adoptedbinning. This is especially relevant for the prediction of galaxy clustering that depends on $p_{\text{d}}^2(z)$ (see Eq. (22)). Instead, we fit for p_d(z) a smooth four-parameter Gram-Charlier series $$p_{\text{d}}(z|\lambda ,\Theta )=\lambda {\text{e}}^{-\frac{x^2}{2}}\left(1+\frac{s}{6}{H}_3(x)+\frac{k}{24}{H}_4(x)\right);x=\frac{z-\overline{z}}{{\sigma }_{\text{z}}},$$(1)

with the Hermite polynomials H₃(x) = x³ − 3x and H₄(x) = x⁴ − 6x² + 3 to a mock sample {z_i: i = 1…n} of n galaxy redshifts; λ is a normalisation constant that depends on the parameter combination $\Theta =(\overline{z},{\sigma }_{\text{z}},s,k)$ and is defined by $${\displaystyle {\int }_0^{\infty }\text{d}}z{p}_{\text{d}}(z|\lambda ,\Theta )=1.$$(2)

For an estimate $\widehat{\Theta }$ of the parameters Θ, we maximise the log-likelihood $$\ln \mathcal{L}(\Theta )={\displaystyle \sum _{i=1}^n\ln }{p}_{\text{d}}(z_i|\lambda ,\Theta )$$(3)

with respect to Θ. This procedure selects the PDF $p_{\text{d}}(z|\lambda ,\widehat{\Theta })$ that is closest to the sample distribution of redshifts z_i in a Kullback-Leibler sense (Knight 1999). The mean $\overline{z}$ and variance ${\sigma }_z^2$ in the fit matches that of the redshift distribution in the mock lens sample. The resulting densities for all our lens samples are shown in the two top panels of Fig. 1.

Fig. 1 Models of the probability densities pd(z) of galaxy redshifts in our lens samples SM1 to SM6, RED and BLUE (two top panels), and the density ps (z) of the source sample (bottom panel).

2.2 Shear catalogues

For mock source catalogues based on the MS data, we construct lensing data by means of multiple-lens-plane ray tracing as described in Hilbert et al. (2009). The ray tracing produces the lensing convergence κ(θ |z_s) and shear distortion γ(θ|z_s) for 4096² line-of-sight directions θ on 64 regular angular grids and a sequence of n_s = 31 source redshifts z_s,i between z_s = 0 and z_s = 2; we denote by Δz_i = z_s,i+1 − z_s,i the difference between neighbouring source redshifts. Each grid covers a solid angle of Ω =4 × 4 deg². For each grid, we then compute the average convergence for sources with redshift PDF p_s (z) by $$\kappa (\theta )=\frac{{\displaystyle {\sum }_{i=1}^{n_{\text{s}}}{p}_{\text{s}}}(z_{\text{s},i})\Delta z_i\kappa (\theta |z_{\text{s}.i})}{{\displaystyle {\sum }_{i=1}^{n_{\text{s}}}{p}_{\text{s}}}(z_{\text{s},i})\Delta z_i},$$(4)

and the average shear γ(θ) from the sequence γ(θ|z_s) accordingly. For p_s(z), we employ the estimated PDF of CFHTLenS sources that is selected through i′ < 24.7 and 0.65 ≤ z_p < 1.2, weighted by their shear-measurement error (SES13; see the bottom panel in Fig. 1). The mean redshift of sources is $\overline{z}\approx 0.93$. To assign source positions on the sky, we uniform-randomly pick a sample {θ_i: i = 1…n} of positions for each grid; the amount of positions is $n=\Omega {\overline{n}}_{\text{s}}$ for a number density of ${\overline{n}}_{\text{s}}=5 {\text{arcmin}}^{-2}$ sources, which roughly equals the effective number density of sources in SES13.

Depending on the type of our lensing analysis, we assign a source at θ_i of one of the following three values for the simulated sheared ellipticity ϵ_i: (i) ϵ_i = γ(θ_i) for source without shape noise; (ii) ϵ_i = A(γ(θ_i), ϵ_s) for noisy sources with shear; and (iii) ${ϵ}_i=A\left(g_i,{ϵ}_{\text{s}}\right)$ for noisy sources with reduced shear g_i = γ(θ_i)∕[1 − κ(θ_i)]. We define here by $A(x,y):=(x+y){(1+x{y}^{\ast })}^{-1}$ the conformal mapping of two complex numbers x and y, and by ϵ_s a random shape noise drawn from a bivariate, truncated Gaussian PDF with zero mean, 1D dispersion σ_ϵ = 0.3, and an exclusion of values beyond |ϵ_s|≥ 1.

2.3 Power spectra

We obtain the true spatial galaxy-galaxy, galaxy-matter, and matter-matter power spectra for all galaxy samples at a given simulation snapshot with fast Fourier transform (FFT) methods. For a choice of pair of tracers (i.e. simulation matter particles or galaxies from different samples) in a snapshot, we compute a series of raw power spectra by “chaining the power” (Smith et al. 2003). We cover the whole simulation volume as well as smaller sub-volumes (by a factor 4³ to 256³, into which the whole box is folded) by regular meshes of 512³ points (providing a spatial resolution from ~1 h⁻¹ Mpc for the coarsest mesh to ~5 h⁻¹ kpc for the finest mesh). We project the tracers onto these meshes using clouds-in-cells (CIC) assignment (Hockney & Eastwood 1981).

We FFT-transform the meshes, record their raw power spectra, apply a shot-noise correction (except for cross-spectra), a deconvolution to correct for the smoothing by the CIC assignment, and an iterative alias correction (similar to what is described in Jing 2005). From these power spectra, we discard small scales beyond half their Nyquist frequency as well as large scales that are already covered by a coarser mesh, and combine them into a single power spectrum covering a range of scales from modes ~ 0.01 h Mpc⁻¹ to modes ~ 100 h Mpc⁻¹.

The composite power spectra are then used as input to estimate alias corrections for the partial power spectra from the individual meshes with different resolutions, and the process is repeated until convergence. From the resulting power spectra, we then compute the true biasing functions, Eq. (10), which we compare to our lensing-based reconstructions in Sect. 7.

3 Projected biasing functions as observed with lensing techniques

The combination of suitable statistics for galaxy clustering, galaxy-galaxy lensing, and cosmic-shear correlations on the sky allows us to infer, without a concrete physical model, the z-averaged spatial biasing functions b(k) and r(k) as projections b_2D(θ_ap) and r_2D (θ_ap) for varying angular scales θ_ap. Later on, we forward-fit templates of spatial biasing functions to these projected functions to perform a stable de-projection. We summarise here the relation between (b(k), r(k)) and the observable ratio-statistics (b_2D(θ_ap), r_2D(θ_ap)). We include corrections to the first-order Born approximation for galaxy-galaxy lensing and galaxy clustering, and corrections for the intrinsic alignment of sources.

3.1 Spatial biasing functions

We define galaxy bias in terms of two biasing functions b(k) and r(k) for a given spatial scale 2π k⁻¹ or wave number k in the following way. Let δ(x) in $\rho (x)=\overline{\rho }[1+\delta (x)]$ be the density fluctuations at position x of a random density field ρ(x), and $\overline{\rho }$ denotes the mean density. A density field is either the matter density ρ_m(x) or the galaxy number density n_g(x) with densitycontrasts δ_m(x) and δ_g(x), respectively. We determine the fluctuation amplitude for a density mode k by the Fourier transform of δ(x), $$\tilde{\delta }(k)={\displaystyle \int {\text{d}}^3}\text{}x\delta (x){\text{e}}^{-\text{i}x\cdot k}.$$(5)

All information on the two-point correlations of $\tilde{\delta }(k)$ is contained in the power spectrum P(k) defined through the second-order correlation function of modes, $$\langle \tilde{\delta }(k)\tilde{\delta }(k^{\prime })\rangle ={(2\pi )}^3{\delta }_{\text{D}}(k+k^{\prime })P(k),$$(6)

where k = |k| is the scalar wave number and δ_D(s) is the Dirac Delta distribution. Specifically, we utilise three kinds of power spectra, $$\begin{array}{lll}\langle {\tilde{\delta }}_{\text{m}}(k){\tilde{\delta }}_{\text{m}}(k^{\prime })\rangle \hfill & =\hfill & {(2\pi )}^3{\delta }_{\text{D}}(k+k^{\prime })P_{\text{m}}(k);\hfill \\ \langle {\tilde{\delta }}_{\text{m}}(k){\tilde{\delta }}_{\text{g}}(k^{\prime })\rangle \hfill & =\hfill & {(2\pi )}^3{\delta }_{\text{D}}(k+k^{\prime })P_{\text{gm}}(k);\hfill \\ \langle {\tilde{\delta }}_{\text{g}}(k){\tilde{\delta }}_{\text{g}}(k^{\prime })\rangle \hfill & =\hfill & {(2\pi )}^3{\delta }_{\text{D}}(k+k^{\prime })\left(P_{\text{g}}(k)+{\overline{n}}_{\text{g}}^{-1}\right),\hfill \end{array}$$(7) (8) (9)

namely the matter power spectrum P_m(k), the galaxy-matter cross-power spectrum P_gm(k), and the galaxy power spectrum P_g(k). The latter subtracts the shot noise ${\overline{n}}_{\text{g}}^{-1}$ from the galaxy power spectrum by definition. In contrast to the smooth matter density, the galaxy number density is subject to shot noise because it consists of a finite number of discrete points that make up the number density field. Traditionally, the definition of P_g (k) assumes a Poisson process for the shot noise in the definition of P_g(k) (Peebles 1980).

The biasing functions (of the second order) express galaxy bias in terms of ratios of the foregoing power spectra, $$b(k):=\sqrt{\frac{P_{\text{g}}(k)}{P_{\text{m}}(k)}} ; r(k):=\frac{P_{\text{gm}}(k)}{\sqrt{P_{\text{g}}(k)P_{\text{m}}(k)}}.$$(10)

Galaxies that sample the matter density by a Poisson process have b(k) = r(k) = 1 for all scales k and are dubbed “unbiased”; for b(k) > 1, we find that galaxies cluster stronger than matter at scale k, and vice versa for b(k) < 1; a decorrelation of r(k)≠1 indicates either stochastic bias, non-linear bias, a sampling process that is non-Poisson, or combinations of these cases (Dekel & Lahav 1999; Guzik & Seljak 2001).

3.2 Aperture statistics and galaxy-bias normalisation

The projected biasing functions b(k) and r(k) are observable by taking ratios of the (co-)variances of the aperture mass and aperture number count of galaxies (van Waerbeke 1998; Schneider 1998). To see this, let ${\kappa }_{\text{g}}(\theta )=N_{\text{g}}(\theta )/{\overline{N}}_{\text{g}}-1$ be the density contrast of the number density of galaxies N_g(θ) on the sky in the direction θ, and ${\overline{N}}_{\text{g}}=\langle N_{\text{g}}(\theta )\rangle$ be their mean number density. We define the aperture number count of N_g (θ) for an angular scale θ_ap at position θ by $$\mathcal{N}({\theta }_{\text{ap}};\theta )={\displaystyle \int {\text{d}}^2}{\theta }^{\prime }U(|{\theta }^{\prime }|;{\theta }_{\text{ap}}){\kappa }_{\text{g}}({\theta }^{\prime }+\theta ),$$(11)

where $$U(\theta ;{\theta }_{\text{ap}})=\frac{1}{{\theta }_{\text{ap}}^2}u(\theta {\theta }_{\text{ap}}^{-1}); u(x)=\frac{9}{\pi }(1-x^2)\left(\frac{1}{3}-x^2\right)\text{\hspace{0.17em}H}(1-x)$$(12)

is the aperture filter of the density field, and H(x) is the Heaviside step function of our polynomial filter profile u(x). The aperture filter is compensated, that is ${\displaystyle {\int }_0^{\infty }\text{d}}xxu(x)=0$. Similarly for the (average) lensing convergence κ(θ) of sources in direction θ, the aperture mass is given by $$M_{\text{ap}}({\theta }_{\text{ap}};\theta )={\displaystyle \int {\text{d}}^2}{\theta }^{\prime }U(|{\theta }^{\prime }|;{\theta }_{\text{ap}})\kappa ({\theta }^{\prime }+\theta ).$$(13)

The aperture statistics consider the variances $\langle {\mathcal{N}}^2\rangle ({\theta }_{\text{ap}})$ and $\langle M_{\text{ap}}^2\rangle ({\theta }_{\text{ap}})$ of $\mathcal{N}({\theta }_{\text{ap}};\theta )$ and M_ap(θ_ap;θ), respectively, across the sky, as well as their co-variance $\langle \mathcal{N}{M}_{\text{ap}}\rangle ({\theta }_{\text{ap}})$ at zero lag.

From these observable aperture statistics, we obtain the galaxy bias factor b_2D (θ_ap) and correlation factor r_2D(θ_ap) through the ratios $$\begin{array}{lll}{b}_{2\text{D}}({\theta }_{\text{ap}})\hfill & =\hfill & \sqrt{\frac{\langle {\mathcal{N}}^2\rangle ({\theta }_{\text{ap}})}{\langle M_{\text{ap}}^2\rangle ({\theta }_{\text{ap}})}}\times f_{\text{b}}({\theta }_{\text{ap}}),\hfill \\ r_{2\text{D}}({\theta }_{\text{ap}})\hfill & =\hfill & \frac{\langle \mathcal{N}{M}_{\text{ap}}\rangle ({\theta }_{\text{ap}})}{\sqrt{\langle {\mathcal{N}}^2\rangle ({\theta }_{\text{ap}})\langle M_{\text{ap}}^2\rangle ({\theta }_{\text{ap}})}}\times f_{\text{r}}({\theta }_{\text{ap}}),\hfill \end{array}$$(14) (15)

where $$\begin{array}{lll}{f}_{\text{b}}({\theta }_{\text{ap}})\hfill & :=\hfill & \sqrt{\frac{{\langle M_{\text{ap}}^2\rangle }_{\text{th}}({\theta }_{\text{ap}})}{{\langle {\mathcal{N}}^2\rangle }_{\text{th}}({\theta }_{\text{ap}};1)}},\hfill \\ f_{\text{r}}({\theta }_{\text{ap}})\hfill & :=\hfill & \frac{\sqrt{{\langle M_{\text{ap}}^2\rangle }_{\text{th}}({\theta }_{\text{ap}}){\langle {\mathcal{N}}^2\rangle }_{\text{th}}({\theta }_{\text{ap}};1)}}{{\langle \mathcal{N}{M}_{\text{ap}}\rangle }_{\text{th}}({\theta }_{\text{ap}};1)}\hfill \end{array}$$(16) (17)

normalise the statistics according to a fiducial cosmology, that means the aperture statistics with subscript “th” as in ${\langle M_{\text{ap}}^2\rangle }_{\text{th}}({\theta }_j)$ denote the expected (co-)variance for a fiducial model. The normalisation is chosen such that we have b_2D (θ_ap) = r_2D(θ_ap) = 1 for unbiased galaxies given the distributions of lenses and sources with distance χ as in the survey, hence the “(θ_ap;1)” in the arguments of the normalisation. The normalisation functions f_f and f_b are typically weakly varying withangular scale θ_ap (Hoekstra et al. 2002). In addition, they depend weakly on the fiducial matter power spectrum P_m (k;z); they are even invariant with respect to an amplitude change P_m(k;z)↦υ P_m(k;z) with some number υ > 0. We explore the dependence on the fiducial cosmology quantitatively in Sect. 7.3.

For this study, we assume that the distance distribution of lenses is sufficiently narrow, which means that the bias evolution in the lens sample is negligible. We therefore skip the argument χ in b(k; χ) and r(k; χ), and we use a b(k) and r(k) independent of χ for average biasing functions instead.

The relation between (b(k), r(k)) and (b_2D(θ_ap), r_2D(θ_ap)) is discussed in the following. Let p_d(χ) dχ and p_s(χ) dχ be the probability to find a lens or source galaxy, respectively, at comoving distance [χ, χ + dχ). The matter power spectrum at distance χ shall be P_m(k;χ), and $k_{\mathcal{l}}^{\chi }:=(\mathcal{l}+0.5)/f_K(\chi )$ is a shorthand for the transverse spatial wave number k at distance χ that corresponds to the angular wave number ℓ. The function f_K(χ) denotes the comoving angular-diameter distance in the given fiducial cosmological model. The additive constant 0.5 in $k_{\mathcal{l}}^{\chi }$ applies a correction to the standard Limber approximation on the flat sky, which gives more accurate results for large angular scales (Kilbinger et al. 2017; Loverde & Afshordi 2008). According to theory, the aperture statistics are then $$\begin{array}{l}{\langle {\mathcal{N}}^2\rangle }_{\text{th}}({\theta }_{\text{ap}};b)=2\pi {\displaystyle \underset{0}{\overset{\infty }{\int }}\text{d}}\mathcal{l}\mathcal{l}{P}_{\text{n}}(\mathcal{l};b){\left[I(\mathcal{l}{\theta }_{\text{ap}})\right]}^2\text{},\\ {\langle \mathcal{N}{M}_{\text{ap}}\rangle }_{\text{th}}({\theta }_{\text{ap}};b,r)=2\pi {\displaystyle \underset{0}{\overset{\infty }{\int }}\text{d}}\mathcal{l}\mathcal{l}{P}_{\text{n}\kappa }(\mathcal{l};b,r){\left[I(\mathcal{l}{\theta }_{\text{ap}})\right]}^2\text{},\\ {\langle M_{\text{ap}}^2\rangle }_{\text{th}}({\theta }_{\text{ap}})=2\pi {\displaystyle \underset{0}{\overset{\infty }{\int }}\text{d}}\mathcal{l}\mathcal{l}{P}_{\kappa }(\mathcal{l}){\left[I(\mathcal{l}{\theta }_{\text{ap}})\right]}^2\text{},\end{array}$$(18) (19) (20)

with the angular band-pass filter $$I(x):={\displaystyle \underset{0}{\overset{\infty }{\int }}\text{d}}ssu(s){\text{J}}_0(sx)=\frac{12}{\pi }\frac{{\text{J}}_4(x)}{x^2},$$(21)

the angular power spectrum of the galaxy clustering $$P_{\text{n}}(\mathcal{l};b)={\displaystyle \underset{0}{\overset{{\chi }_{\text{h}}}{\int }}\frac{\text{d}\chi p_{\text{d}}^2(\chi )}{f_K^2(\chi )}}{b}^2(k_{\mathcal{l}}^{\chi })P_{\text{m}}\left(k_{\mathcal{l}}^{\chi };\chi \right),$$(22)

the galaxy-convergence cross-power $$\begin{array}{l}{P}_{\text{n}\kappa }(\mathcal{l};b,r)\\ =\frac{3H_0^2{\Omega }_{\text{m}}}{2c^2}{\displaystyle \underset{0}{\overset{{\chi }_{\text{h}}}{\int }}\frac{\text{d}\chi p_{\text{d}}(\chi )g_{\text{s}}(\chi )}{a(\chi )f_K(\chi )}}b(k_{\mathcal{l}}^{\chi })r(k_{\mathcal{l}}^{\chi })P_{\text{m}}\left(k_{\mathcal{l}}^{\chi };\chi \right),\end{array}$$(23)

and the convergence power-spectrum $$P_{\kappa }(\mathcal{l})=\frac{9H_0^4{\Omega }_{\text{m}}^2}{4c^4}{\displaystyle \underset{0}{\overset{{\chi }_{\text{h}}}{\int }}\frac{\text{d}\chi g_{\text{s}}^2(\chi )}{a^2(\chi )}}{P}_{\text{m}}\left(k_{\mathcal{l}}^{\chi };\chi \right),$$(24)

all in the Born and Limber approximation. In the integrals, we use the lensing kernel $$g_{\text{s}}(\chi )={\displaystyle \underset{\chi }{\overset{{\chi }_{\text{h}}}{\int }}\text{d}}{\chi }^{\prime }{p}_{\text{s}}({\chi }^{\prime })\frac{f_K({\chi }^{\prime }-\chi )}{f_K({\chi }^{\prime })},$$(25)

the scale factor a(χ) at distance χ, the maximum distance χ_h of a source, and the nth-order Bessel function J_n(x) of the first kind. By c we denote the vacuum speed of light. The power spectra and aperture statistics depend on specific biasing functions as indicated by the b and r in the arguments. For given biasing functions b(k) and r(k), we obtain the normalised galaxy bias inside apertures therefore through $$\begin{array}{l}{b}_{2\text{D}}({\theta }_{\text{ap}};b)=\sqrt{\frac{{\langle {\mathcal{N}}^2\rangle }_{\text{th}}({\theta }_{\text{ap}};b)}{{\langle {\mathcal{N}}^2\rangle }_{\text{th}}({\theta }_{\text{ap}};1)}},\\ r_{2\text{D}}({\theta }_{\text{ap}};b,r)=\frac{1}{b_{2\text{D}}({\theta }_{\text{ap}};b)}\frac{{\langle \mathcal{N}{M}_{\text{ap}}\rangle }_{\text{th}}({\theta }_{\text{ap}};b,r)}{{\langle \mathcal{N}{M}_{\text{ap}}\rangle }_{\text{th}}({\theta }_{\text{ap}};1)},\end{array}$$(26) (27)

which can be compared to measurements of Eqs. (14) and (15).

3.3 Intrinsic alignment of sources

Recent studies of cosmic shear find evidence for an alignment of intrinsic source ellipticities that contribute to the shear-correlation functions (Hildebrandt et al. 2017; Abbott et al. 2016; Troxel & Ishak 2015; Heymans et al. 2013; Joachimi et al. 2011; Mandelbaum et al. 2006). These contributions produce systematic errors in the reconstruction of b(k) and r(k) if not included in their normalisation f_b and f_r. Relevant are “II”-correlations between intrinsic shapes of sources in $\langle M_{\text{ap}}^2\rangle$ and “GI”-correlations between shear and intrinsic shapes in both $\langle \mathcal{N}{M}_{\text{ap}}\rangle$ and $\langle M_{\text{ap}}^2\rangle$. The GI term in $\langle \mathcal{N}{M}_{\text{ap}}\rangle$ can be suppressed by minimising the redshift overlap between lenses and sources. Likewise, the II term is suppressed by a broad redshift distribution of sources which, however, increases the GI amplitude. The amplitudes of II and GI also vary with galaxy type and luminosity of the sources (Joachimi et al. 2011).

An intrinsic alignment (IA) of sources has an impact on the ratio statistics b_2D (θ_ap) and r_2D (θ_ap), Eqs. (14) and (15), mainly through $\langle M_{\text{ap}}^2\rangle ({\theta }_{\text{ap}})$ if we separate sources and lenses in redshift. The impact can be mitigated by using an appropriate model for ${\langle M_{\text{ap}}^2\rangle }_{\text{th}}({\theta }_{\text{ap}})$ and ${\langle \mathcal{N}{M}_{\text{ap}}\rangle }_{\text{th}}({\theta }_{\text{ap}})$ in the normalisation of the measurements. For this study, we do not include II or GI correlations in our synthetic mock data but, instead, predict the amplitude of potential systematic errors when ignoring the intrinsic alignment for future applications in Sect. 7.3.

For a reasonable prediction of the GI and II contributions to ${\langle M_{\text{ap}}^2\rangle }_{\text{th}}({\theta }_{\text{ap}})$, we use the recent non-evolution model utilised in Hildebrandt et al. (2017). This model is implemented by using $$P_{\kappa }^{\prime }(\mathcal{l})=P_{\kappa }(\mathcal{l})+P_{\kappa }^{\text{II}}(\mathcal{l})+P_{\kappa }^{\text{GI}}(\mathcal{l})$$(28)

instead of (24) in Eq. (20). The new II and GI terms are given by $$\begin{array}{l}{P}_{\kappa }^{\text{II}}(\mathcal{l})={\displaystyle {\int }_0^{{\chi }_{\text{h}}}\frac{\text{d}\chi p_{\text{s}}^2(\chi )}{f_K^2(\chi )}}{F}_{\text{ia}}^2(\chi )P_{\text{m}}(k_{\mathcal{l}}^{\chi };\chi );\\ P_{\kappa }^{\text{GI}}(\mathcal{l})=\frac{3H_0^2{\Omega }_{\text{m}}}{c^2}{\displaystyle {\int }_0^{{\chi }_{\text{h}}}\frac{\text{d}\chi p_{\text{s}}(\chi )g_{\text{s}}(\chi )}{a(\chi )f_K(\chi )}}{F}_{\text{ia}}(\chi )P_{\text{m}}(k_{\mathcal{l}}^{\chi };\chi ),\end{array}$$(29) (30)

where $$\begin{array}{ll}{F}_{\text{ia}}(\chi )\hfill & :=-A_{\text{ia}}{\text{C}}_1{\overline{\rho }}_{\text{crit}}\frac{{\Omega }_{\text{m}}}{D_{+}(\chi )}\hfill \\ \hfill & \approx 2.4\times {10}^{-2}\left(\frac{A_{\text{ia}}}{3.0}\right)\left(\frac{{\Omega }_{\text{m}}}{0.3}\right){\left(\frac{D_{+}(\chi )}{0.5}\right)}^{-1}\hfill \end{array}$$(31)

controls the correlation amplitude in the so-called “non-linear linear” model (see Hirata & Seljak 2004; Bridle & King 2007; or Joachimi et al. 2011; for details). The factor A_ia scales the amplitude; it broadly falls within A_ia ∈ [−3, 3] for recent cosmic-shear surveys and is consistent with A_ia ≈ 2 for sources in the Kilo-Degree Survey (Joudaki et al. 2018; Hildebrandt et al. 2017; Heymans et al. 2013). For the normalisation of F_ia (χ), we use ${\text{C}}_1=5\times {10}^{-14}{h}^{-2}{M}_{\odot }^{-1}{\text{Mpc}}^3$, and the linear structure-growth factor D₊(χ), normalised to unity for χ = 0 (Peebles 1980). By comparing $P_{\kappa }^{\text{II}}(\mathcal{l})$ and P_n(ℓ) in Eq. (22) we see that II contributions are essentially the clustering of sources on the sky (times a small factor). Likewise, $P_{\kappa }^{\text{GI}}(\mathcal{l})$ is essentially the correlation between source positions and their shear on the sky (cf. Eq. (23)). In this IA model, we assume a scale-independent galaxy bias for sources in the IA modelling since F_ia (χ) does not depend on k.

In Fig. 2, we plot the predicted levels of II and GI terms in the observed $\langle M_{\text{ap}}^2\rangle$ for varying values of A_ia as black and red lines for our MS cosmology and the p_s(z) in our mock survey. The corresponding value of A_ia is shown asa number in the figure key. We use negative values of A_ia for GI to produce positive correlations for the plot; the corresponding predictions for − A_ia have the same amplitude as those for A_ia but with opposite sign. The II terms, on the other hand, are invariant with respect to a sign flip of A_ia. All curves in the plot use a matter power spectrum P_m(k;χ) computed with Halofit (Smith et al. 2003) and the update in Takahashi et al. (2012). For comparison, we plot as blue line GG the theoretical $\langle M_{\text{ap}}^2\rangle$ without GI and II terms. For |A_ia|≈ 3, GI terms can reach levels up to 10% to 20% of the correlation signal for θ_ap ≳ 1′, whereas II terms are typically below 10%. The GI and II terms partly cancel each other for A_ia > 0 so that the contamination is worse for negative A_ia.

Moreover, we quantify the GI term in $\langle \mathcal{N}{M}_{\text{ap}}\rangle$ by using in (19) the modified power spectrum $$P_{\text{n}\kappa }^{\prime }(\mathcal{l};b,r)=P_{\text{n}\kappa }(\mathcal{l};b,r)+P_{\text{n}\kappa }^{\text{GI}}(\mathcal{l};b,r)$$(32)

with $$P_{\text{n}\kappa }^{\text{GI}}(\mathcal{l};b,r)={\displaystyle {\int }_0^{{\chi }_{\text{h}}}\frac{\text{d}\chi p_{\text{s}}(\chi )p_{\text{d}}(\chi )}{f_K^2(\chi )}}b(k_{\mathcal{l}}^{\chi })r(k_{\mathcal{l}}^{\chi })F_{\text{ia}}(\chi )P_{\text{m}}(k_{\mathcal{l}}^{\chi };\chi ).$$(33)

This is the model in Joudaki et al. (2018; see their Eq. (11)) with an additional term $r(k_{\mathcal{l}}^{\chi })$ that accounts for a decorrelation of the lens galaxies. This GI model is essentially the relative clustering between lenses and unbiased sources on the sky and therefore vanishes in the absence of an overlap between the lens and source distributions, which means that ∫ dχ p_s(χ) p_d(χ) = 0. In Fig. 3, we quantify the relative change in $\langle \mathcal{N}{M}_{\text{ap}}\rangle$ owing to the GI term for different values of A_ia. Since the change is very similar for all galaxy samples in the same redshift bin, we plot only the results for SM4. The overlap between sources and lenses is only around 4% for low-z samples and, therefore, the change stays within 2% for all angular scales considered here (SES13). On the other hand, for high-z samples where we have roughly 14% overlap between the distributions, the change can amount to almost 10% for A_ia ≈±2 and could have a significant impact on the normalisation.

Fig. 2 Levels of GI and II contributions to $\langle M_{\text{ap}}^2\rangle$ for different values of Aia (red and black lines labelled “II ± Aia” and “GI ± Aia”). The line “GG” is the theoretical $\langle M_{\text{ap}}^2\rangle$ without GI and II terms; the data points are measurements on the mocks for sources with shear and shape noise (MS γ+n), reduced shear and shape noise (MS g+n), and shear without shape noise (MS γ). The error bars indicate jackknife errors inflated by a factor of five for clarity (Appendix B).

Fig. 3 Relative change of $\langle \mathcal{N}{M}_{\text{ap}}\rangle$ for present GI correlations with different amplitudes Aia as indicated by “GI ± Aia”. The figure uses SM4 as fiducial lens sample; the results for other samples are similar. The thin lines within ± 2% are for the low-z sample, and the thick lines are for the high-z sample.

Fig. 4 Relative errors in the aperture statistics due to magnification bias of the lenses. Left: errors for $\langle \mathcal{N}{M}_{\text{ap}}\rangle$ where different line styles distinguish the galaxy samples. Larger errors for the same sample correspond to the high-z bin, smaller errors to low-z. Right: percentage errors for $\langle {\mathcal{N}}^2\rangle$ where larger errors for the same line style are the high-z bias.

3.4 Higher-order corrections

Corrections to the (first-order) Born approximation or for the magnification of the lenses cannot always be neglected as done in Eq. (23) (e.g. Ziour & Hui 2008; Hilbert et al. 2009; Hartlap 2009). This uncorrected equation over-predicts the power spectrum P_nκ(ℓ) by up to 10% depending on the galaxy selection and the mean redshift of the lens sample; the effect is smaller in a flux-limited survey but also more elaborate to predict as it depends on the luminosity function of the lenses. Hilbert et al. (2009) tests this for the tangential shear around lenses by comparing (23) to the full-ray-tracing results in the MS data, which account for contributions from lens-lens couplings and the magnification of the angular number density of lenses.

For a volume-limited lens sample, Hartlap (2009, H09 hereafter), derives the second-order correction (in our notation) $$P_{\text{n}\kappa }^{(2)}(\mathcal{l})=-\frac{9H_0^4{\Omega }_{\text{m}}}{2c^4}{\displaystyle {\int }_0^{{\chi }_{\text{h}}}\text{d}}\chi \frac{g_{\text{s}}(\chi )g_{\text{d}}(\chi )}{a^2(\chi )}{P}_{\text{m}}(k_{\mathcal{l}}^{\chi };\chi ),$$(34)

where $$g_{\text{d}}(\chi )={\displaystyle {\int }_{\chi }^{{\chi }_{\text{h}}}\text{d}}{\chi }^{\prime }{p}_{\text{d}}({\chi }^{\prime })\frac{f_K({\chi }^{\prime }-\chi )}{f_K({\chi }^{\prime })},$$(35)

for a more accurate power spectrum $P_{\text{n}\kappa }^{\prime }(\mathcal{l})=P_{\text{n}\kappa }(\mathcal{l})+P_{\text{n}\kappa }^{(2)}(\mathcal{l})$ that correctly describes the correlations in the MS. Physically, this correction accounts for the magnification of the projected number density of lens galaxies by matter in the foreground. We find that the thereby corrected ${\langle \mathcal{N}{M}_{\text{ap}}\rangle }_{\text{th}}({\theta }_{\text{ap}};b,r)$ can be different to the uncorrected aperture statistic by up to a few % (see left-hand panel in Fig. 4). This directly affects the normalisation of r_2D (θ_ap): the measured, normalised correlation r_2D(θ_ap) would be systematically low. We obtain Fig. 4 by comparing the uncorrected to the corrected ${\langle \mathcal{N}{M}_{\text{ap}}\rangle }_{\text{th}}$ for each of our lens-galaxy samples. In accordance with H09, we find that the systematic error is not negligible for some lens sample, and we therefore include this correction by employing $P_{\text{n}\kappa }^{\prime }(\mathcal{l})$ instead of P_{n κ} (ℓ) in the normalisation f_r (θ_ap) and in the prediction r_2D (θ_ap;b, r). This improves the accuracy of the lensing reconstruction of r(k) by up to a few %, most notably the sample blue high-z, especially around k ≈ 1 h Mpc⁻¹, which corresponds to θ_ap ≈ 10′.

Additional second-order terms for P_nκ(ℓ) arise due to a flux limit of the survey (Eqs. (3.129) and (3.130) in H09), but they require a detailed model of the luminosity function for the lenses. We ignore these contributions here because our mock lens samples, selected in redshift bins of Δ z ≈ 0.2 and for stellar masses greater than 5 × 10⁹ M_⊙, are approximately volume limited because of the lower limit of stellar masses and the redshift binning (see Sect. 4.1 in Simon et al. 2017, which uses our lens samples).

Similarly, by $$\begin{array}{l}{P}_{\text{n}}^{(2)}(\mathcal{l};b,r)=\frac{9H_0^2{\Omega }_{\text{m}}^2}{c^2}{\displaystyle {\int }_0^{{\chi }_{\text{h}}}\text{d}}\chi \frac{g_{\text{d}}^2(\chi )}{a^2(\chi )}{P}_{\text{m}}(k_{\mathcal{l}}^{\chi };\chi )-\frac{6H_0^2{\Omega }_{\text{m}}}{c^2}\\ \times {\displaystyle {\int }_0^{{\chi }_{\text{h}}}\text{d}}\chi \frac{p_{\text{d}}(\chi )g_{\text{d}}(\chi )}{a(\chi )f_K(\chi )}b(k_{\mathcal{l}}^{\chi })r(k_{\mathcal{l}}^{\chi })P_{\text{m}}(k_{\mathcal{l}}^{\chi };\chi )\end{array}$$(36)

H09 gives a second-order correction for P_n(ℓ;b) in addition to more corrections for flux-limited surveys (Eqs. (3.140)–(3.143)). We include $P_{\text{n}}^{(2)}(\mathcal{l};b,r)$ by using $P_{\text{n}}^{\prime }(\mathcal{l};b,r)=P_{\text{n}}(\mathcal{l};b)+P_{\text{n}}^{(2)}(\mathcal{l};b,r)$ instead of Eq. (22) in the following for f_b(θ_ap) and b_2D(θ_ap;b, r), although this correction is typically below half a % here (see the right-hand panel in Fig. 4).

4 Model templates of biasing functions

Apart from the galaxy-bias normalisation, the ratio statistics b_2D and r_2D are model-free observables of the spatial biasing functions, averaged for the radial distribution of lenses. The de-projection of the ratio statistics into (an average) b(k) and r(k) is not straightforward due to the radial and transverse smoothing in the projection. Therefore, for a de-projection we construct a parametric family of templates that we forward-fit to the ratio statistics. In principle, this family could be any generic function but we find that physical templates that can be extrapolated to scales unconstrained by the observations result in a more stable de-projection. To this end, we pick a template prescription that is motivated by the halo-model approach but with more freedom than is commonly devised (Cooray & Sheth 2002, for a review). Notably, we derive explicit expressions for b(k) and r(k) in a halo-model framework.

4.1 Separation of small and large scales

Before we outline the details of our version of a halo model, used to construct model templates, we point out that any halo model splits the power spectra P_m(k), P_gm(k), and P_g (k) into one- and two halo terms, $$P(k)=P^{1\text{h}}(k)+P^{2\text{h}}(k).$$(37)

The one-halo term P^1h(k) dominates at small scales, quantifying the correlations between density fluctuations within the same halo, whereas the two-halo term P^2h (k) dominates the power spectrum at large scales where correlations between fluctuations in different halos and the clustering of halos become dominant.

We exploit this split to distinguish between galaxy bias on small scales (one-halo terms) and galaxy bias on large scales (two-halo terms), namely $$b^{1\text{h}}(k):=\sqrt{\frac{P_{\text{g}}^{1\text{h}}(k)}{P_{\text{m}}^{1\text{h}}(k)}};b^{2\text{h}}(k):=\sqrt{\frac{P_{\text{g}}^{2\text{h}}(k)}{P_{\text{m}}^{2\text{h}}(k)}}$$(38)

and $$r^{1\text{h}}(k):=\frac{P_{\text{gm}}^{1\text{h}}(k)}{\sqrt{P_{\text{g}}^{1\text{h}}(k)P_{\text{m}}^{1\text{h}}(k)}};r^{2\text{h}}(k):=\frac{P_{\text{gm}}^{2\text{h}}(k)}{\sqrt{P_{\text{g}}^{2\text{h}}(k)P_{\text{m}}^{2\text{h}}(k)}},$$(39)

and we derive approximations for both regimes separately. We find that the two-halo biasing functions are essentially constants, and the one-halo biasing functions are only determined by the relation between matter and galaxy density inside halos.

To patch together both approximations of the biasing functions in the one-halo and two-halo regime, we then do the following. Based on Eq. (10), the function b² (k) is a weighted mean of b^1h(k) and b^2h(k): $$\begin{array}{lll}{b}^2(k)\hfill & =\hfill & \frac{P_{\text{g}}^{1\text{h}}(k)+P_{\text{g}}^{2\text{h}}(k)}{P_{\text{m}}(k)}\hfill \\ \hfill & =\hfill & \frac{P_{\text{m}}^{1\text{h}}(k){[b^{1\text{h}}(k)]}^2}{P_{\text{m}}(k)}+\frac{P_{\text{m}}^{2\text{h}}(k){[b^{2\text{h}}(k)]}^2}{P_{\text{m}}(k)}\hfill \\ \hfill & =\hfill & (1-W_{\text{m}}(k)){[b^{1\text{h}}(k)]}^2+W_{\text{m}}(k){[b^{2\text{h}}(k)]}^2,\hfill \end{array}$$(40)

where the weight $$W_{\text{m}}(k):=\frac{P_{\text{m}}^{2\text{h}}(k)}{P_{\text{m}}(k)}$$(41)

is the amplitude of the two-halo matter power spectrum relative to the total matter power spectrum. Deep in the one-halo regime we have W_m(k) ≈ 0 but W_m(k) ≈ 1 in the two-halo regime. Since the two-halo biasing is approximately constant, the scale dependence of galaxy bias is mainly a result of the galaxy physics inside halos and the shape of W_m(k).

Once the weight W_m(k) is determined for a fiducial cosmology, it does not rely on galaxy physics and we can use it for any model of b^1h (k) and b^2h (k). In principle, the weight W_m(k) could be accurately measured from a cosmological simulation by correlating only the matter density from different halos for $P_{\text{m}}^{2\text{h}}(k)$, which is then normalised by the full power spectrum P_m(k) in the simulation. We, however, determine W_m(k) by computing the one-halo and two-halo term of P_m(k) with the setup of Simon et al. (2009). Our results for W_m(k) at different redshifts are plotted in Fig. 5. There we find that the transition between the one-halo and two-halo regime, W_m ~ 0.5, is at k ~ 0.3 h Mpc⁻¹ for z = 0, whereas the transition point moves to k ~ 1 h Mpc⁻¹ for z ~ 1.

Similar to b(k), we can expand the correlation function r(k) in terms of its one-halo and two-halo biasing functions. To this end, let $$W_{\text{g}}(k):=\frac{P_{\text{g}}^{2\text{h}}(k)}{P_{\text{g}}(k)}={\left(\frac{b^{2\text{h}}(k)}{b(k)}\right)}^2{W}_{\text{m}}(k)$$(42)

be a weight by analogy to W_m(k). For unbiased galaxies, that is b^2h(k) = b(k) = 1, we simply have W_g(k) = W_m(k). Using the definition of r(k) in Eq. (10) and Eq. (39), we generally find $$\begin{array}{ll}r(k)\hfill & =\sqrt{(1-W_{\text{m}}(k))(1-W_{\text{g}}(k))}{r}^{1\text{h}}(k)\hfill \\ \hfill & +\sqrt{W_{\text{m}}(k)W_{\text{g}}(k)}{r}^{2\text{h}}(k).\hfill \end{array}$$(43)

Fig. 5 Weight Wm(k) of the two-halo term in the matter-power spectrum for varying redshifts z.

4.2 Halo-model definitions

For approximations of the biasing functions in the one- and two-halo regime, we apply the formalism in Seljak (2000) and briefly summarise it here. All halo-related quantities depend on redshift. In the fits with the model later on, we use for this the mean redshift of the lens galaxies.

We shall denote by n(m) dm the (comoving) number density of halos within the halo-mass range [m, m + dm); ⟨N|m⟩ is the mean number of galaxies inside a halo of mass m; ⟨N(N − 1)|m⟩ is the mean number of galaxy pairs inside a halo of mass m. Let u(r, m) be the radial profile of the matter density inside a halo or the galaxy density profile. Also let $$\tilde{u}(k,m)=\frac{{\displaystyle {\int }_0^{\infty }\text{d}}ss{k}^{-1}u(s,m)\sin (ks)}{{\displaystyle {\int }_0^{\infty }\text{d}}s{s}^{2}u(s,m)}$$(44)

be its normalised Fourier transform. Owing to this normalisation, profiles obey ũ (k, m) = 1 at k = 0. To assert a well-defined normalisation of halos, we truncate them at their virial radius r_vir, which we define by the over density ${\Omega }_{\text{m}}{\overline{\rho }}_{\text{crit}}{\Delta }_{\text{vir}}(z)$ within the distance r_vir from the halo centre and by Δ_vir(z) as in Bullock et al. (2001). Furthermore, the mean matter and galaxy number density (comoving) are $${\overline{\rho }}_{\text{m}}={\displaystyle \int \text{d}}mn(m)m; {\overline{n}}_{\text{g}}={\displaystyle \int \text{d}}mn(m)\langle N|m\rangle .$$(45)

The one-halo terms of the galaxy power spectrum P_g(k), the matter power spectrum P_m(k), and the galaxy-matter cross-power spectrum P_gm(k) are $$\begin{array}{lll}{P}_{\text{g}}^{1\text{h}}(k)\hfill & =\hfill & {\displaystyle {\int }_0^{\infty }\frac{\text{d}mn(m)}{{\overline{n}}_{\text{g}}^2}}{\tilde{u}}_{\text{g}}^{2p}(k,m)\langle N(N-1)|m\rangle ;\hfill \\ P_{\text{m}}^{1\text{h}}(k)\hfill & =\hfill & {\displaystyle {\int }_0^{\infty }\frac{\text{d}mn(m)m^2}{{\overline{\rho }}_{\text{m}}^2}}{\tilde{u}}_{\text{m}}^2(k,m);\hfill \\ P_{\text{gm}}^{1\text{h}}(k)\hfill & =\hfill & {\displaystyle {\int }_0^{\infty }\frac{\text{d}mn(m)m}{{\overline{\rho }}_{\text{m}}{\overline{n}}_{\text{g}}}}{\tilde{u}}_{\text{m}}(k,m){\tilde{u}}_{\text{g}}^q(k,m)\langle N|m\rangle .\hfill \end{array}$$(46) (47) (48)

In these equations, the exponents p and q are modifiers of the statistics for central galaxies, which are accounted for in the following simplistic way: central galaxies are by definition at the halo centre r = 0; one galaxy inside a halo is always a central galaxy; their impact on galaxy power spectra is assumed to be only significant for halos that contain few galaxies. Depending on whether a halo contains few galaxies or not, the factors (p, q) switch on or off a statistics dominated by central galaxies through $$\begin{array}{lll}p\hfill & =\hfill & \{\begin{array}{ll}1,\hfill & \text{\hspace{0.17em}for} \langle N(N-1)|m\rangle >1\text{}\hfill \\ 1/2,\hfill & \text{\hspace{0.17em}otherwise}\hfill \end{array};\hfill \\ q\hfill & =\hfill & \{\begin{array}{ll}1,\hfill & \text{\hspace{1em}\hspace{0.17em}for} \langle N|m\rangle >1\text{}\hfill \\ 0,\hfill & \text{\hspace{1em}\hspace{0.17em}otherwise}\hfill \end{array}.\hfill \end{array}$$(49)

We note that p and q are functions of the halo mass m. Later in Sect. 4.6, we consider also more general models where there can be a fraction of halos that contains only satellite galaxies. We achieve this by mixing (46)–(48) with power spectra in a pure-satellite scenario, this means a scenario where always p ≡ q ≡ 1.

We now turn to the two-halo terms in this halo model. We approximate the clustering power of centres of halos with mass m by $b_{\text{h}}^2(m)P_{\text{lin}}(k)$, where P_lin(k) denotes the linear matter power spectrum, and b_h(m) is the halo bias factor on linear scales; the clustering of halos is thus linear and deterministic in this description. Likewise, this model approximates the cross-correlation power-spectrum of halos with the masses m₁ and m₂ by b_h (m₁) b_h(m₂) P_lin(k). The resulting two-halo terms are then $$\begin{array}{lll}{P}_{\text{g}}^{2\text{h}}(k)\hfill & =\hfill & \frac{P_{\text{lin}}(k)}{{\overline{n}}_{\text{g}}^2}{\left({\displaystyle {\int }_0^{\infty }\text{d}}mn(m)\langle N|m\rangle b_{\text{h}}(m){\tilde{u}}_{\text{g}}(k,m)\right)}^2\text{};\hfill \\ P_{\text{m}}^{2\text{h}}(k)\hfill & =\hfill & \frac{P_{\text{lin}}(k)}{{\overline{\rho }}_{\text{m}}^2}{\left({\displaystyle {\int }_0^{\infty }\text{d}}mn(m)m{b}_{\text{h}}(m){\tilde{u}}_{\text{m}}(k,m)\right)}^2\text{};\hfill \\ P_{\text{gm}}^{2\text{h}}(k)\hfill & =\hfill & \frac{P_{\text{lin}}(k)}{{\overline{n}}_{\text{g}}{\overline{\rho }}_{\text{m}}}{\displaystyle {\int }_0^{\infty }\text{d}}mn(m)\langle N|m\rangle b_{\text{h}}(m){\tilde{u}}_{\text{g}}(k,m)\hfill \\ \hfill & \hfill & \times {\displaystyle {\int }_0^{\infty }\text{d}}mn(m)m{b}_{\text{h}}(m){\tilde{u}}_{\text{m}}(k,m).\hfill \end{array}$$(50) (51) (52)

The two-halo terms ignore power from central galaxies because it is negligible in the two-halo regime.

4.3 A toy model for the small-scale galaxy bias

We first consider an insightful toy model of b(k) and r(k) at small scales. In this model, both the matter and the galaxy distribution shall be completely dominated by halos of mass m₀, such that we find an effective halo-mass function n(m) ∝ δ_D(m − m₀); its normalisation is irrelevant for the galaxy bias. In addition, the halos of the toy model shall not cluster so that the two-halo terms of the power spectra vanish entirely. The toy model has practical relevance in what follows later because the one-halo biasing functions that we derive afterwards are weighted averages of toy models with different m₀. For this reason, most of the features can already be understood here, albeit not all, and it already elucidates biasing functions on small scales.

Let us define the variance ${\sigma }_N^2(m_0)=\langle N^2|m_0\rangle -{\langle N|m_0\rangle }^2$ of the halo-occupation distribution (HOD) in excess of a Poisson variance ⟨N|m₀⟩ by $$\Delta {\sigma }_N^2(m_0)={\sigma }_N^2(m_0)-\langle N|m_0\rangle .$$(53)

If the model galaxies obey Poisson statistics, they have $\Delta {\sigma }_N^2(m_0)=0$. We can nowwrite the mean number of galaxy pairs as $$\begin{array}{ll}\langle N(N-1)|m_0\rangle \hfill & =\langle N^2|m_0\rangle -\langle N|m_0\rangle \hfill \\ \hfill & ={\langle N|m_0\rangle }^2\left(1+\frac{\Delta {\sigma }_N^2(m_0)}{{\langle N|m_0\rangle }^2}\right).\hfill \end{array}$$(54)

By using Eqs. (46)–(48) with n(m) ∝ δ_D(m − m₀), the correlation factor reads $$\begin{array}{ll}R(k,m_0)\hfill & =\frac{{\tilde{u}}_{\text{g}}^{q-p}(k,m_0)\langle N|m_0\rangle }{\sqrt{\langle N(N-1)|m_0\rangle }}\hfill \\ \hfill & ={\tilde{u}}_{\text{g}}^{q-p}(k,m_0){\left(1+\frac{\Delta {\sigma }_N^2(m_0)}{{\langle N|m_0\rangle }^2}\right)}^{-1/2}\text{},\hfill \end{array}$$(55)

and the bias factor is $$\begin{array}{ll}B(k,m_0)\hfill & =\frac{{\tilde{u}}_{\text{g}}^p(k,m_0)\sqrt{\langle N(N-1)|m_0\rangle }}{{\tilde{u}}_{\text{m}}(k,m_0)\langle N|m_0\rangle }\hfill \\ \hfill & =\frac{{\tilde{u}}_{\text{g}}^q(k,m_0)}{{\tilde{u}}_{\text{m}}(k,m_0)}\frac{1}{R(k,m_0)}.\hfill \end{array}$$(56)

To avoid ambiguities in the following, we use capital letters for the biasing functions in the toy model.

We dub galaxies “faithful tracers” of the matter density if they have both (i) ũ_g (k, m) = ũ_m(k, m) and (ii) no central galaxies (p = q = 1). Halos with relatively small numbers of galaxies, that is ⟨N|m₀⟩, ⟨N(N − 1)|m₀⟩≲ 1, are called “low-occupancy halos” in the following. This toy model then illustrates the following points.

• Owing to galaxy discreteness, faithful tracers are biased if they not obey Poisson statistics. Namely, for a sub-Poisson variance, $\Delta {\sigma }_N^2(m_0)<0$, they produce opposite trends R(k, m₀) > 1 and B(k, m₀) < 1 with k, and vice versa for a super-Poisson sampling, but generally we find the relation R(k, m₀) × B(k, m₀) = 1.

• Nevertheless faithful tracers obey B(k, m₀), R(k, m₀) ≈ 1 if the excess variance becomes negligible, that is if $\Delta {\sigma }_N^2(m_0)\ll {\langle N|m_0\rangle }^2$. The discreteness of galaxies therefore becomes only relevant in low-occupancy halos.

• A value of R(k, m₀) > 1 occurs once central galaxies are present (p, q < 1). As a central galaxy is always placed at the centre, central galaxies produce a non-Poisson sampling of the profile u_m(r, m₀). In contrast to faithful galaxies with a non-Poisson HOD, we then find agreeing trends with scale k for R(k, m₀) and B(k, m₀) if $\Delta {\sigma }_{\text{N}}^2(m_0)=0$. Again, this effect is strong only in low-occupancy halos.

• The biasing functions in the toy model are only scale-dependent if galaxies are not faithful tracers. The bias function B(k, m₀) varies with k if either ũ_m(k, m₀)≠ũ_g(k, m₀) or for central galaxies (p≠1). The correlation function R(k, m₀) is scale-dependent only for central galaxies, that is p − q≠0, which then obeys $R(k,m_0)\propto {\tilde{u}}_{\text{g}}^{-1/2}(k,m_0)$. Variations with k become small for both functions, however, if ũ_m(k, m₀), ũ_g(k, m₀) ≈ 1, which is on scales larger than the size r_vir of a halo.

We stress again that a counter-intuitive r(k) > 1 is a result of the definition of P_g(k) relative to Poisson shot-noise and the actual presence of non-Poisson galaxy noise. One may wonder here if r > 1 is also allowed for biasing parameters defined in terms of spatial correlations rather than the power spectra. That this is indeed the case is shown in Appendix A for completeness.

4.4 Galaxy biasing at small scales

Compared to the foregoing toy model, no single halo mass scale dominates the galaxy bias at any wave number k for realistic galaxies. Nevertheless, we can express the realistic biasing functions b^1h (k) and r^1h (k) in the one-halo regime as weighted averages of the toy model B(k, m) and R(k, m) with modifications.

To this end, we introduce by $$b(m)=\frac{\langle N|m\rangle }{m}\frac{{\overline{\rho }}_{\text{m}}}{{\overline{n}}_{\text{g}}}$$(57)

the “mean biasing function”, which is the mean number of halo galaxies ⟨N|m⟩ per halo mass m in units of the cosmic average ${\overline{n}}_{\text{g}}/{\overline{\rho }}_{\text{m}}$ (Cacciato et al. 2012). If galaxy numbers linearly scale with halo mass, that means ⟨N|m⟩∝ m, we find a mean biasing function of b(m) = 1, while halos masses devoid of galaxies have b(m) = 0. For convenience, we make useof ⟨N|m⟩∝ m b(m) instead of ⟨N|m⟩ in the following equations because we typically find ⟨N|m⟩∝ m^β with β ≈ 1: b(m) is therefore usually not too different from unity.

Using Eqs. (46) and (47) we then find $${[b^{1\text{h}}(k)]}^2=\frac{P_{\text{g}}^{1\text{h}}(k)}{P_{\text{m}}^{1\text{h}}(k)}={\displaystyle {\int }_0^{\infty }\text{d}}m{w}_{20}^{1\text{h}}(k,m)b^2(m)B^2(k,m),$$(58)

with $w_{20}^{1\text{h}}(k,m)$ being one case in a family of (one-halo) weights, $$w_{ij}^{1\text{h}}(k,m):=\frac{n(m)m^2{\tilde{u}}_{\text{m}}^i(k,m){[b(m){\tilde{u}}_{\text{g}}(k,m)]}^j}{{\displaystyle {\int }_0^{\infty }\text{}}\text{d}mn(m)m^2{\tilde{u}}_{\text{m}}^i(k,m){[b(m){\tilde{u}}_{\text{g}}(k,m)]}^j}.$$(59)

This family and the following weights w(k, m) are normalised, which means that ∫ dm w(k, m) = 1. The introduction of these weight functions underlines that the biasing functions are essentially weighted averages across the halo mass spectrum as, for example, ${[b^{1\text{h}}(k)]}^2$, which is the weighted average of b²(m) B²(k, m).

The effect of $w_{20}^{1\text{h}}(k,m)$ is to down-weight large halo masses in the bias function because $w_{20}^{1\text{h}}(k,m_1)/w_{20}^{1\text{h}}(k,m_2)\propto {\tilde{u}}_{\text{m}}^2(k,m_1)/{\tilde{u}}_{\text{m}}^2(k,m_2)$ decreases with m₁ for a fixed m₂ < m₁ and k. Additionally, the relative weight of a halo with mass m decreases towards larger k because ũ_m(k, m) tends to decrease with k. As a result, at a given scale k only halos below a typical mass essentially contribute to the biasing functions (Seljak 2000).

We move on to the correlation factor r^1h(k) in the one-halo regime. Using Eqs. (46)–(48) and the relations $$\begin{array}{l}\langle N(N-1)|m\rangle =R^{-2}(k,m){\tilde{u}}_{\text{g}}^{2q-2p}(k,m)\langle N|m\rangle ;\\ {\tilde{\rho }}_{\text{m}}(k,m)=m{\tilde{u}}_{\text{m}}(k,m);\\ {\tilde{\rho }}_{\text{g}}(k,m)=\langle N|m\rangle {\tilde{u}}_{\text{g}}(k,m),\end{array}$$(60) (61) (62)

we write $$r^{1\text{h}}(k)=\frac{P_{\text{gm}}^{1\text{h}}(k)}{\sqrt{P_{\text{g}}^{1\text{h}}(k)P_{\text{m}}^{1\text{h}}(k)}} =: {\zeta }_{\text{sat}}(k){\zeta }_{\text{cen}}(k){\zeta }_{\Delta \sigma }(k)$$(63)

as product of the three separate factors $$\begin{array}{lll}{\zeta }_{\text{sat}}(k)\hfill & :=\hfill & \frac{{\displaystyle {\int }_0^{\infty }\text{d}}mn(m){\tilde{\rho }}_{\text{m}}(k,m){\tilde{\rho }}_{\text{g}}(k,m)}{{\left({\displaystyle {\int }_0^{\infty }\text{d}}mn(m){\tilde{\rho }}_{\text{g}}^2(k,m){\displaystyle {\int }_0^{\infty }\text{d}}mn(m){\tilde{\rho }}_{\text{m}}^2(k,m)\right)}^{1/2}};\hfill \\ {\zeta }_{\text{cen}}(k)\hfill & :=\hfill & {\displaystyle {\int }_0^{\infty }\text{d}}m{w}_{11}^{1\text{h}}(k,m){\tilde{u}}_{\text{g}}^{q-1}(k,m);\hfill \\ {\zeta }_{\Delta \sigma }(k)\hfill & :=\hfill & {\left({\displaystyle {\int }_0^{\infty }\text{d}}m{w}_{02}^{1\text{h}}(k,m){\tilde{u}}_{\text{g}}^{2q-2}(k,m)R^{-2}(k,m)\right)}^{-1/2}\hfill \end{array}$$(64) (65) (66)

with the following meaning:

• The first factor ζ_sat(k) quantifies, at spatial scale k, the correlation between the radial profiles of the matter density ${\tilde{\rho }}_{\text{m}}(k,m)$ and the (average) number density of satellite galaxies ${\tilde{\rho }}_{\text{g}}(k,m)$ across the halo mass spectrum n(m). As upperbound we always have |ζ_sat(k)|≤ 1 because of the Cauchy-Schwarz inequality when applied to the nominator of Eq. (64). Thus ζ_sat(k) probably reflects best what we intuitively understand by a correlation factor between galaxies and matter densities inside a halo. Since it only involves the average satellite profile, the satellite shot-noise owing to a HOD variance is irrelevant at this point. The next two factors can be seen as corrections to ζ_sat(k) owing tocentral galaxies or a non-Poisson HOD variance.

• The second factor ζ_cen(k) is only relevant in the sense of ζ_cen(k)≠1 through low-occupancy halos with central galaxies (q≠1). It has the lower limit ζ_cen(k) ≥ 1 because of ũ_g(k, m) ≤ 1 and hence ${\tilde{u}}_{\text{g}}^{q-1}(k,m)\ge 1$. This correction factor can therefore at most increase the correlation r^1h(k).

• The third factor ζ_Δσ(k) is the only one that is sensitive to an excess variance $\Delta {\sigma }_{\text{N}}^2(m)\ne 0$ of the HOD, namely through R(k, m). In the absence of central galaxies, that means for p ≡ q ≡ 1, ${\zeta }_{\Delta \sigma }^2(k)$ is the (weighted) harmonic mean of R²(k, m), or the harmonic mean of the reduced ${[{\tilde{u}}_{\text{g}}(k,m)R(k,m)]}^2\le R^2(k,m)$ otherwise.