© ESO 2021

1. Introduction

The broad diversity of galaxy types and their large-scale clustering requires complex and numerically challenging physics that in models is approximated or empirically calibrated. This leaves us today with a variety of modelling techniques that produce a range of theoretical, not always consistent predictions for galaxy statistics (e.g. Vogelsberger et al. 2020; Springel et al. 2018; Lacey et al. 2016; Henriques et al. 2017; Campbell et al. 2018; Saghiha et al. 2017). Of major interest for model tests, due to the key role of dark matter halos for galaxy physics, is the predicted relative clustering of matter and galaxies: the galaxy bias.

Galaxy bias varies with galaxy type (Mo et al. 2010, and references therein), rendering it a type-dependent footprint for the galaxy-matter connection. Motivated by this idea, as a test we compare measurements of these footprints in the Canada-France-Hawaii Lensing Survey (CFHTLenS) data (Erben et al. 2013) to predictions by two semi-analytic models (SAMs) for various galaxy types, redshifts, and spatial scales from the linear down to the highly non-linear regime at k ≈ 10 h Mpc⁻¹. They are the Munich model ‘H15’ by Henriques et al. (2015) and the Durham model ‘L12’ by Lagos et al. (2012). These SAMs are different implementations of a technique that populates dark matter halos in cosmological simulations, here the Millennium Simulation (Springel et al. 2005), with model galaxies by emulating known galaxy physics (Cole et al. 2000, for a general overview). Their predictions for the excess surface matter-density around galaxy pairs were already studied in Saghiha et al. (2017) and Simon et al. (2019). Here we focus on galaxy biasing in these SAMs.

Among the many ways to define galaxy bias (see e.g. Desjacques et al. 2018), our model tests employ the definition by Tegmark & Bromley (1999) in terms of galaxy and matter power spectra for the (comoving) Fourier modes k,

$$\begin{array}{c}\hfill b(k):=\sqrt{\frac{P_{\mathrm{g}}(k)}{P_{\mathrm{m}}(k)}};r(k):=\frac{P_{\mathrm{gm}}(k)}{\sqrt{P_{\mathrm{m}}(k)P_{\mathrm{g}}(k)}}.\end{array}$$(1)

These two biasing functions are measured for a specific galaxy population and redshift. They incorporate the power spectra of fluctuations in the matter density field, P_m(k); in the galaxy density field, P_g(k); and their cross power P_gm(k). In the definition of P_g(k), a Poisson sampling shot-noise is subtracted, which may produce r(k) > 1 (Peebles 1980; Guzik & Seljak 2001). The Eq. (1) are sensitive to a stochastic bias; however, restricted to two-point statistics, they do not discriminate between a non-linear deterministic bias and a stochastic bias. The only exceptions are linear scales where k ≪ 1 h Mpc⁻¹ and density fluctuations $\stackrel{\sim }{\delta }(\mathit{k})$ are small (Somerville et al. 2001; Dekel & Lahav 1999). For example, on linear scales r(k) < 1 indicates a stochastic scatter ϵ_k in the relation ${\stackrel{\sim }{\delta }}_{\mathrm{g}}(\mathit{k})=b(k){\stackrel{\sim }{\delta }}_{\mathrm{m}}(\mathit{k})+{ϵ}_k$ between galaxy density contrast, ${\stackrel{\sim }{\delta }}_{\mathrm{g}}$, and matter density contrast, ${\stackrel{\sim }{\delta }}_{\mathrm{m}}$, whereas for large density fluctuations r(k) < 1 can also indicate ${\stackrel{\sim }{\delta }}_{\mathrm{g}}(\mathit{k})=b_k[{\stackrel{\sim }{\delta }}_{\mathrm{m}}(\mathit{k})]+{ϵ}_k$ with a non-linear mapping b_k[x] and smaller (or vanishing) ϵ_k. Studies of galaxy clustering show that both non-linearity and stochasticity increase towards smaller scales, and galaxy bias is essentially linear and deterministic (ϵ_k = 0) at 2πk⁻¹ ≳ 30 h⁻¹ Mpc (e.g. Friedrich et al. 2018; Gruen et al. 2018; Dvornik et al. 2018; Di Porto et al. 2016; Patej & Eisenstein 2016; Kovač et al. 2011; Swanson et al. 2008; Wild et al. 2005; Tegmark & Bromley 1999). Despite the limitations, these measures of galaxy biasing highlight differences between galaxy populations and their relation to the matter density field, they are easily measured in simulations (Springel et al. 2018; Weinberg et al. 2004), and are hence a useful quantity to test galaxy models.

Perfectly suited to measuring scale-dependent galaxy bias is the analysis of the weak gravitational-lensing effect (Schneider et al. 2006). This exploits small coherent distortions of distant galaxy images (sources) by the large-scale matter-density field in the foreground. For second-order statistics the analysis correlates shear distortions with lens galaxy positions to obtain, apart from an amplitude normalisation, galaxy bias in projection. This can be normalised or deprojected if the redshift distributions of sources and lens galaxies, as well as the background cosmology, mainly Ω_m, is known (van Waerbeke 1998; Schneider 1998; Hoekstra et al. 2001). Several previous studies have successfully applied this (or a similar) lensing analysis (Dvornik et al. 2018; Prat et al. 2018; Comparat et al. 2013; Jullo et al. 2012; Simon et al. 2007; Pen et al. 2003; Hoekstra et al. 2001), but they neglect the effects of lensing magnification and source intrinsic alignments (e.g. Hirata & Seljak 2004; Unruh et al. 2020). In particular, intrinsic alignments (IA) have to be included in the lensing signal to reach an accuracy of several per cent for b(k) and r(k), as demonstrated in Simon & Hilbert (2018, hereafter SH18) by using mock galaxy and shear catalogues that are similar to ours here. We apply the SH18 procedure with some improvements to measure the functions in Eq. (1).

This lensing technique probes galaxy biasing with a minimum of assumptions and is also valid in the non-linear regime. This sets it apart from other techniques that, for instance, only apply in the weakly linear regime; require a theoretical reference model for the clustering of matter as a function of scale; or depend on a detailed physical model for the clustering of matter and galaxies (e.g. Cacciato et al. 2013; Marulli et al. 2013; Dvornik et al. 2018; Heymans et al. 2021). However, our technique employs flexible model functions for b(k) and r(k) to stabilise the deprojection of the projected galaxy bias. Although based on a halo model prescription, their purpose notably is not an accurate interpretation of the projected galaxy bias in terms of a physical model, but rather an accurate deprojection of the model-free biasing functions. The completely generic templates in SH18 could be used equally well, albeit with less accuracy. To highlight our focus on the biasing functions, we refer to the stabilising model functions in the reconstruction as ‘templates’.

The outline of this paper is as follows. Section 2 describes our CFHTLenS shear catalogue and the selection of galaxy (foreground) samples for the analysis. For these we estimate the redshift distributions (lenses and sources) and the mean number density of the lenses; both are used in the analysis. For the model predictions, we select similar galaxy samples from the SAMs in Sect. 3, accounting for the stellar-mass errors in CFHTLenS. Section 4 summarises the template-based deprojection technique for b(k) and r(k). Our new analysis improves on that of SH18 with regard to the transformation bias in the aperture statistics (Sect. 4.4), the prior for the template parameters of the biasing functions (Appendix E), and an allowed mismatch between matter and satellite galaxy distribution inside dark matter halos (Sect. 4.6). We present raw (i.e. not normalised) CFHTLenS measurements for the projected galaxy bias already in this section. We compare our main results for b(k) and r(k) to the H15 and L12 predictions in Sect. 5, along with some physical interpretation, and discuss the results in Sect. 6.

The fiducial background cosmology defines the projection of galaxy bias and the lensing kernel; a small additional dependence is introduced by the shape of a fiducial matter power spectrum that impacts the smoothing in the projection. Our results adopt a flat Λ cold dark matter (ΛCDM) cosmology with a matter density parameter of Ω_m = 0.306, a Hubble constant of H₀ = h 100 km s⁻¹ Mpc⁻¹ and h = 0.683, a baryon matter density of Ω_b = 0.0479, and a primordial spectral index of n_s = 0.965, in agreement with Planck Collaboration XIII (2016) (Sect. 4.1). Systematic errors in the amplitude of the CFHTLenS biasing functions due to uncertainties in the fiducial cosmological model and galaxy redshift distributions are explored in Sect. 4.9.

2. Canada-France-Hawaii Lensing Survey

The CFHTLenS¹ is a multi-colour lensing survey based on the Wide part of the Canada-France-Hawaii Survey (CFHTLS) taken with the 2k × 4k MegaPrime camera (MegaCam) between 2003 and 2008 (Erben et al. 2013). The pixel size of the camera CCD chips is 0.187 arcsec. The survey consists of four patches W1 (∼ 63.8 deg²), W2 (∼ 22.6 deg²), W3 (∼ 44.2 deg²), and W4 (∼ 23.3 deg²) with 154 deg² total sky area. The patches are mosaics of 171 1 × 1 deg² MegaCam fields, observed with the five filters u^*g′r′i′z′ to the limiting AB magnitudes of z′ ≲ 23.5 mag and between 24.5 and 25.5 mag for the other filters (5σ detection limit in a 2 arcsec aperture). All filters were obtained under sub-arcsec seeing conditions of which i′, used for shear estimates, has a median seeing of 0.68 arcsec. After a break of the i′ filter in 2008, 33 fields were observed with the slightly different y′ filter. The survey combines weak-lensing data processing with THELI (Erben et al. 2013), shear measurement with lensfit (Miller et al. 2013), and photometric redshift measurement with PSF-matched photometry (Hildebrandt et al. 2012). A full systematic-error analysis of the shear measurements in combination with the photometric redshifts is presented in Heymans et al. (2012), with additional error analyses of the photometric redshift measurements presented in Benjamin et al. (2013).

The survey data product is the object and lensing catalogue v1.1², for which galaxies were selected from the i′-band stacks with the SExtractor software (Bertin & Arnouts 1996). The revised photometry of these objects yields a root mean square (RMS) uncertainty of σ ≈ 0.01 − 0.03 mag for the photometric calibration in all five bands (Hildebrandt et al. 2012). From these data, photometric redshifts z_p as galaxy distance estimators were inferred with the Bayesian Photometric Redshift Code (BPZ; Benítez 2000).

From the catalogue we used only the 129 science-ready fields that are considered suitable for cosmic shear studies in Heymans et al. (2012). To estimate the effective unmasked area of this subset of fields, we counted the number of pixels in the mask files of the co-added images with mask = 0; only source and lens galaxies from these regions are selected for the analysis. Table 1 lists the effective survey area of all patches W1 to W4, amounting in total to an effective area of 92.12 deg². We summarise here the details of our lens and source samples.

2.1. Galaxy samples

We measure galaxy bias for foreground galaxies that were binned in stellar mass (SM1 to SM6) or colour (RED or BLUE), and additionally in redshift (low-z or high-z). Table 2 is an overview of all 2 × 8 = 16 samples and their selection criteria.

All samples select only bright objects (17.5 ≤ i′< 22.5) within regions of mask = 0. At this brightness the error on the redshift estimators is typically σ_z = (1 + z) 0.04 with a few per cent outliers. We used stellar-mass estimates in the CFHTLenS catalogue. These estimators apply stellar-synthesis modelling to the photometric data in the five bands, assuming a Chabrier (2003) initial stellar mass-function, the stellar population synthesis model by Bruzual & Charlot (2003), the dust extinction law in Calzetti et al. (2000), and a flat ΛCDM model with Hubble parameter of h = 0.7 and Ω_m = 0.27 for the luminosity distance, which is similar to our fiducial cosmology in Table 4. For the stellar masses, we assume the RMS error of 0.27 dex in Velander et al. (2014).

As in SH18 we split the galaxy samples into six stellar mass bins between 5 × 10⁹ and $3.2\times {10}^{11}{h}_{70}^{-2}{\mathrm{M}}_{\odot }$, and into two colour bins along the line u^* − r′ = 1.93 z_p + 1.85 for the full stellar mass range. We further divided these eight samples in photometric redshift to obtain low-z samples with 0.2 ≤ z_p <  0.44 ($\overline{z}\approx 0.35$) and high-z samples with 0.44 ≤ z_p <  0.6 ($\overline{z}\approx 0.51$). For each sample, the total number of objects is in the range 10³ − 10⁵, which corresponds to a mean number density of 0.003 − 0.3 arcmin⁻² on the sky.

2.2. Shear catalogues

As background sources for the lensing analysis, we use resolved galaxies with reliable photo-z estimates and i′< 24.7. Each source i has an estimated ϵ_i = ϵ_i1 + i ϵ_i2 gravitational shear γ_i in the source line-of-sight direction and a statistical weight w_i >  0, which depends on the variance of its posterior ellipticity and the ellipticity distribution in the entire source galaxy population. The additional calibration factor m_i >  −1, determined with simulated CFHTLenS images, addresses the multiplicative shear bias. We accounted for the shear calibration in our analysis, equivalent to the procedure described in Miller et al. (2013), by rescaling the source ellipticities ϵ_i → ϵ_i (1 + m_i)⁻¹ and weights w_i → w_i (1 + m_i). Furthermore, we subtracted from ϵ_i2 the additive shear bias c₂, provided by the catalogue.

We then constructed six source samples according to Table 3. For the full analysis, sources are in two samples, BACK low-z and BACK high-z, which have little redshift overlap with the low-z and high-z foreground samples in Sect. 2.1, but aim at a large number density of sources on the sky. For low-z we applied the cut 0.65 ≤ z_p <  1.3, and for high-z the cut 0.75 ≤ z_p <  1.3. Additionally, for robustness tests of our results later on, we subdivided each source bin into the two complementary photo-z bins BACKa and BACKb. In total, there are roughly 2 × 10⁶ sources in each sample BACK low-z and high-z, which corresponds to a number density of 6 arcmin⁻². However, due to the source weights, the effective source number $N_{\text{eff}}:={({\sum }_i{w}_i)}^2/{\sum }_i{w}_i^2$ is roughly 20% smaller. The mean, w-weighted redshift of BACK low-z (high-z) is $\overline{z}=0.97(1.03)$, and the weighted ellipticity dispersion (i.e. the source shape noise) in all samples is ⟨|ϵ|²⟩^1/2 ≈ 0.41, or 0.29 for each component individually.

2.3. Redshift distributions

The average BPZ redshift posterior of the sample galaxies is the estimated redshift distribution of a sample (which is the sum of all posteriors followed by a normalisation). According to Benjamin et al. (2013), these estimates provide a sufficiently accurate model for the actual redshift distribution p(z) of CFHTLenS galaxies at z <  1.3.

However, the BPZ posterior is binned with bin width Δz = 0.05 and therefore artificially broadened. While this has a small effect for cosmic shear or galaxy-galaxy lensing studies, which are mostly sensitive to the mean of the distribution, it is potentially problematic for predictions of the angular clustering of galaxies and for the galaxy-bias normalisation. As correction of the broadening, we deconvolved the binned, average posterior for the lens samples as described in Appendix B. For the source samples, where the broadening is not relevant, we directly used the BPZ posteriors.

The two top panels in Fig. 1 show the deconvolved p(z) for all our lens samples. Despite the same photo-z cuts, the differences in the redshift distributions between the lens samples is evident: bluer samples or samples with less stellar mass tend to be more spread out in redshift. The bottom panel plots the average BPZ posteriors (weighted) of all source samples.

Fig. 1. Redshift distributions p(z) of all galaxy samples. Two top panels: lens samples, low-z and high-z, bottom panel: two source samples (solid lines) and their two subsamples (dotted or dashed lines).

2.4. Average number densities

We obtain the mean (comoving) number density of lens galaxies with our modified V_max-estimator, described in Appendix C. Originally, V_max assumed an exactly determined maximum (comoving) volume V_i inside which a galaxy i is observable under the adopted selection criteria; the estimator of the number density was then ${\overline{n}}_{\mathrm{g}}=\sum _{i=1}^n{V}_i^{-1}$ for all n galaxies in the sample (Takeuchi et al. 2000; Schmidt 1968). In our case with photometric redshifts, however, the volume V_i is uncertain. Our modified V_max estimator therefore accounts for this uncertainty, and we find values for ${\overline{n}}_{\mathrm{g}}$ that are then a few per cent larger than the standard V_max based on the best-fit photo-z. Our estimates safely ignore a density correction for inhomogeneous galaxy samples because our broad redshift bins at z >  0.1 in combination with the ∼ 100 deg² survey area span a sufficiently large and homogeneous volume (e.g. Baldry et al. 2012).

Using the V_max densities for CFHTLenS in Fig. 2 (black data points with error bars) as priors improves the convergence of our Monte Carlo analysis because they break parameter degeneracies in the galaxy-bias templates. The priors are log-normal distributions with RMS errors for ${log}_{10}{\overline{n}}_{\mathrm{g}}$ based on jackknife re-sampling (Knight 1999); we removed one MegaCam field at a time for each jackknife sample. The number densities typically range between 2 × 10⁻⁴ − 4 × 10⁻³ h³ Mpc⁻³ and have a clear trend towards fewer galaxies for higher stellar masses.

Fig. 2. Comoving number density of galaxies for CFHTLenS (black data points) and the SAMs (green data points) in our stellar-mass bins for low-z (filled symbols) and high-z (open symbols). The sample labels SM1 to SM6 are indicated at the bottom; the small numbers at the top are the mean stellar masses inside the bins in units of $h_{70}^{-2}{\mathrm{M}}_{\odot }$. The coloured lines are the galaxy SMFs by Wright et al. (2018, GAMA/G10-COSMOS/3D-HST), Bielby et al. (2012, WIRDS/CFHTLS), Davidzon et al. (2013, VIPERS/CFHTLS), and Baldry et al. (2012, GAMA). The SMFs are, before integrating over our stellar-mass bins, convolved with a stellar-mass error of 0.27 dex to account for the Eddington bias in our data.

Figure 2 also compares our estimates for ${\overline{n}}_{\mathrm{g}}$ to measured stellar-mass functions (SMFs) from selected studies, which are shown as lines. For this comparison, we convolved the best-fitting SMFs in Wright et al. (2018), Bielby et al. (2012), Davidzon et al. (2013), and Baldry et al. (2012) with our RMS error of 0.27 dex for stellar-mass estimates (log-normal error model), and we then integrated the convolved functions over the stellar-mass ranges in Table 2. The convolution accounts for the Eddington bias, which softens the drop in galaxy counts beyond the characteristic mass in the Schechter function around $5\times {10}^{10}{h}_{70}^{-2}{\mathrm{M}}_{\odot }$ and which suppresses the counts below this cutoff. The top axis of the figure indicates, averaged for low-z and high-z, the mean stellar masses of our samples. All SMFs in this comparison make assumptions for the stellar population synthesis and the stellar initial mass-function that are consistent with the CFHTLenS stellar masses, and they employ, with the exception of Wright et al. (2018) and partly Davidzon et al. (2013), the same dust extinction law.

That there are still differences between the surveys is probably a sign of inconsistent stellar-mass estimates or galaxy selections. Closest to our selection is the VIPERS study by Davidzon et al. (2013), which identifies galaxies in the patches W1 and W4 of CFHTLS (at that time ∼ 15 deg²) with our flux limit i_AB <  22.5. This sample was, superior to our photometric redshift estimates, followed up with the VIMOS multi-object spectrograph at the ESO VLT to obtain precise estimates of stellar mass and redshift with HYPERZMASS (Bolzonella et al. 2010, 2000). In addition to the i_AB-band flux limit, they applied a colour cut to exclude low-redshift objects (z <  0.5). The resulting SMF in their lowest redshift bin z ∈ [0.5, 0.6] agrees well with our high-z galaxy counts SM5 and SM6, but is ∼40% lower for the less massive galaxies SM1 to SM3. In another study, Bielby et al. (2012) combine deep CFHTLS (∼ 2 deg²) photometry with near-infrared observations by the WIRcam Deep Survey (WIRDS) to attain the galaxy SMF at various redshifts, in particular also around our low-z and high-z. Again, there is consistency for SM5 and SM6, but our galaxy counts are lower (higher) for low-z (high-z) in the other stellar-mass bins. This discrepancy for less massive galaxies may be related to their near-infrared colour cuts, which were used to separate stars from galaxies, and the fainter flux-limit of i = 25.5. Finally, the studies by Wright et al. (2018) and Baldry et al. (2012) are based on the Galaxy And Mass Assembly (GAMA) survey and, in the former case, additionally on the deeper observations in the Cosmos Evolution Survey (COSMOS) and 3D-HST. They consistently find higher galaxy counts for low-z, whereas now high-z SM1 to SM4 are in good agreement with our results.

In conclusion, our measurements for ${log}_{10}{\overline{n}}_{\mathrm{g}}$ are broadly supported within 0.2 − 0.3 dex by other studies. Mismatches are presumably a mixture of inconsistent selections and stellar mass errors. A clear explanation probably requires a completely consistent selection of galaxies, beyond the scope of this study.

3. Semi-analytic models of galaxies

Semi-analytic galaxy models (SAMs) populate matter halos in purely dark matter simulations with galaxies, reducing galaxies and halos to entities with a few global properties that follow a simplified set of physical and empirical laws. By applying prescriptions for a stellar population synthesis and dust extinction the SAM results are, with increased uncertainty, converted into galaxy spectral energy-distributions, and thus mock galaxy catalogues with multi-band photometry. A detailed description of the semi-analytic approach can be found in Lacey et al. (2016), among others. We use the Munich model H15 and the Durham model L12.

3.1. Simulation and galaxy sample selection

Our galaxy mocks are based on the Millennium Simulation (Springel et al. 2005). This N-body simulation of the dark matter density field has a comoving spatial resolution of ∼ 5 h⁻¹ kpc over a cubic volume with 500 h⁻¹ Mpc side length, sampled by ∼10¹⁰ mass particles. The fiducial cosmology of the simulation (Ω_m = 1 − Ω_Λ = 0.25, Ω_b = 0.045, σ₈ = 0.9, h = 0.73, n_s = 1.0) is different from our fiducial Planck cosmology in Table 4. As shown in Henriques et al. (2017), however, this has only a mild impact on the galaxy clustering because the combination of a higher σ₈ and a lower Ω_m gives a similar structure growth. We therefore do not expect large differences between the galaxy models and the CFHTLenS galaxy bias from this side.

To compute predictions of the galaxy biasing functions, we combine the galaxy and dark matter power spectra, and their cross-power spectra at the Millennium Simulation snapshots z = 0.3623,  0.5086 for low-z and high-z, exactly as described in Sect. 2.3 of SH18. To prepare the various lens samples in these snapshots, we selected SAM galaxies according to Table 2, after adding a random RMS error of 0.27 dex to the SAM stellar masses log₁₀M_* and after converting SDSS magnitudes of the SAMs to CFHTLenS magnitudes according to Sect. 3.2 in Saghiha et al. (2017). By SAM stellar mass we specifically mean that we do not apply BPZ to the SAM photometry to emulate CFHTLenS, rather we utilise the randomised SAM stellar masses. For the colour cut, we identify z_p with the snapshot z.

3.2. Model differences

The H15 and L12 models follow somewhat differing SAM approaches and hence produce distinct predictions for the galaxy bias. A major difference is their calibration: H15 fits stellar-mass functions and the red-galaxy fractions as a function of M_* to observations for z ≤ 1; L12 calibrates the model with local luminosity functions (near-infrared and blue bands) and the observed correlation between galactic bulge mass and central black hole mass. In addition, the identification of subhalos in the Millennium Simulation differs and so does the treatment of subhalos and satellite galaxies. L12 satellites, for instance, can merge with a central galaxy while they still have resolved subhalos, whereas H15 satellites have to disrupt their halos first. This directly affects the galaxy clustering on small scales and the radial distribution of satellites inside galaxy clusters (Contreras et al. 2013). Another difference, affecting the photometry of the model galaxies, are the initial mass functions, the implemented models for the stellar population synthesis, and the dust extinction law.

More can be seen in the stellar mass function (SMF) of the models. The green data points in Fig. 2 show the number density of the SAM galaxy samples for the various stellar masses (filled data points: low-z; open data points: high-z). We obtained these data by dividing the total number of galaxies at the snapshot redshifts by the comoving simulation volume of (500 h⁻¹ Mpc)³. Within 0.1 dex these data points closely follow those of CFHTLenS with two exceptions. Firstly, the L12 high-z samples (open triangles) have up to 0.5 dex fewer galaxies for SM1 and SM2, due to a strong evolution of the stellar-mass function at low stellar masses. This strong evolution is not present in H15. Secondly, the SAMs predict approximately 0.2 dex more SM6 galaxies for both redshifts.

For the colour-selected samples we find values broadly compatible with CFHTLenS, but also differences between H15 and L12: for RED low-z in H15 (L12), we have ${log}_{10}({\overline{n}}_{\mathrm{g}}/h^3{\mathrm{Mpc}}^{-3})=-2.20(-2.11)$ and for RED high-z −2.37(−2.37). The corresponding figures for BLUE low-z and high-z are −2.02(−2.37) and −1.99(−2.51), respectively. L12 thus leans towards a smaller number of blue (SM1 and SM2) galaxies, especially at high-z.

4. Aperture statistics and reconstruction of biasing functions

This section briefly summarises the SH18 formalism for galaxy-bias reconstruction and reports the CFHTLenS aperture statistics to which we apply this technique. The basic idea is to measure the fluctuations of galaxy numbers and lensing convergence inside circular apertures for a range of angular scales. The ratios of their (co-)variances directly probe the projected scale-dependent galaxy bias on the sky. By then forward-fitting templates to the observed ratio statistics we reconstruct the average biasing functions b(k) and r(k) for a given fiducial cosmology and galaxy redshift distributions.

4.1. Fiducial cosmology

The baseline cosmology in the analysis are the parameters listed in Table 4, taken from the recent cosmic microwave background (CMB) analysis by the Planck Collaboration XIII (2016). In this baseline, we changed the normalisation σ₈ of the matter power spectrum to the smaller value of σ₈ = 0.714 in order to be consistent with cosmic-shear constraints σ₈(Ω_m/0.27)^0.5 ≈ 0.76 in CFHTLenS (Joudaki et al. 2017; Fu et al. 2014; Heymans et al. 2013). For the following analysis of the ratio aperture-statistics, however, this smaller value has no practical relevance. This is reflected by the insensitivity of the ratio statistics with respect to any rescaling of the matter power spectrum amplitude.

Nevertheless, the ratio statistics depends weakly on the shape of the matter power spectrum, if it is not a power law in k (van Waerbeke 1998). In our analysis, the power spectrum P_m(k, χ) at comoving distance χ is the updated Halofit (Takahashi et al. 2012; Smith et al. 2003), as it is implemented in the latest version of Nicea³ (Kilbinger et al. 2009). In addition, we modify the power spectrum P_m(k, χ) ↦ F(k, χ) P_m(k, χ) in the non-linear regime by the transfer function F(k, χ) in Chisari et al. (2018) to correct for AGN feedback, star formation, and adiabatic contraction of matter halos due to gas cooling (Schneider & Teyssier 2015). Since F(k, χ) is not well known, we randomly vary F(k, z)−1 later on to include systematic errors by baryons.

4.2. Aperture statistics

This section is a summary of the relation between power spectra and aperture statistics. For a circular aperture of size θ_ap at position θ, the density fluctuation of galaxy numbers inside the aperture is

$$\begin{array}{c}\hfill \mathcal{N}({\theta }_{\mathrm{ap}};\mathit{\theta })={\displaystyle \int }{\mathrm{d}}^2{\theta }^{\prime }U(|{\mathit{\theta }}^{\prime }|;{\theta }_{\mathrm{ap}}){\kappa }_{\mathrm{g}}({\mathit{\theta }}^{\prime }+\mathit{\theta }).\end{array}$$(2)

This aperture number count convolves the galaxy number-density contrast ${\kappa }_{\mathrm{g}}(\mathit{\theta })=N(\mathit{\theta })/{\overline{N}}_{\mathrm{g}}-1$ of galaxy density N(θ) and mean density ${\overline{N}}_{\mathrm{g}}$ on the sky with the compensated polynomial filter

$$\begin{array}{cc}& U(\theta ;{\theta }_{\mathrm{ap}})=\frac{1}{{\theta }_{\mathrm{ap}}^2}u(\theta {\theta }_{\mathrm{ap}}^{-1});\hfill \\ \hfill & u(x)=\frac{9}{\pi }(1-x^2)(\frac{1}{3}-x^2)\mathrm{H}(1-x),\hfill \end{array}$$(3)

where H(x) is the Heaviside step function. The discreteness of galaxies at positions θ_i (i.e. N(θ) = ∑_iδ_D(θ − θ_i) for Dirac delta functions δ_D(θ)) gives rise to shot-noise in 𝒩(θ_ap;θ), which is corrected for in the statistical moments of 𝒩 below. The aperture number count 𝒩 compares to the aperture mass (i.e. the fluctuation of the lensing convergence κ(θ) inside the same aperture):

$$\begin{array}{c}\hfill M_{\mathrm{ap}}({\theta }_{\mathrm{ap}};\mathit{\theta })={\displaystyle \int }{\mathrm{d}}^2{\theta }^{\prime }U(|{\mathit{\theta }}^{\prime }|;{\theta }_{\mathrm{ap}})\kappa ({\mathit{\theta }}^{\prime }+\mathit{\theta }).\end{array}$$(4)

Our statistical analysis utilises second-order moments of the joint distribution of 𝒩 and M_ap over all positions θ. These second-order moments are measures of angular power spectra P(ℓ),

$$\begin{array}{cc}& {⟨{\mathcal{N}}^2⟩}_{\mathrm{th}}({\theta }_{\mathrm{ap}};b)=2\pi {\displaystyle \underset{0}{\overset{\infty }{\int }}\mathrm{d}\ell \ell P_{\mathrm{n}}(\ell ;b){[I(\ell {\theta }_{\mathrm{ap}})]}^2,}\hfill \end{array}$$(5)

$$\begin{array}{cc}\hfill & {⟨\mathcal{N}{M}_{\mathrm{ap}}⟩}_{\mathrm{th}}({\theta }_{\mathrm{ap}};b,r)=2\pi {\displaystyle \underset{0}{\overset{\infty }{\int }}\mathrm{d}\ell \ell P_{\mathrm{n}\kappa }(\ell ;b,r){[I(\ell {\theta }_{\mathrm{ap}})]}^2,}\hfill \end{array}$$(6)

$$\begin{array}{cc}\hfill & {⟨M_{\mathrm{ap}}^2⟩}_{\mathrm{th}}({\theta }_{\mathrm{ap}})=2\pi {\displaystyle \underset{0}{\overset{\infty }{\int }}\mathrm{d}\ell \ell P_{\kappa }(\ell ){[I(\ell {\theta }_{\mathrm{ap}})]}^2,}\hfill \end{array}$$(7)

for the angular band-pass filter

$$\begin{array}{c}\hfill I(x)=\frac{12}{\pi }\frac{{\mathrm{J}}_4(x)}{x^2}.\end{array}$$(8)

We estimate these statistics from observations (see Sect. 4.4), and denote them just by ⟨𝒩²(θ_ap)⟩, ⟨𝒩M_ap(θ_ap)⟩, and $\langle M_{\text{ap}}^2({\theta }_{\text{ap}})\rangle$ without the subscript ‘th’ (for theory) and arguments ‘b’ or ‘b, r’.

To model the aperture statistics, we use Eqs. (5)–(7) and the following three power-spectra models: P_n(ℓ; b), P_nκ(ℓ; b, r), and P_κ(ℓ). In the models the sources have the probability distribution p_s(χ) dχ in comoving distance χ, and lens galaxies have the distribution p_d(χ) dχ. The specific biasing functions of the lens galaxies are given by b(k) and r(k) as indicated by the arguments ‘b’ or ‘b, r’. There is the angular galaxy power-spectrum

$$\begin{array}{c}\hfill P_n(\ell ;b)={\displaystyle \underset{0}{\overset{{\chi }_{\mathrm{h}}}{\int }}\frac{\mathrm{d}\chi p_{\mathrm{d}}^2(\chi )}{f_K^2(\chi )}{b}^2(k_{\ell }^{\chi })P_{\mathrm{m}}(k_{\ell }^{\chi };\chi ),}\end{array}$$(9)

the galaxy-convergence cross-power

$$\begin{array}{cc}& P_{n\kappa }(\ell ;b,r)=\frac{3H_0^2{\mathrm{\Omega }}_{\mathrm{m}}}{2c^2}\hfill \\ \hfill & {\displaystyle \underset{0}{\overset{{\chi }_{\mathrm{h}}}{\int }}\frac{\mathrm{d}\chi p_{\mathrm{d}}(\chi )g_{\mathrm{s}}(\chi )}{a(\chi )f_K(\chi )}b(k_{\ell }^{\chi })r(k_{\ell }^{\chi })P_{\mathrm{m}}(k_{\ell }^{\chi };\chi ),}\hfill \end{array}$$(10)

and the convergence power-spectrum

$$\begin{array}{c}\hfill P_{\kappa }(\ell )=\frac{9H_0^4{\mathrm{\Omega }}_{\mathrm{m}}^2}{4c^4}{\displaystyle \underset{0}{\overset{{\chi }_{\mathrm{h}}}{\int }}\frac{\mathrm{d}\chi g_{\mathrm{s}}^2(\chi )}{a^2(\chi )}{P}_{\mathrm{m}}(k_{\ell }^{\chi };\chi ),}\end{array}$$(11)

where f_K(χ) is the comoving angular diameter distance,

$$\begin{array}{c}\hfill g_{\mathrm{s}}(\chi ):={\displaystyle \underset{\chi }{\overset{{\chi }_{\mathrm{h}}}{\int }}\mathrm{d}{\chi }^{\prime }{p}_{\mathrm{s}}({\chi }^{\prime })\frac{f_K({\chi }^{\prime }-\chi )}{f_K({\chi }^{\prime })}}\end{array}$$(12)

is a lensing kernel, $k_{\mathcal{l}}^{\chi }:=(\mathcal{l}+0.5)/f_K(\chi )$ is a wave number at χ, a(χ) is the scale factor, χ_h is the maximum distance to a source, and J_n(x) is the nth-order Bessel function of the first kind.

4.3. Intrinsic alignment and magnification bias

We add corrections to the theoretical power spectra, Eqs. (5)–(7), to account for the intrinsic source alignments (IAs) and the magnification of lenses. To this end, we employ the approximate ‘non-linear linear’ model for IA (Hirata & Seljak 2004; Bridle & King 2007; Joachimi et al. 2011). Specifically, instead of P_κ(ℓ) in the integral Eq. (7), we integrate over

$$\begin{array}{c}\hfill P_{\kappa }^{\prime }(\ell )=P_{\kappa }(\ell )+P_{\kappa }^{\mathrm{II}}(\ell )+P_{\kappa }^{\mathrm{GI}}(\ell ),\end{array}$$(13)

where the additional intrinsic-intrinsic (II) and shear-intrinsic (GI) terms are⁴

$$\begin{array}{cc}\hfill P_{\kappa }^{\mathrm{II}}(\ell )& ={\displaystyle {\int }_0^{{\chi }_{\mathrm{h}}}\frac{\mathrm{d}\chi p_{\mathrm{s}}^2(\chi )}{f_K^2(\chi )}{F}_{\mathrm{ia}}^2(\chi )P_{\mathrm{m}}(k_{\ell }^{\chi };\chi );}\hfill \end{array}$$(14)

$$\begin{array}{cc}\hfill P_{\kappa }^{\mathrm{GI}}(\ell )& =\frac{3H_0^2{\mathrm{\Omega }}_{\mathrm{m}}}{c^2}{\int }_0^{{\chi }_{\mathrm{h}}}\frac{\mathrm{d}\chi p_{\mathrm{s}}(\chi )g_{\mathrm{s}}(\chi )}{a(\chi )f_K(\chi )}{F}_{\mathrm{ia}}(\chi )P_{\mathrm{m}}(k_{\ell }^{\chi };\chi ),\hfill \end{array}$$(15)

and

$$\begin{array}{c}\hfill F_{\mathrm{ia}}(\chi )\approx -2.4\times {10}^{-2}\left(\frac{A_{\mathrm{ia}}}{3.0}\right)\left(\frac{{\mathrm{\Omega }}_{\mathrm{m}}}{0.3}\right){\left(\frac{D_{+}(\chi )}{0.5}\right)}^{-1}\end{array}$$(16)

uses the linear growth factor of structure D₊(χ), normalised by D₊(χ = 0) (Peebles 1980). The parameter A_ia defines the amplitude of the intrinsic alignments in the model. In our Bayesian analysis we use as prior $A_{\mathrm{ia}}=-0.{48}_{-0.87}^{+0.75}$ the constraints from Heymans et al. (2013) (68% credible interval, hereafter CI).

Likewise, we correct for the magnification of the lens number density by foreground structure in the power spectrum P_nκ(ℓ; b, r), Eq. (6), through

$$\begin{array}{c}\hfill P_{n\kappa }^{\prime }(\ell ;b,r)=P_{n\kappa }(\ell ;b,r)+P_{n\kappa }^{\mathrm{GI}}(\ell ;b,r)+P_{n\kappa }^{\mathrm{mb}}(\ell )\end{array}$$(17)

with

$$\begin{array}{c}\hfill P_{n\kappa }^{\mathrm{GI}}(\ell ;b,r)={\displaystyle {\int }_0^{{\chi }_{\mathrm{h}}}\frac{\mathrm{d}\chi p_{\mathrm{s}}(\chi )p_{\mathrm{d}}(\chi )}{f_K^2(\chi )}b(k_{\ell }^{\chi })r(k_{\ell }^{\chi })F_{\mathrm{ia}}(\chi )P_{\mathrm{m}}(k_{\ell }^{\chi };\chi )}\end{array}$$(18)

and

$$\begin{array}{c}\hfill P_{n\kappa }^{\mathrm{mb}}(\ell )=-\frac{9H_0^4{\mathrm{\Omega }}_{\mathrm{m}}^2}{2c^4}{\displaystyle {\int }_0^{{\chi }_{\mathrm{h}}}\mathrm{d}\chi \frac{g_{\mathrm{s}}(\chi )g_{\mathrm{d}}(\chi )}{a^2(\chi )}{P}_{\mathrm{m}}(k_{\ell }^{\chi };\chi ),}\end{array}$$(19)

where

$$\begin{array}{c}\hfill g_{\mathrm{d}}(\chi ):={\displaystyle {\int }_{\chi }^{{\chi }_{\mathrm{h}}}\mathrm{d}{\chi }^{\prime }{p}_{\mathrm{d}}({\chi }^{\prime })\frac{f_K({\chi }^{\prime }-\chi )}{f_K({\chi }^{\prime })}}\end{array}$$(20)

is the lensing kernel of the foreground structure (Unruh et al. 2020; Hilbert et al. 2009; Hartlap 2009). The expression $P_{n\kappa }^{\mathrm{mb}}(\ell )$ assumes a volume-complete sample of lens galaxies which is approximately correct for our stellar-mass limited samples. The flux limit of i′< 22.5 mostly affects the samples SM1 and SM2 for low-z and high-z (Sect. 3.2 in Simon et al. 2013). The estimated magnification bias for ⟨𝒩M_ap⟩ can reach levels of several per cent for lenses with low clustering amplitude and high redshifts (z ∼ 0.5) and should therefore be modelled. The impact of magnification on ⟨𝒩²⟩, on other hand, is negligible but nevertheless included here (using Eq. (36) in SH18).

4.4. Estimators and data covariance

The second-order moments of the aperture statistics are linear transformations of the angular correlations functions for the lens clustering ω(θ), mean tangential shear around lenses ${\overline{\gamma }}_{\mathrm{t}}(\theta )$, and the cosmic shear ξ_±(θ), which all can be obtained relatively straightforwardly from observational data (e.g. Simon et al. 2007). To perform this transformation we measure with Athena⁵ the correlation functions for angular separations between 9 arcsec and 3 deg in 8000 bins in the four patches W1 to W4 and combine the measurements (Kilbinger et al. 2014). The result is statistics for n = 8 aperture radii θ_ap between 2 and 90 arcmin, equally spaced on a logarithmic scale.

The largest aperture size is our choice to retain a reasonably large number of jackknife samples; the resampling scheme removes subfields that stretch the full angular range of the correlation functions. The aim of the smallest aperture size, on the other hand, is to keep the bias in the estimated statistics below 10% at a few arcmin. The bias arises due to the blending of galaxy images which compromises lensfit below 9 arcsec, our smallest galaxy separation, and hence produces a gap in the correlation functions (Miller et al. 2013; Heymans et al. 2012). Previous studies set the correlation signal to zero inside the gap and produced a ‘transformation bias’ in the aperture statistics. We follow a revised, even though not perfect, approach with less transformation bias.

Specifically, we extrapolate the correlations below 9 arcsec in two different ways, depending on which correlation function is used. For ${\overline{\gamma }}_{\mathrm{t}}(\theta )$ and ξ_±(θ), which are respectively converted into ⟨𝒩M_ap⟩ and $\langle M_{\text{ap}}^2\rangle$, we extend the lowest angular bin to 1 arcsec, which means that we assume a constant signal equal to that at 9 arcsec. For ω(θ), which is converted into ⟨𝒩²⟩, we fit a power law ω(θ)∝θ^−δ to the data points between 9 arcsec ≤ θ ≤ 1 arcmin, and we then extrapolate the best fit to 0.1 arcsec ≤ θ ≤ 9 arcsec in the transformation. A power-law correlation function for galaxy clustering is a sensible approximation on small scales (e.g. Zehavi et al. 2005). To quantify the remaining transformation bias in the estimated aperture statistics after application of these extrapolations, we use mock galaxy (SM1-SM6, RED, and BLUE low-z and high-z) and shear catalogues that were constructed from the H15 data and previously used for the work in SH18. With these mocks we compare the estimated aperture statistics (with extrapolation) to the true aperture statistics (with no truncation at 9 arcsec). In summary, the transformation bias for the aperture moments ⟨𝒩M_ap⟩ and $\langle M_{\text{ap}}^2\rangle$ is typically below 1% at θ_ap ≈ 2 arcmin. For ⟨𝒩²⟩ it is about 10% at θ_ap ≈ 2 arcmin (biased low), but quickly falls below a few per cent beyond several arcmin.

For estimates of the errors in the aperture statistics, we perform jackknife resampling of the data (Singh et al. 2017; Norberg et al. 2009; Knight 1999). The resulting jackknife covariance matrix,

$$\begin{array}{c}\hfill \mathtt{C}=\frac{n_{\mathrm{jk}}-1}{n_{\mathrm{jk}}}{\displaystyle \sum _{i=1}^{n_{\mathrm{jk}}}({\mathit{d}}_i-\overline{\mathit{d}}){({\mathit{d}}_i-\overline{\mathit{d}})}^t;\overline{\mathit{d}}=\frac{1}{n_{\mathrm{jk}}}\sum _{i=1}^{n_{\mathrm{jk}}}{\mathit{d}}_i,}\end{array}$$(21)

for i = 1…n_jk jackknife samples,

$$\begin{array}{c}\hfill {\mathit{d}}_i=({\mathit{d}}_i^{(1)},{\mathit{d}}_i^{(2)},{\mathit{d}}_i^{(3)}),\end{array}$$(22)

comprises the error covariances between different angular scales and types of the aperture statistics, namely $d_{\mathit{ij}}^{(1)}={⟨{\mathcal{N}}^2⟩}_i({\theta }_{\mathrm{ap},j})$, $d_{\mathit{ij}}^{(2)}={⟨\mathcal{N}{M}_{\mathrm{ap}}⟩}_i({\theta }_{\mathrm{ap},j})$, and $d_{\mathit{ij}}^{(3)}={⟨M_{\mathrm{ap}}^2⟩}_i({\theta }_{\mathrm{ap},j})$. For each jackknife sample, we remove two adjacent MegaCam pointings from one field, redo the measurement, and finally combine the results from all four patches W1 to W4 into d_i. In three samples, where the fields have an odd number of pointings, we remove three pointings. In total, we thereby obtain n_jk = 63 jackknife samples to estimate the error covariance.

Figure 3 shows our measurements for $\langle M_{\text{ap}}^2\rangle$ and the jackknife error; the top panel is the E-mode signal, the bottom panel the B-mode signal, which is consistent with zero, relative to the E-mode. In our analysis the measurements for BACK low-z and BACK high-z are plotted as circles. The noisier data of the source-sample subdivisions a and b are depicted as crosses and stars. Importantly, the circle data points are the weighted average of crosses, stars, and the cross-correlation signal of the subdivisions (not shown here). The dashed and dotted lines are predictions based the cosmology in Table 4.

Fig. 3. Variance $\langle M_{\text{ap}}^2\rangle$ as function of aperture size for our source samples. Top panel: E-mode and the predictions (lines) for the BACK low-z and BACK high-z data points (solid and open circles). The crosses and stars correspond to the (noisier) signals of the four subdivisions a and b of the source samples. Bottom panel: B-mode relative to the E-mode signal (in percentage). This plot uses seven evenly spaced aperture radii instead of the eight used in the analysis.

4.5. Ratio statistics and projected galaxy bias

The ratio statistics

$$\begin{array}{cc}& {\widehat{b}}_{2\mathrm{D}}({\theta }_{\mathrm{ap}})=\sqrt{\frac{⟨{\mathcal{N}}^2⟩({\theta }_{\mathrm{ap}})}{⟨M_{\mathrm{ap}}^2⟩({\theta }_{\mathrm{ap}})}};\hfill \\ \hfill & {\widehat{r}}_{2\mathrm{D}}({\theta }_{\mathrm{ap}})=\frac{⟨\mathcal{N}{M}_{\mathrm{ap}}⟩({\theta }_{\mathrm{ap}})}{\sqrt{⟨{\mathcal{N}}^2⟩({\theta }_{\mathrm{ap}})⟨M_{\mathrm{ap}}^2⟩({\theta }_{\mathrm{ap}})}}\hfill \end{array}$$(23)

is the observed projection of the biasing functions b(k) and r(k) on the sky. These raw data are plotted in Fig. 4, as function of aperture size, for CFHTLenS with the bias factor ${\widehat{b}}_{2\mathrm{D}}$ in the top panels A and B, and the correlation coefficient ${\widehat{r}}_{2\mathrm{D}}$ in the bottom panels C and D. The error bars indicate symmetric jackknife errors, whereas the analysis below accounts for asymmetric errors in the ratio statistics through a more complex likelihood model. The likelihood model also accounts for the correlation between the data points at different aperture sizes, encoded in the jackknife covariance. Moreover, statistical errors for different lens samples are correlated because they share the same shear data.

Fig. 4. Observed ratio statistics with jackknife errors for the CFHTLenS lens samples in the low-z bin (panels A and C) and the high-z bin (panels B and D). The source samples are BACK low-z and BACK high-z, respectively. The data points are slightly shifted along the x-axis for better visibility. The dashed lines are predictions for unbiased galaxies (the galaxy-bias normalisation) with SM3 pd(z). The shaded regions bracket the prediction variations using the pd(z) of the other lens samples. Top panels: relative clustering ${\widehat{b}}_{2\mathrm{D}}({\theta }_{\mathrm{ap}})$ of galaxies and matter. Bottom panels: correlation ${\widehat{r}}_{2\mathrm{D}}({\theta }_{\mathrm{ap}})$ between the (projected) galaxy and matter distribution. This plot uses seven evenly spaced aperture sizes instead of the eight used in the analysis.

The lines and shaded areas inside the panels are predictions for unbiased galaxies and the fiducial cosmology in Table 4. They are traditionally used as normalisation for the projected galaxy bias to correct, among other things, for the distinct radial projection kernels of galaxies and matter (Hoekstra et al. 2002). Specifically, the dashed lines are the normalisation for SM3 low-z (left panels) and high-z (right panels). The normalisation varies weakly with aperture size, and therefore essentially rescales the amplitude of the inferred galaxy bias. This is also the case for the other samples. Across all lens samples in the same photo-z bin, however, the normalisation amplitude varies by about 10% within the shaded regions owing to differences among the p_d(z) (Fig. 1). The predictable uncertainties in the normalisation amplitude for a specific galaxy sample are estimators for the systematic error in the reconstructed biasing functions. This is exploited below.

Coming back to the CFHTLenS data points, our samples are usually more strongly clustered than unbiased galaxies (i.e. ${\widehat{b}}_{2\mathrm{D}}$ is larger). In addition, our galaxy samples are strongly correlated with the matter-density field, which means that ${\widehat{r}}_{2\mathrm{D}}$ is consistent with unbiased galaxies. Inside large apertures the correlation is similar for all lens samples, whereas there is some up-scatter below 20 arcmin. This hints at sample differences in the galaxy biasing at small scales. To investigate this in detail we deproject the aperture statistics.

4.6. Template models for deprojection

To reconstruct spatial biasing functions from the projected ratio statistics in Fig. 4 we fit the models

$$\begin{array}{cc}& {\widehat{b}}_{2\mathrm{D}}({\theta }_{\mathrm{ap}};b)=\sqrt{\frac{{⟨{\mathcal{N}}^2⟩}_{\mathrm{th}}({\theta }_{\mathrm{ap}};b)}{{⟨M_{\mathrm{ap}}^2⟩}_{\mathrm{th}}({\theta }_{\mathrm{ap}})}};\hfill \\ \hfill & {\widehat{r}}_{2\mathrm{D}}({\theta }_{\mathrm{ap}};b,r)=\frac{{⟨\mathcal{N}{M}_{\mathrm{ap}}⟩}_{\mathrm{th}}({\theta }_{\mathrm{ap}};b,r)}{\sqrt{{⟨{\mathcal{N}}^2⟩}_{\mathrm{th}}({\theta }_{\mathrm{ap}};b){⟨M_{\mathrm{ap}}^2⟩}_{\mathrm{th}}({\theta }_{\mathrm{ap}})}}\hfill \end{array}$$(24)

to the measurements by using templates for b(k) and r(k). Templates are useful here because, due to the smoothing in the projections (24), a deprojection into b(k) and r(k) for a given ${\widehat{b}}_{2\mathrm{D}}$ and ${\widehat{r}}_{2\mathrm{D}}$ without regularisation is unstable. We denote the templates b(k; Θ) and r(k; Θ) to indicate the set Θ of parameters they are based on.

We employ templates that are thoroughly tested with mock galaxy survey data, similar to our CFHTLenS data, in SH18. They are based on the halo model in Seljak (2000) but with more freedom for the halo-occupation distribution (HOD) of galaxies and a model-free galaxy bias on linear scales. This halo model ignores relevant physics, such as halo exclusion or assembly bias or galactic subhalos, and it has a crude treatment of central galaxies. Nevertheless, the templates produce an accurate deprojection of the biasing functions which is our primary interest here. We outline the template construction in Appendix D and restrict ourselves here to a brief description of the template parameters, listed in Table 5.

The templates have 49 free parameters Θ = (…,b_i, …, V_i, …, ζ, m_piv, b_ls, r_ls, f_cen). The large number of parameters is due to the fine sampling of the HOD mean number of galaxies, ⟨N|m⟩, and galaxy pairs, ⟨N(N − 1)|m⟩, with 22 interpolation points each that cover the broad halo mass range 10⁴ ≤ m ≤ 10¹⁶ h⁻¹ M_⊙. Our intention is to keep the HOD as flexible as possible and to avoid a model bias in the deprojection with a too restrictive HOD. The meaning of the parameters is as follows:

• The two parameters b_ls = b(k → 0) and r_ls = r(k → 0) quantify the linear stochastic bias on large (linear) scales. They are completely independent from the galaxy bias on small scales, i.e. independent from the template HOD;

• The parameter f_cen ∈ [0, 1] is the fraction of halos that permit a central galaxy in the selected galaxy population. This fraction is not to be confused with the fraction of central galaxies in a sample;

• A concentration parameter ζ ≠ 1 indicates biased satellite galaxies that have, averaged over all halos, a different distribution from the dark matter inside halos. For ζ >  1 satellites have a more diffuse distribution, and for ζ <  1 a more concentrated distribution;

• The pivotal mass m_piv defines the halo mass scale, where ⟨N|m_piv⟩ = 1 and N comprises both central and satellite galaxies;

• The mean number of galaxies in a halo of mass m, ⟨N|m⟩=b(m) m/(b(m_piv) m_piv), is expressed by the mean biasing function b(m) and its value at m_piv. The b(m) is linearly interpolated on a logarithmic m-scale with 22 interpolation points, and it has to satisfy the normalisation condition in Eq. (D.12).

• The mean number of galaxy pairs, expressed by the normalised excess variance V(m) via ⟨N(N − 1)|m⟩=⟨N|m⟩ [⟨N|m⟩+V(m)], is sensitive to the variance of galaxy numbers inside halos, which means for V(m) = 0 that galaxy numbers have a Poisson variance, for V(m) > 0 a super-Poisson variance, and for V(m) < 0 a sub-Poisson variance. Similar to b(m), V(m) is interpolated with 22 points. To be physically meaningful, the function V(m) has to be larger than max{−1, −⟨N|m⟩}, which is accounted for in the template fitting.

Unbiased galaxies have b_ls = r_ls = ζ = 1, f_cen = 0, b(m) = 1, and V(m) = 0 for all halo masses m.

The template parameter space is highly degenerate, which may pose numerical problems when sampling the space during our deprojection technique. Therefore, we also fit the relation

$$\begin{array}{c}\hfill {\overline{n}}_{\mathrm{g}}(\mathbf{\Theta })={\displaystyle {\int }_0^{\infty }\mathrm{d}mn(m)⟨N|m⟩=\frac{{\mathrm{\Omega }}_{\mathrm{m}}{\overline{\rho }}_{\mathrm{crit}}}{m_{\mathrm{piv}}b(m_{\mathrm{piv}})}}\end{array}$$(25)

to our observed ${\overline{n}}_{\mathrm{g}}$ in Fig. 2 to better confine m_piv and to break parameter degeneracies⁶. Here, ${\overline{\rho }}_{\mathrm{crit}}=3H_0^2/(8\pi G_{\mathrm{N}})$ denotes the critical matter density at z = 0. Equation (25) follows from the definition of b(m) and the condition Eq. (D.12).

4.7. Bayesian reconstruction of biasing functions

We forward-fit in a Bayesian analysis the templates b(k; Θ) and r(k; Θ) to the ratio statistics (Eq. (23)); the vector Θ comprises all template parameters in Table 5. For the likelihood function of the data we assume a Gaussian multivariate error distribution of the data,

$$\begin{array}{c}\hfill \mathit{d}=(\dots ,⟨\mathcal{N}⟩({\theta }_{\mathrm{ap},i}),\cdots ,⟨\mathcal{N}{M}_{\mathrm{ap}}⟩({\theta }_{\mathrm{ap},i}),\cdots ,⟨M_{\mathrm{ap}}^2⟩({\theta }_{\mathrm{ap},i}),\dots ),\end{array}$$(26)

and a lognormal error distribution of the galaxy number density ${\overline{n}}_{\mathrm{g}}$. Specifically, the probability distribution function p_δ of statistical errors δd and δ_n in d and ${log}_{10}{\overline{n}}_{\mathrm{g}}$, respectively, is

$$\begin{array}{c}\hfill -2ln{p}_{\delta }(\delta \mathit{d},{\delta }_n)=\mathrm{const}+n_{\mathrm{jk}}ln(1+\frac{\delta {\mathit{d}}^{\mathrm{T}}{\mathtt{C}}^{-1}\delta \mathit{d}}{n_{\mathrm{jk}}-1})+\frac{{\delta }_n^2}{{\sigma }_{\mathrm{logn}}^2},\end{array}$$(27)

for the error covariance C of d (Sect. 4.4) and an error variance σ_{log n} of ${log}_{10}{\overline{n}}_{\mathrm{g}}$ (error bars in Fig. 2). This noise model also marginalises over statistical errors in our estimate of C, which is constructed from n_jk jackknife samples (Sellentin & Heavens 2016). As discussed in SH18, the posterior density of the template parameters given the ratio statistics in Eq. (23) is then a marginalisation of p_δ with respect to $\mathit{x}:=(\dots ,\langle M_{\text{ap}}^2\rangle ({\theta }_{\text{ap},i}),\dots )$,

$$\begin{array}{cc}& p(\mathbf{\Theta }|\mathit{d},{\overline{n}}_{\mathrm{g}})\propto p_{\mathrm{\Theta }}(\mathbf{\Theta }){\displaystyle \int }\mathrm{d}\mathit{x}{p}_{\delta }(\mathit{d}-\mathit{m}(\mathbf{\Theta },\mathit{x}),\hfill \\ \hfill & {log}_{10}\frac{{\overline{n}}_{\mathrm{g}}}{{\overline{n}}_{\mathrm{g}}(\mathbf{\Theta })})p_x(\mathit{x}),\hfill \end{array}$$(28)

where m(Θ, x) = (m⁽¹⁾,m⁽²⁾,m⁽³⁾) and

$$\begin{array}{cc}& m_i^{(1)}={\widehat{b}}_{2\mathrm{D}}^2({\theta }_{\mathrm{ap},i};b)x_i;\hfill \end{array}$$(29)

$$\begin{array}{cc}\hfill & m_i^{(2)}={\widehat{b}}_{2\mathrm{D}}({\theta }_{\mathrm{ap},i};b){\widehat{r}}_{2\mathrm{D}}({\theta }_{\mathrm{ap},i};b,r)x_i;\hfill \end{array}$$(30)

$$\begin{array}{cc}\hfill & m_i^{(3)}=x_i.\hfill \end{array}$$(31)

The integral kernel p_δ(…) inside (28) is the likelihood of d given Θ, while p_Θ(Θ) and p_x(x) are prior densities for Θ and x_i. The prior densities are detailed in the following section.

Since we seek constraints on the biasing functions rather than Θ, we convert the posterior density of Θ into a posterior distribution for b(k) and r(k) by drawing ${\mathbf{\Theta }}_i\sim p(\mathbf{\Theta }|\mathit{d},{\overline{n}}_{\mathrm{g}})$ from (28) in a Monte Carlo process and by computing the mean and variance of {b(k; Θ_i)} and {r(k; Θ_i)} for a range of k.

4.8. Bayesian priors

For the prior densities p_Θ(Θ) = ∏_ip_Θi(Θ_i) and p_x(x) = ∏_ip_xi(x_i) in Eq. (28), we assert, similar to the analysis in SH18, uniform top-hat priors and the boundaries V(m)∈[max{ − 1, −⟨N|m⟩}, 1], ζ ∈ [0.1, 2], ${log}_{10}(m_{\text{piv}}h{\text{M}}_{\odot }^{-1})\in [4,16]$, b_ls ∈ [0, 5], r_ls ∈ [0, 2], and f_cen ∈ [0, 1].

An exception here is the prior density for b_i := b(m_i), which compared to SH18 now has the less informative density p_bi(b_i)∝1/b_i or, equivalently, $p_{b_i}(ln{b}_i)=\mathrm{const}$. This change is needed to avoid a bias in the r(k) reconstruction for k ≳ 3 h⁻¹ Mpc⁻¹ due to the noisier CFHTLenS data (for more details on this problem, see Appendix E). Likewise, the prior density for x_i is now p_x(x)∝1/∏_ix_i.

4.9. Systematic error of the galaxy-bias amplitude

Errors in the fiducial cosmology or the redshift distributions of lenses and sources (i.e. in the projection kernel) cause systematic errors in the amplitude of the reconstructed biasing functions b(k) and r(k). Since changes in the projection kernel produce variations in the amplitudes of ${\widehat{b}}_{2\mathrm{D}}$ and ${\widehat{r}}_{2\mathrm{D}}$ with only weak dependence on angular scale, the k-dependence of errors in the deprojected biasing functions is also weak, and the amplitude errors in the deprojection are of similar magnitude. We therefore estimate errors in b(k) and r(k) by propagating variations in the lensing kernel parameters to ${\widehat{b}}_{2\mathrm{D}}$ and ${\widehat{r}}_{2\mathrm{D}}$ for unbiased galaxies (i.e. the normalisation for the projected bias), and we average their RMS variations over 1′≲θ_ap ≲ 2°.

As baseline for the error model we assume the RMS values in Table 6 for (i) the mean and width of the redshift distributions, (ii) the set of fiducial cosmological parameters, (iii) the baryon physics in the non-linear matter power spectrum, and (iv) the amplitude $A_{\mathrm{ia}}=-0.{48}_{-0.87}^{+0.75}$ of the IA, based on the CFHTLenS constraints in Heymans et al. (2013). The error distribution of cosmological parameters for (ii) is the posterior distribution of the Planck TT, TE, and EE cosmology constraints combined with data from BAO experiments; the reference cosmology is a flat wCDM model (Planck Collaboration XIII 2016)⁷. For the fiducial redshift distribution of the lenses in (i), we take that of SM1, SM4, SM6, RED, and BLUE, and we then average the results. The redshift distribution of the sources are either BACK low-z or BACK high-z.

We perform thousands of random variations in (i) to (iv) and average for all variations the quadratic means in the interval 1′≲θ_ap ≲ 2° for $\delta b={\widehat{b}}_{2\mathrm{D}}/{\widehat{b}}_{2\mathrm{D}}^{\mathrm{fid}}-1$ and $\delta r={\widehat{r}}_{2\mathrm{D}}/{\widehat{r}}_{2\mathrm{D}}^{\mathrm{fid}}-1$ relative to our fiducial models ${\widehat{b}}_{2\mathrm{D}}^{\mathrm{fid}}$ and ${\widehat{r}}_{2\mathrm{D}}^{\mathrm{fid}}$ in Table 4 and the redshift distributions in Fig. 1. For every new variation in the cosmological model we draw a set of cosmological parameters from the published Planck MCMC posterior, thereby accounting for the correlations of their errors. Variations in the mean (δz) and RMS width (δσ) of redshift distributions of lenses and sources, the IA amplitude (δA_ia) and the baryon transfer function (δF) are independently drawn from a Gaussian distribution with the variances given in Table 6. To change the mean of a redshift distribution p(z), we map p(z) → p(z[1 + δz]), whereas p(z) → p(z)^{1/(1 − δσ)2} varies the width of the distribution (Sect. 7.3 in SH18). We perform this analysis for a scenario where all errors (i)–(iv) are present simultaneously but are statistically independent, and scenarios where only one error source (i)–(iv) is switched on at a time, giving the isolated RMS error σ_i. The detailed error distributions of the bias amplitudes can be found in Appendix F.

In summary, the RMS error σ_tot for b(k) ranges between 6.1% and 11.6%, and between 4.3% and 6.9% for r(k). For these values, the lower limits adopt uncorrelated systematic errors (i) to (iv) (the optimistic scenario, i.e. ${\sigma }_{\text{tot}}^2={\sum }_i{\sigma }_i^2$), while the upper limits assume that systematic errors (i) to (iv) conspire to maximise the amplitude error (the conservative scenario, i.e. σ_tot = ∑_iσ_i).

However, these are not all the systematic errors that need to be considered. There may be additional systematic errors from the methodology in Sect. 4, for instance model bias by the templates or the likelihood, and its specific implementation in a computer code. According to SH18, these errors fall between 3% and 7% for b(k) and between 3% and 5% for r(k). These estimates were obtained from galaxy-bias reconstructions with mock lensing data (H15 based) that have exactly known galaxy bias, redshift distributions, and cosmological and IA parameters. Since they are related to methodology, they are uncorrelated to uncertainties in the galaxy bias normalisation.

Combining all the systematic errors above, the accuracy for b(k) ranges somewhere between 6.8% (optimistic) and 13.5% (conservative), and for r(k) between 5.2% and 8.5%. A systematic error would be similar on all scales, and therefore offset the inferred biasing functions.

5. Results

5.1. Biasing functions

Our b(k) and r(k) results are shown on the following pages, after the references, as Figs. A.1–A.4 in the Appendix. The orange regions depict the 68% and 95% CIs of the reconstructed biasing functions around the median (dashed lines). The reconstructions assume the fiducial parameters in Table 4 and the galaxy redshift distributions in Fig. 1. There is one panel for each galaxy sample, indicated by the in-panel labels at the top left. For most CFHTLenS samples, the bias factor b(k) clearly changes with k, while the high galaxy-matter correlation, r(k) ≈ 1, has a moderate scale dependence.

The green data points show the SAM predictions in comparison: filled diamonds are for H15, and filled triangles are those for L12. In comparison to the models, the relative changes of the CFHTLenS b(k) with k is broadly similar to the models, but real galaxies are clearly more biased for k ≲ 0.5 h Mpc⁻¹. Towards the small spatial scales at higher k, the biasing of galaxies in both CFHTLenS and the models increases. Due to a steeper increase of L12 for SM2 to SM5, however, the Durham model typically agrees better with CFHTLenS than H15. Still, the models have a too small bias for SM1 and BLUE on all scales. With regard to the matter-galaxy correlation, the value of r(k) evidently increases for k ≳ 1 h Mpc⁻¹ in the CFHTLenS samples SM3 to SM6. This reflects, despite the noise, the scale-dependence in the SAMs. The model agreement is somewhat worse for high-z where the high value of r(k) at k ≈ 10 h Mpc⁻¹ is not supported by the CFHTLenS samples SM1, SM2, RED, and BLUE.

For a test of signal robustness in our data we split the source catalogues into two subsamples of distinct mean redshifts $\overline{z}$: BACK low-z ($\overline{z}\approx 0.97$) into BACKa low-z ($\overline{z}\approx 0.86$) and BACKb low-z ($\overline{z}\approx 1.12$); BACK high-z ($\overline{z}\approx 1.03$) into BACKa high-z ($\overline{z}\approx 0.93$) and BACKb high-z ($\overline{z}\approx 1.14$). Then we repeat every galaxy-bias reconstruction twice using those subsamples. The overlaid transparent grey regions in the figures are the 68% posterior CIs about the median for these reconstructions. Since they constrain the same biasing functions as the combined source sample, they should overlap with the filled dark-orange region, and they should overlap with each other. This is indeed the case for most reconstructions. A few exceptions are found near k ≈ 2 h Mpc⁻¹ for b(k) or at k ≲ 0.1 h Mpc⁻¹ for r(k), such as for RED low-z in both cases. Another disagreement is visible at small scales for SM4 high-z where the two extra reconstructions prefer smaller b(k) values. On the whole, however, and with some disagreement to be expected due to random chance, we consider the reconstruction overall robust. The next section tests the adequacy of the templates in describing the ratio statistics.

5.2. Posterior predictive check of bias reconstruction

A posterior predictive check validates the template models in our Bayesian deprojection (Gelman et al. 2003). We construct the posterior predictive, plotted in Figs. 5 and 6, by drawing values $\mathbf{\Theta }\sim p(\mathbf{\Theta }|\mathit{d},{\overline{n}}_{\mathrm{g}})$ from the posterior distribution of template parameters in Eq. (28) and by comparing the scatter of their b_2D(θ_ap; b) (left panels) and r_2D(θ_ap; b, r) (right panels) to the CFHTLenS data points of the ratio statistics. The ratio statistics is normalised by the prediction for unbiased galaxies using the fiducial cosmology in Table 4. The 68% and 95% CIs are indicated by the shaded regions. We plot, for reference, b_2D = r_2D = 1 as dotted lines. With the exception of r_2D at θ_ap ≈ 6 arcmin for high-z, the replicated and CFHTLenS data agree very well for the (normalised) ratio statistics b_2D and r_2D. The origin of the r_2D misfit is unclear, but it may be a bias of an overly smooth template or a systematic error in the shear data at this particular angular scale.

Fig. 5. Posterior predictive checks for the low-z CFHTLenS samples. The data points are normalised bias measurements b2D (left panel) and r2D (right panel) with jackknife errors as a function of aperture size. The orange regions are the posterior predictive distributions of our templates with 68% and 95% CIs; the dashed lines are the median. The green solid lines depict the normalised bias and 68% confidence intervals for H15 mock data in SH18 without shape noise.

Fig. 6. As in Fig. 5, but for the high-z samples. There are signs of a model misfit to r2D at around 6 arcmin, especially for RED and SM5. There is also weak evidence in favour of a stochastic bias, r2D <  1, for SM1 and BLUE at degree scales (bottom right panels).

The projected galaxy bias in Figs. 5 and 6 nicely illustrates that the conflict with the SAM galaxy bias is already present before our deprojection of the ratio statistics; it has nothing to do with the specific templates in our analysis. For comparison with the CFHTLenS data points, we plot the normalised statistics and their 68% confidence intervals for the H15-based mock galaxy and shear catalogues without shape noise, used in SH18, as green solid lines (cosmological parameters as in the Millennium Simulation). Throughout, the H15 data for b_2D are systematically lower than CFHTLenS on all angular scales. The correlation factor, r_2D, on the other hand, is broadly similar to CFHTLenS, with the notable exception of SM1 high-z. Similar data for L12 are not available to us. However, we do not expect their predictions to be significantly different to H15 due to their similar biasing functions for k ≲ 5 h Mpc⁻¹ (cf. green data points in Figs. A.1–A.4).

5.3. Galaxy bias on linear scales

Our technique separates galaxy bias in the non-linear (one-halo) regime from the linear stochastic bias b_ls and r_ls on large scales. The CFHTLenS results and model predictions are plotted in Fig. 7. Diamonds with error bars in the left panel are CFHTLenS posterior medians of b_ls with 68% CIs or r_ls in the right panel. Filled data points denote the low-z samples and open data points the high-z samples. The corresponding galaxy samples are shown at the bottom of each panel and the mean stellar masses at the top. The shaded rectangular regions indicate the 68% intervals of the average bias, obtained from merged samples SM1-SM6 in an extra analysis.

Fig. 7. Large-scale bias bls (left panel) and correlation rls (right panel) for the galaxy populations indicated on the x-axis. The black data points are CFHTLenS results (medians and 68% CIs), the green data points in the left panel are SAM predictions (diamonds for H15, triangles for L12, filled symbols for low-z, and open symbols for high-z). The shaded areas are 68% CIs for CFHTLenS, merging SM1 to SM6 in an extra analysis. The blue squares are VIPERS bls for stellar-mass binned galaxies with iAB <  22.5 and median $\overline{z}=0.61$ by Marulli et al. (2013).

On the one hand, for b_ls there are clear differences from the model predictions, which are the green diamonds (H15) and triangles (L12): model galaxies are less clustered relative to matter. On the other hand, the relative bias b_i/b_j for any two samples i and j is similar for CFHTLenS and the models. This can be seen when multiplying the H15 or L12 values by roughly 1.5, bringing them into broad agreement with CFHTLenS. Furthermore, the models predict an evolution of biasing from low-z (open symbols) to high-z (filled symbols) that is too small for the CFHTLenS errors. Nevertheless, our merged samples SM1-SM6 weakly indicate an evolution but with the opposite trend to in the models, namely a decline in b_ls (68% credibility) from $1.{98}_{-0.18}^{+0.20}$ at low-z to $1.{65}_{-0.12}^{+0.17}$ at high-z.

The blue squares in the left panel are non-lensing measurements of b_ls by Marulli et al. (2013) in the VIMOS Public Extragalactic Redshift Survey (VIPERS) for galaxy samples similar to ours at high-z (Sect. 2.4). For their estimate of b_ls, the VIPERS samples were divided into three stellar-mass bins, encompassing ${10}^{9.8}{h}_{70}^{-2}{\mathrm{M}}_{\odot }$ to ${10}^{11.8}{h}_{70}^{-2}{\mathrm{M}}_{\odot }$, and the galaxy-clustering correlation function between 1 and 10 h⁻¹ Mpc was normalised by the theoretical clustering of dark matter for a flat ΛCDM universe (Ω_m = 0.25, σ₈ = 0.8, n_s = 1). To make these results comparable to ours, we interpolate the VIPERS data points (their Table 4) to the mean stellar masses of SM2 to SM6 by using a best-fitting second-order polynomial. Then, except for our higher SM2 data point, the VIPERS results agree well with our lensing-based high-z measurements (open diamonds) and conflict with the SAMs for the highest stellar masses SM5 and SM6.

The linear correlation factors in the right panel are 68% consistent with a deterministic bias r_ls = 1 (and the SAMs), except for two cases: BLUE and, in particular, SM1 high-z significantly fall below this value. Since SM1 high-z is an interesting case, Fig. 8 (right panel) plots its ${\widehat{r}}_{2\mathrm{D}}$ as function of aperture size. The lines are predictions for unbiased galaxies with IA amplitudes varying between −3 ≤ A_ia ≤ 3. Importantly, all CFHTLenS galaxies with r_ls = 1, bias or unbiased, should for large apertures be consistent with our fiducial value A_ia = −0.48. While this is clearly the case for SM1 low-z (left panel), the data points of SM1 high-z in the right panel fall even below A_ia = −3 for θ_ap ≳ 30 arcmin, reflecting our low value for r_ls in the Bayesian interpretation of the aperture statistics.

Fig. 8. Ratio statistics ${\widehat{r}}_{2\mathrm{D}}$ as a function of aperture size. The data points with error bars (jackknife 68% errors) are for SM1 low-z (left panel A) and SM1 high-z (right panel B). The lines are predictions for unbiased galaxies with varying values of Aia ∈ [ − 3, +3]. Our fiducial value is Aia = −0.48. The solid magenta line ‘no MB’ is a model for Aia = 0 and no magnification.

5.4. Template parameters

Our main results are the reconstructed biasing functions for CFHTLenS in Figs. A.1–A.4. Nevertheless, marginal constraints for the template parameters in the one-halo regime (i.e. b(m), V(m), f_cen, ζ, and m_piv) are also reported in Appendix G and are briefly summarised here. These inform us about physical galaxy parameters, albeit with significantly lower accuracy compared to the reconstructed b(k), r(k), and their asymptotic values b_ls and r_ls on large scales.

Typically, we find only weak constraints for HOD-related template parameters in the CFHTLenS data, yet differences to the SAMs are visible for b(m): CFHTLenS galaxies of same stellar mass seemingly prefer to populate halos of higher mass when compared to the true b(m) of the SAMs. However, this shift of the inferred CFHTLenS b(m) towards higher halo masses is at least partly an artefact of the crude template modelling because it is also visible in the simulated reconstruction of biasing functions in SH18, using H15-based mock data. We also note here that the specific results for b(m) (and V(m)) rely on the adopted halo-mass spectrum n(m) which differs for CFHTLenS and the Millennium Simulation due to different fiducial cosmologies. Furthermore, the normalised excess variance V(m) is mostly consistent with a Poisson variance of galaxy numbers inside halos (V(m) = 0) or possibly a sub-Poisson variance (V(m) < 0) around 10¹³ ≲ m_piv ≲ 10¹⁴ h⁻¹ M_⊙ in some cases, SM6 high-z for instance. In comparison, the V(m) for the SAMs also indicates a sub-Poisson variance but again offset to somewhat lower halo masses. For high halo masses m ≳ 10¹⁴ h⁻¹ M_⊙ and stellar masses below SM5 the variance of galaxy numbers N(m) becomes super-Poisson (V(m) > 0) for H15. This is neither seen in L12 nor visible in our CFHTLenS results. Then the distribution of satellites is consistent with that of matter inside halos (ζ ∼ 1) within typical errors of roundabout ±0.45 (68% CI). And, finally, the mixing parameter of templates with and without central galaxies, f_cen, falls broadly within 20% and 80%, including the high values of 70%≲f_cen ≲ 90% in the SAMs.

6. Discussion

Our results for the galaxy bias in Figs. 7 and A.1–A.4 (Appendix) clearly show a conflict with the SAM predictions by H15 (Munich model) and L12 (Durham model). This conflict hints at inaccuracies in the SAM modelling of galaxy physics through its connection with galaxy bias (e.g. Weinberg et al. 2004; Blanton et al. 1999). There is no clear preference for any model by CFHTLenS, despite the inter-model differences on small scales, most prominent for b(k) in SM3 to SM5. We discuss our findings and possible systematic errors below.

The linear bias factor b_ls is overall 50% higher relative to the models due to a stronger clustering of galaxies in CFHTLenS. This may highlight modelling problems for z ≳ 0.3 that affect all stellar masses. We measure for b(k) and b_ls the amplitude ratio of galaxy and matter clustering in the data. Comparing the amplitudes of ⟨𝒩²⟩ in CFHTLenS on scales 4′≤θ_ap ≤ 40′ to the H15 mock data in SH18, we find typical values in CFHTLenS that are 20% to 40% higher. This conflict is further increased in b_ls by a weaker clustering of matter in CFHTLenS, ${\mathrm{\Sigma }}_8={\sigma }_8\sqrt{{\mathrm{\Omega }}_{\mathrm{m}}/0.27}\approx 0.76$, relative to that in the Millennium Simulation, Σ₈ ≈ 0.87. Problems with galaxy clustering in the Durham and Munich models are known for low-mass galaxies at z ≲ 0.1 on small scales of ∼ 1 h⁻¹ Mpc (Henriques et al. 2017; van Daalen et al. 2016; Farrow et al. 2015; Kang 2014). According to Fig. 7, supported by the VIPERS data for SM3-SM6, this problem deepens in the range 0.3 ≲ z ≲ 0.6 and then also affects linear scales and stellar masses as high as $M_{*}\sim {10}^{11}{h}_{70}^{-2}{\mathrm{M}}_{\odot }$. In a halo-model picture the lower b_ls suggests that H15 or L12 galaxies tend to reside inside matter halos of systematically low mass, decreasing the clustering amplitude on linear scales (Berlind & Weinberg 2002). At first sight this seems to be supported by the constraints of our template parameters for b(m)∝⟨N|m⟩/m (Fig. G.2), but we caution here that the templates, which are employed for a different purpose, have not been tested for their ability to truthfully reconstruct b(m). Clearly visible in Fig. 2, however, is that the SMF of the SAMs in the two redshift bins closely matches that in CFHTLenS. First of all, this supports a selection of galaxy populations compatible with CFHTLenS; secondly, if there are differences in the SAMs’ ⟨N|m⟩, then they have to be such that we nevertheless obtain similar galaxy number densities ${\overline{n}}_{\mathrm{g}}={\int }_0^{\infty }\mathrm{d}mn(m)⟨N|m⟩$. In other words, variations δ⟨N|m⟩ between the SAMs and CFHTLenS have to satisfy ${\int }_0^{\infty }\text{d}mn(m)\delta \langle N|m\rangle \approx 0$. The large-scale bias factor, on the other hand, can be changed by $\delta b_{\mathrm{ls}}={\overline{n}}_{\mathrm{g}}^{-1}{\int }_0^{\infty }\mathrm{d}mn(m)b_{\mathrm{h}}(m)\delta ⟨N|m⟩$ while keeping ${\overline{n}}_{\mathrm{g}}$ fixed. To identify HOD variations between CFHTLenS and the SAMs, permitted by this δ⟨N|m⟩ constraint, let us crudely assume that all galaxies in a sample live inside halos of characteristic mass m_*, hence exhibiting the large-scale bias b_ls = b_h(m_*). In Fig. 9, we explore the relative m_* differences needed to explain the b_ls of observational data and models, using the same reference cosmology for both to make them comparable (see that of the Millennium Simulation, Sect. 3.1). For the halo bias factor b_h(m), we employ Tinker et al. (2005) at the mean redshift of a galaxy sample. On average, the m_* in CFHTLenS is broadly a factor of 6(3) higher than the average SAM at low-z (high-z), as indicated by the solid (dashed) line; in the SAMs m_* increases from around 4 × 10¹² h⁻¹ M_⊙ for SM1 to 2 × 10¹³ h⁻¹ M_⊙ for SM6. In conclusion, our results hint at problems in H15 and L12 for z ≳ 0.3. These problems may be related to inaccuracies in their HOD, preferring parent halos of 3 to 6 times lower mass than in reality. This view, however, assumes well-controlled systematic errors in our galaxy selection and the reconstruction of galaxy bias.

Fig. 9. Characteristic halo mass m* that corresponds to the large-scale bias bls in Fig. 7 for the halo bias model by Tinker et al. (2005) and the Millennium Simulation cosmology as reference (Sect. 3.1). Data points with error bars indicate the CFHTLenS galaxy samples (x-axis) and the green points are the SAM predictions. The green lines average the SAM values and rescale them by a factor of 6 (solid line for low-z) or 3 (dashed line for high-z).

In this respect it is conceivable that the conflicting results might also indicate large systematic errors in our galaxy bias amplitude combined with an inconsistent selection of SAM galaxies, such as systematically low CFHTLenS stellar masses. The estimated systematic RMS error in b_ls is somewhere within [6.8%,13.5%] (Sect. 4.9). This is well below the b_ls conflict of ∼50%, and therefore cosmology or redshift errors alone probably do not explain the model conflict. But two other systematic errors in the galaxy sample selection can also bias our results. Firstly, our flux limit of i′< 22.5 might remove too many faint weakly clustered galaxies in CFHTLenS if the SAM i′-band magnitudes are too faint relative to CFHTLenS (the flux limit is also applied to the SAMs). Even then the SAM conflict for the massive (luminous) galaxies $M_{*}\gtrsim {10}^{11}{h}_{70}^{-2}{\mathrm{M}}_{\odot }$ is hard to explain where the flux limit is irrelevant. In addition, the SAM luminosity functions give good fits to observations for z ≲ 1 (e.g. Bower et al. 2006); the Durham model even uses luminosity functions in two bands for the model calibration. We note that the low-mass samples SM1 and SM2 are especially affected by the flux limit, implying a higher value for b_ls compared to volume-limited stellar-mass samples because the bias factor increases with galaxy luminosity (e.g. Zehavi et al. 2005; Norberg et al. 2002). More relevant could be a second selection bias: stellar-mass estimates might be systematically low in CFHTLenS. Shifting the data points in the left panel of Fig. 7 by Δ log₁₀M_* = +0.6 dex gives a better match to the SAM predictions for b_ls (e.g. SM1 to SM3, SM2 to SM4). On the other hand, a stellar-mass error that large is hard to justify considering previous error estimates (Sect. 2.1 of Velander et al. 2014), and considering that the VIPERS results for both b_ls and ${\overline{n}}_{\mathrm{g}}$, using full spectroscopic information, are comparable to ours in the high-z bin (Marulli et al. in Fig. 7; Davidzon et al. in Fig. 2). In addition, the stellar-mass functions of the SAMs and CFHTLenS in Fig. 2 are good matches, which is unlikely for systematically different CFHTLenS stellar masses. Nevertheless, a large yet realistic stellar-mass error of ∼ 0.3 dex combined with a very large systematic error of ∼25% for b_ls could in principle resolve the SAM conflict for the large-scale galaxy clustering. A systematic error of such magnitude in b_ls, however, would imply that all our results for r_ls ≈ 1 in Fig. 7 have to be offset by about 15%, which would then mildly conflict with a deterministic galaxy bias on large scales; systematic errors in r_ls are typically 50% to 70% of those in b_ls (Appendix F).

For systematic errors within the estimated range, our results for r_ls agree with the SAM picture of a deterministic galaxy bias. Only the high-z samples SM1 and BLUE are weak evidence for a stochastic linear bias, which will be tested with upcoming survey data. Although the biasing functions b(k) and r(k) generally do not distinguish non-linear deterministic bias from a stochastic bias, the linear regime is an exception, where ${\stackrel{\sim }{\delta }}_{\mathrm{g}}\approx b_{\mathrm{ls}}{\stackrel{\sim }{\delta }}_{\mathrm{m}}+ϵ$ and r_ls is only sensitive to the shot-noise corrected stochastic scatter ϵ (Tegmark & Peebles 1998; Dekel & Lahav 1999). For this regime, and after the IA (and magnification) correction, the data points in the right panel of Fig. 7 establish, in agreement with the SAMs, a deterministic bias r_ls = 1 for CFHTLenS galaxies on linear scales. More interestingly, and contrary to the SAMs, the samples SM1 and BLUE high-z prefer an average stochastic bias r_ls ≈ 0.75 ± 0.14 (stat.) ± 0.06 (sys.) at 68% credibility. Some stochasticity of blue galaxies, r ≈ 0.9 on scales of 15 h⁻¹ Mpc, is supported by other observations, but it is unclear if this is still present on linear scales and for our redshift regime (Patej & Eisenstein 2016; Swanson et al. 2008; Wild et al. 2005; Blanton 2000). At least for the CMB lensing results of i-band limited galaxies by Giannantonio et al. (2016), a linear stochastic bias is one interpretation (Sect. 7.4). Unfortunately, our new evidence for stochasticity (2σ) is further weakened by the systematic amplitude uncertainty of [5.2%,8.5%], quoted above as ‘sys.’ for the pessimistic 8.5% (Sect. 4.9). In the near future we anticipate better constraints from the Kilo Degree Survey (de Jong et al. 2013), the Hyper Suprime-Cam Survey (Aihara et al. 2018), or the Dark Energy Survey (Abbott et al. 2016) which will have ten times the survey area of CFHTLenS and a better control of galaxy redshift distributions, which are, in addition to IA, the major source for the systematic error in r_ls (bottom panels in Fig. F.1).

Turning to the non-linear regime, CFHTLenS galaxy biasing clearly varies with scale, stellar mass, and colour, but differently from the SAMs and with no clear preference for either H15 or L12. Nevertheless, there is a qualitative agreement with the models insofar as the bias factor b(k) increases on non-linear scales k ∼ 10 h Mpc⁻¹ relative to linear scales by up to a factor of two (Figs. A.1 and A.2). The relative increase of r(k) is weaker, typically a factor of 1.2 (Figs. A.3 and A.4). These relative changes are insensitive to systematic errors in our method which offset the biasing functions as a whole. Physically the b(k) increase is probably related to central galaxies inside halos that have ⟨N|m⟩≲1 and dominate the signal for sufficiently large k (Sect. 5.3 in SH18). Central galaxies affect r(k) in a similar way, like probably the stellar mass samples of CFHTLenS at low-z. For other samples, however, where the r(k) increase is small or absent in contrast to the SAMs, this interpretation is less clear, RED and BLUE high-z for instance. For a physical interpretation of these observations, an accurate halo-model fit (e.g. Dvornik et al. 2018; Seljak et al. 2005) or a study of b(k) and r(k) variations with respect to SAM parameters is needed.

Finally, in a broader context, the possibility of a stochastic galaxy bias, IA distortions, or magnification effects on linear scales is also relevant for the popular cosmological probe E_G for gravity (Reyes et al. 2010). This probe relies on a deterministic bias, whereas a stochastic bias changes the value of E_G to r_ls E_G. Additionally, lensing magnification changes E_G through the galaxy-galaxy lensing (GGL) amplitude by up to several per cent, depending on type and mean redshift of lens galaxies (Fig. 4 in SH18). Likewise, GI correlations bias the GGL signal if the distributions of lenses and sources overlap in radial direction (Hirata & Seljak 2004). SH18 estimate the GI bias of GGL, and thus of E_G, to up to 5% for −2 ≲ A_ia ≲ −1 and lens, source samples similar to ours (Fig. 3 in SH18). While this GI contamination can in principle be reduced to zero by avoiding any overlap, the magnification bias is always present and has to be accounted for by applying the method in Unruh et al. (2020), for example. Importantly, our work suggests not to use BLUE (star-forming) or SM1-like galaxies as lenses for E_G due to the possibility of a stochastic bias.

At the CDS, we provide Monte Carlo realisations of the posterior distributions of b(k) and r(k) for every galaxy sample in our analysis. The scatter among the realisations reflects the distribution of their statistical errors, and the full set of realisations can be utilised to propagate this uncertainty in future studies.

Acknowledgments

We thank Christopher Bonnett for fruitful discussions at an early stage of this project, and we thank Peter Schneider and Catherine Heymans for comments on the paper manuscript. This work has been supported by Collaborative Research Center TR33 ‘The Dark Universe’ and by the Deutsche Forschungsgemeinschaft through the project SI 1769/1-1. This work is based on observations obtained with MegaPrime/MegaCam, a joint project of CFHT and CEA/DAPNIA, at the Canada-France-Hawaii Telescope (CFHT) which is operated by the National Research Council (NRC) of Canada, the Institut National des Sciences de l’Univers of the Centre National de la Recherche Scientifique (CNRS) of France, and the University of Hawaii. This research used the facilities of the Canadian Astronomy Data Centre operated by the National Research Council of Canada with the support of the Canadian Space Agency. CFHTLenS data processing was made possible thanks to significant computing support from the NSERC Research Tools and Instruments grant program.

References

Appendix A: Reconstructed biasing functions

Fig. A.1. Posterior scale-dependent bias b(k) for galaxy samples at redshift $\overline{z}\approx 0.35$. The galaxy sample is indicated in the top left corner. The orange filled regions and the dashed lines show the median, 68%, and 95% credible region for the CFHTLenS measurement using the full source sample BACK low-z. The overlapping transparent grey regions show 68% regions that are based on the source subsamples BACKa and BACKb low-z. These noisier measurements at mean source redshifts $\overline{z}\approx 0.86$ and $\overline{z}\approx 1.12$ test the robustness of the measurement. The filled green diamonds and triangles indicate the SAM predictions by H15 and L12, respectively.

Fig. A.2. As in Fig. A.1, but for the scale-dependent bias b(k) of the high-z galaxy samples ($\overline{z}\approx 0.52$). The grey regions depict 68% credible regions, based on the source subsamples BACKa and BACKb high-z ($\overline{z}\approx 0.93,1.14$). The green data points are SAM predictions.

Fig. A.3. As in Fig. A.1, but for the scale-dependent correlation factor r(k) of the low-z galaxy samples ($\overline{z}\approx 0.35$). The grey 68% regions are based on two source subsamples. The green data points are SAM predictions.

Fig. A.4. As in Fig. A.1, but for the scale-dependent correlation factor r(k) of the high-z samples ($\overline{z}\approx 0.52$). The grey 68% regions are based on two source subsamples. The green data points are SAM predictions.

Appendix B: Deconvolution of binned redshift distributions

To estimate the redshift distribution of our lens galaxies from the binned BPZ posterior distribution we discuss here a maximum-likelihood technique. Generally, the binning of a probability density p(z) dz into i = 0…(N_bin − 1) bins,

$$\begin{array}{c}\hfill p_i:={\displaystyle \underset{i\mathrm{\Delta }z}{\overset{(i+1)\mathrm{\Delta }z}{\int }}\mathrm{d}zp(z),}\end{array}$$(B.1)

is not invertible with respect to p(z). By restricting the inversion for z >  0 to a family of densities of four parameters $\mathbf{\Theta }=(\overline{z},{\sigma }_z,s,k)$,

$$\begin{array}{c}\hfill p(z|\lambda ,\mathbf{\Theta })=\lambda {\mathrm{e}}^{-\frac{ϵ}{z}}{\mathrm{e}}^{-\frac{x^2}{2}}(1+\frac{s}{6}H{e}_3(x)+\frac{k}{24}H{e}_4(x));x=\frac{z-\overline{z}}{{\sigma }_z},\end{array}$$(B.2)

however, the problem is tractable for our lens samples; we set p(z|λ, Θ) = 0 for z ≤ 0. The normalisation λ is defined through ${\int }_0^{\infty }\text{d}z{p}_{\text{d}}(z|\lambda ,\Theta )=1$, and the functions He₃(x) = x³ − 3x and He₄(x) = x⁴ − 6x² + 3 are Hermite polynomials. This family of models (B.2) is a Gram-Charlier series that approximates the narrow redshift distributions of samples inside photometric redshift bins by a normal distribution and corrections for a skewness or kurtosis. In addition, the exponential function with ϵ = 10⁻² asserts p_d(z|λ, Θ) → 0 for z → 0. This correction influences the distribution function only very close to z = 0. The aim of our technique is to determine the parameters Θ for given p_i.

To determine the best-fitting parameters in (B.2) for a given BPZ posterior, we find the binned distribution

$$\begin{array}{c}\hfill {\widehat{p}}_i(\lambda ,\mathbf{\Theta }):={\displaystyle \underset{i\mathrm{\Delta }z}{\overset{(i+1)\mathrm{\Delta }z}{\int }}\mathrm{d}zp(z|\lambda ,\mathbf{\Theta })}\end{array}$$(B.3)

that has the minimal Kullback-Leibler divergence with respect to p_i, which means finding the Θ_ml that minimise

$$\begin{array}{cc}& d_{\mathrm{KL}}(\mathbf{\Theta })=-{\displaystyle \sum _{i=0}^{N_{\mathrm{bin}}-1}{p}_{i}ln\left(\frac{{\widehat{p}}_i(\lambda ,\mathbf{\Theta })}{p_i}\right)}\hfill \\ \hfill & =\mathrm{const}-{\displaystyle \sum _{i=0}^{N_{\mathrm{bin}}-1}{p}_{i}ln{\widehat{p}}_i(\lambda ,\mathbf{\Theta }).}\hfill \end{array}$$(B.4)

These values are the maximum-likelihood parameters in the family of model distributions (Knight 1999).

Figure B.1 shows two examples of this fitting procedure for a galaxy sample SM1-SM6 in the two photo-z bins low-z and high-z. In this figure the dashed histograms are the binned BPZ posteriors p_i, the solid histograms are the best-fitting binned models ${\widehat{p}}_i$, and the solid orange lines are the best-fitting p(z|λ, ⃗Θ) used in the analysis.

Fig. B.1. Two examples of maximum-likelihood fits of models of p(z) to the combined binned posterior of the BPZ redshift. The left-hand peak is for the low-z bin of lenses SM1 to SM6 combined; the right-hand peak corresponds to the high-z bin. The dashed histograms are the average BPZ posteriors of lenses; the solid histograms are the binned best fits. The smooth curves are the unbinned best fits as used in the analysis.

Appendix C: V_max-estimator for photo-z selected galaxies

For our modification of the V_max-estimator, we consider galaxies that are selected in magnitude m_min ≤ m <  m_max and from a photo-z interval z_min ≤ z_ph <  z_max. Due to the latter there is a non-negligible probability that galaxies below z ≤ z_min or above z ≥ z_max are also selected; the selection thus does not clearly define a survey volume. To implement the volume uncertainty in an estimator of the galaxy number density, we assume that we know of every galaxy i the probability density p_i(z|z_ph) of its true redshift z given its photometric redshift, and that we know the probability density p_i(z_ph|z) of photometric redshifts given the true galaxy redshift. We first study how galaxies are selected.

To define the selection of galaxies, let p_i(z, θ) dV = C_i(z) S(θ) dV be the probability that the ith galaxy can be observed inside a comoving volume element dV in the line-of-sight direction θ at redshift z. The volume element is in a cosmological context

$$\begin{array}{c}\hfill \mathrm{d}V=f_K^2(z)\frac{\mathrm{d}\chi (z)}{\mathrm{d}z}\mathrm{d}\mathrm{\Omega }\mathrm{d}z,\end{array}$$(C.1)