PROPAGATING RESIDUAL BIASES IN COSMIC SHEAR POWER SPECTRA

In this paper we derive a full expression for the propagation of multiplicative and additive shape measurement biases into the cosmic shear power spectrum. In doing so we identify several new terms that are associated with selection effects, as well as cross-correlation terms between the multiplicative and additive biases and the shear field. The computation of the resulting bias in the shear power spectrum scales as the fifth power of the maximum multipole considered. Consequently the calculation is unfeasible for large l-modes, and the only tractable way to assess the full impact of shape measurement biases on cosmic shear power spectrum is through forward modelling of the effects. To linear order in bias parameters the shear power spectrum is only affected by the mean of the multiplicative bias field over a survey and the cross correlation between the additive bias field and the shear field. If the mean multiplicative bias is zero then second order convolutive terms are expected to be orders of magnitude smaller.


INTRODUCTION
The statistical properties of the large-scale matter distribution over cosmic time encodes key information about the late time evolution of the Universe, and also allows us to improve constraints on the initial conditions. Thanks to technological advances we can now efficiently survey larger and larger areas of sky, but the interpretation of galaxy redshift surveys is hampered by the fact that galaxies are biased tracers of the underlying dark matter distribution. Fortunately, the distortion of space-time by matter results in correlations in the ellipticities of distant galaxies that are the result of the differential deflection of light rays, a phenomenon called gravitational lensing. The statistics of these correlations can be directly related to those of the large-scale structure. This in turn enables us to constrain the nature of dark energy and to test gravity on cosmological scales.
The cosmological lensing signal has now been robustly measured using large ground-based imaging surveys (e.g. Hildebrandt et al. 2018;Troxel et al. 2018). However to shed light on the nature of dark energy, the precision needs to increase significantly. This is the objective of a number of planned projects that will commence soon. Euclid (Laureijs et al. 2011) aims to image 15 000 deg 2 of extragalactic sky from space, while the Large Synoptic Survey Telescope (LSST) will survey a similar area from the ground. To exploit fully the potential of these data for cosmology, it is essential that astrophysical and t.kitching@ucl.ac.uk instrumental sources of biases are accounted for at levels that are small compared to the statistical uncertainties on measured cosmological parameters. Accurate measurements of the shapes of small, faint galaxies are therefore essential.
The observed ellipticities of galaxies used in weak lensing studies are typically biased with respect to the true ellipticities that would have been measured given ideal data and an ideal measurement algorithm. The dominant sources of bias are a result of the convolution by the point spread function (PSF) and noise in the images. For this reason the performance of shape measurements has been studied extensively (Heymans et al. 2006;Massey et al. 2007;Bridle et al. 2010;Kitching et al. 2012;Mandelbaum et al. 2015).
To first order the biases can be separated into multiplicative and additive functions that act on the true ellipticity. Additive biases arise from anisotropies in the data, such as an anisotropic PSF or detector effects. These do not only affect the measurement of the galaxy shape, but also the detection and selection of sources. Characterizing and correcting for these sources of bias is essential, but residual spurious alignments might still be removed through empirical corrections. For instance, the mean shear when averaged in the coordinate frame defined by the detector should vanish. The detection of a coherent signal would thus indicate an imperfect correction, but that signal could also be fitted for in the cosmological analysis. However we note that this is only partially effective because sources of additive bias are expected to introduce multiplicative biases with a similar amplitude. Unfortunately, multiplicative bias cannot be readily inferred from the imaging data directly. Instead image simulations are used to calibrate the biases in the shape measurement algorithms (e.g. Hoekstra et al. 2015;Kannawadi et al. 2018), although we note that alternative approaches have been recently proposed (Huff & Mandelbaum 2017a).
The desired accuracy in cosmological parameter estimates determines the level at which shape measurement biases can be tolerated. The propagation of biases, or residual biases after calibration, into the weak lensing power spectra (or 'cosmic shear' power spectrum) is not straightforward. Some studies approximate the full expression (Taylor & Kitching 2016;Kitching et al. 2012;Massey et al. 2013) but, as we show in this paper, these results do not capture the spatially varying sources of biases correctly. This is of particular importance because the theoretical propagation of such biases into power spectrum residuals drives the design requirements for experiments that use weak lensing as a cosmological probe (Cropper et al. 2013).
In this paper we show how multiplicative and additive biases in shape measurement propagate through the cosmic shear power spectra, discuss how this formalism relates to previous studies, and discuss the implications for the assessment of shape measurement biases on cosmological parameter performance verification of experiments. In Section 2 we present the formalism, in Section 3 we present some simple simulations that demonstrate the accuracy of the formalism, and in Section 4 we examine the implications of this study; conclusions are presented in Section 5.

METHOD
We begin with the expression for the measured shear in real (angular) space where γ(Ω) is the measured shear as a function of angle Ω = (θ, φ) where θ and φ are R.A. and declination respectively, m 0 is a constant multiplicative bias, m(Ω) = m R (Ω) + im I (Ω) is a position-dependent multiplicative bias term, γ(Ω) = γ 1 (Ω) + iγ 2 (Ω) is the true shear, c 1,0 and c 2,0 are constant additive biases that contribute to the real and imaginary parts of the additive field, and c(Ω) is a position-dependent additive bias. We assume no non-local terms, e.g. m(Ω ′ )γ(Ω), since such terms could be always re-written as a local per galaxy m(Ω)γ(Ω) term.
We discuss the choice of the multiplicative bias expression in Appendix A. We note that the choice to express the multiplicative effect as a product of complex numbers makes the spin-preserving assumption that any effect only changes the amplitude and/or angle of the observed ellipse relative to the unbiased case. There are more general expressions that can be used to capture biases that can occur in the case of anisotropic systematic effects, however image simulations suggest that our adopted approach is accurate for residual systematic effects after calibration with image simulations. In the case where no rotational change is present this reduces to a multiplication by a single scalar field m(Ω). We will revisit these assumptions later in the analysis.

Spherical Harmonic Representation
We now determine the spherical harmonic representation of a biased shear field. We adapt the methodology from CMB pseudo-C ℓ analysis here for the general bias case; in particular we follow Lewis et al. (2002); Brown et al. (2005); Zaldarriaga & Seljak (1997); Grain et al. (2012) but we generalize their formalism further to include general spin-2 bias functions and additive terms. In Kitching et al. (2012) a similar adaption was made, but under simplifying assumptions that did not capture the general case.
Since for cosmic shear the cosmological information is contained within the E-mode (gradient) component of the field, and not in the B-mode (curl) component, we work on spherical harmonic coefficients γ E ℓm and γ B ℓm . We can write down the E-and B-mode coefficients in terms of the shear field as where ℓ and m are angular wavenumbers1, 2 Y ℓm (Ω) is the standard spin-weighted spherical harmonic function for a spin-2 field and * denotes a complex conjugate. This expression is exact for an all-sky unbiased measurement of the shear. This formalism could be generalised to the pure-mode case (Grain et al. 2012), that would become important in the presence of masks, but we leave this masked data generalisation for future work.
To compute the effect of the biases we now replace γ(Ω) in equation (2) with γ(Ω) from equation (1). This results in the following expressions Here we expanded the additive term as c( and In this derivation we have expressed the real and imaginary parts of the multiplicative bias as m(Ω) = m R (Ω) + im I (Ω). We note that when considering residual systematic effects, i.e. when any amplitude and rotational changes caused by multiplicative systematic effects are small (see Appendix A), that m R (Ω) ≃ m(Ω) and m I (Ω) ≃ 02.
Already from the expressions in equation (3) it can be seen that multiplicative biases in general mix E and B-modes together, both from the underlying shear field and the multiplicative bias field, and the propagation of such terms is in the form of a convolution represented as a sum over wavenumbers. Furthermore the window function caused by multiplicative biases is ℓ and m-mode dependent since in general these are not isotropic on the celestial sphere.
We note that in this case the constant additive biases c 1,0 and c 2,0 do not appear in equation (3). This is because a constant term only affects the ℓ = 0 mode, but shear is a spin-2 field where the spherical harmonic transform is not defined for ℓ < 2; because s Y ℓm (Ω) = 0 for ℓ < |s|. Therefore a constant additive bias cannot affect the cosmic shear power spectrum.

Biased Cosmic Shear Power Spectra
We now compute the expressions for the biased cosmic shear power spectra by taking the correlation of the expressions in equation (3). The full expression can be written as a series of terms that pertain to multiplicative, additive and cross-terms, and depend on the true EE, EB and BB power spectra. The power spectra estimates are computed by taking the correlation of the spherical harmonic coefficients from equation (3) where We provide the full expanded expression for the biased power spectra in Appendix B. If we assume that C E B ℓ = 0, which is the case in all but the most exotic dark energy models, then the three estimated power spectra (EE, BB 2 We note that m R (Ω) ≃ m(Ω) and m I (Ω) ≃ 0 implies that, if one expresses the multiplicative biases as m 1 (Ω)γ 1 (Ω) + im 2 (Ω)γ 2 (Ω) (where 1 and 2 denote the ellipticity components measured parallel to Cartesian axes in a measurement frame, and measured at 45 degrees to these axes), m(Ω) = m 1 (Ω) = m 2 (Ω). This is found to be the case in state-of-the-art methods, e.g. Pujol et al. (2019). When considering residual systematic effects, after calibration with simulations, this is also expected to be the case. and EB) are: The various terms in the full expression are The power spectra in the full expression are labelled in their superscripts as either correlations between shear coefficients (E E, E B, BB), correlations between the addi-

or cross correlations between shear and additive bias terms (c
It can be shown that the N terms are simply the mean of the spatially varying multiplicative bias field. If we consider N + ℓ we find We can simplify this expression further by using the generalised addition theorem for spin-weighted spherical harmonics Grain et al. (2012) where (α, β, γ) are Euler angles between Ω and Ω ′ which in our case are zero, and D ℓ ss ′ (α, β, γ) are the Wigner ro-tation matrices which for D ℓ (0, 0, 0) = δ K ss . This leads to and similarly N − ℓ = m I (Ω) . We discuss further simplifications of these expressions below. In Appendix C we show the generalisation of this to the case of multiple tomographic bins.

Discussion of the terms
We can now discuss each term in the full expression and its physical meaning.
• m 0 and m 2 0 : These terms are the normal contribution from the constant multiplicative bias terms. These arise from the limitations with which shape measurement algorithms can be calibrated (see e.g. Hoekstra et al. 2017). As we are concerned with the residual biases after such a calibration, m 2 0 ≪ m 0 as m 0 2 × 10 −3 (Cropper et al. 2013) • N: These terms represent multiplicative biases, and are the mean of the multiplicative bias field.
• M: These terms represent multiplicative terms of order m 2 . The rows and columns show how the E and B-mode power have mixed terms, where the ++ terms pick up contributions from the real part of the multiplicative bias field, the −− terms pick up contributions from the imaginary part, and +− or −+ are mixed terms.
• B: These terms represent third-order, bispectrumlike, correlations between the position-dependent multiplicative bias m(Ω)γ(Ω) and the positiondependent additive bias c(Ω). Such effects are likely since areas in a survey, or particular pointings, that have detector, telescope or background effects that cause additive biases will also lead to multiplicative biases. This is because any anisotropic change in the quadrupole moments will modify the size, and thus the multiplicative bias.
• c E E, c E B, and c B B: these terms capture the correlations between the underlying shear field and the additive bias terms. Such terms are expected to be caused by selection effects in a real survey, where for example blending in high shear regions (e.g. about clusters) could cause an additive bias contribution.
• c E c E , c E c B , and c B c B : These are the power spectra of the position-dependent additive biases. We note again that constant additive bias terms do not contribute to cosmic shear power spectra.

Linear Expressions
Here we show the linearised expressions of the biased cosmic shear power, that only include terms that are linear in the bias parameters. We find that We have included B-mode power since as shown in Schneider et al. (2002) source redshift clustering can cause a small B-mode component. We see that the impact of spatially varying biases will be, to linear order, captured by the mean of the multiplicative bias and the additive-shear cross correlation power spectrum, but in the presence of intrinsic B-modes, C E B ℓ now includes a term 2 m I (Ω) C BB ℓ . We note that if m R (Ω) = m I (Ω) , then (twice) the EB power spectrum could be added to the EE power spectrum to cancel out multiplicative effects; however this is not expected to be the case in general, or for small biases.

SIMPLE SIMULATIONS
To test that the above formalism can indeed capture the propagation of general position-dependent multiplicative and additive bias terms into the cosmic shear power spectrum we generate several toy examples and investigate the contributions of each term to the overall change. For each case we define a multiplicative constant and field, m 0 and m(Ω), and an additive constant and field, c 0 = c 1,0 = c 2,0 and c(Ω), although the choice for c 0 has no impact on cosmic shear power spectra by definition. We normalise these fields such that m 0 + m(Ω) = 2 × 10 −3 and c 0 + c(Ω) = 1 × 10 −4 , which represent the overall requirements for a Euclid -like experiment (Cropper et al. 2013). For each case we compare the computation of the analytic expression in equation (7) and a numerical case where we compute the real space shear field γ(Ω) = [1+m 0 +m(Ω)]γ(Ω)+[c 0 +c(Ω)] and then compute the power spectra of this directly using a spherical harmonic transform. In all cases we compute the original γ(Ω) field using a Gaussian random field generated by using a cosmic shear power spectrum based on the Planck ΛCDM cosmology (Planck Collaboration et al. 2018), using the massmappy code (Wallis et al. 2017). In all cases we use SSHT McEwen et al. (2013) to compute the spin-weighted spherical harmonics, which sample the sphere using the sampling scheme of McEwen & Wiaux (2011).
The cases we consider are shown below. Note that we express these in terms of an arbitrary amplitude A since these are all normalised to have m 0 + m(Ω) = 2 × 10 −3 and c 0 + c(Ω) = 1 × 10 −4 . The cases are simple examples but nonetheless are approximations of realistic spatial variations that could occur: 1. Case 1, Simple Galactic Plane: • m R (Ω) = Ai, where i is an iterative pixel number count, which is reset when i = 10, m I (Ω) = 0, • c R (Ω) = c I (Ω) = Ai 2 , where i is an iterative pixel number count, which is reset when i = 10, • c 0 = A, m 0 = A. The first case approximates a Galactic plane dependency, the second case approximates a patch-dependent systematic effect, and the third case is a non-analytic case that approximates a scanning sequence of exposures. To demonstrate the complexity of the spatial variation of the cases we show in Figure 1 the real part of the multiplicative field for each of the cases (we do not show all the fields associated with the systematic effects since they are largely similar in form).
In Figure 2 we show the residual power spectra δC ℓ = C E E ℓ − C E E ℓ for each of the cases considered. We note that in all cases we use a maximum angular multipole of L = 32. This is because the calculations are particularly numerically demanding. The W ± ℓℓ ′ mm ′ calculations have dimension L 4 , and for each of these spin-weighted spherical harmonic functions must be computed each of which scale like L 2 log L at best (McEwen et al. 2013). This point is discussed further in Section 4.1. In all cases the analytic formula given in equation (7) accurately captured the form of the residual power spectrum; the very small differences are due to the numerical stability of the spin-weighted spherical harmonic transform calculations.
With regard to the different terms we find in all cases that the NC E E ℓ term is dominant, which is expected since it is of linear order in m(Ω), followed by the C c E E ℓ crosscorrelation term and the (m 0 +m 2 0 )C E E ℓ terms. The convolutive M, B and C c E c E ℓ terms are all at least an order of magnitude lower in all cases. Therefore the linearised expression in equation (11) is a good approximation in these simple examples. In the case that the mean of m(Ω) and c(Ω) are both zero, all of the remaining terms at second and third order in bias would become important at approximately the same level.

DISCUSSION
We have shown in general how constant and positiondependent shape measurement biases propagate through to cosmic shear power spectra.
The multiplicative bias terms shown here are similar to those that result in CMB polarisation pseudo-C ℓ analyses, where masking of the data results in expressions that also include temperature power spectra (Lewis et al. 2002;Zaldarriaga & Seljak 1997;Grain et al. 2012;Brown et al. 2005). The difference here is that instead of a mask we have a multiplicative bias field that is in general spin-dependent. We note that this formalism equally applies to the case of masked cosmic shear data where m(Ω) may be zero in some regions, and that such a case would lead to further mode mixing.
In Kitching et al. (2012) a pseudo-C ℓ formalism was used to assess position-dependent shape measurement errors. However in that study the linear terms N, bispectrum and additive terms were not included and the nonlinear convolution term M used a simpler form based on a flat-sky approximation. In Taylor & Kitching (2016) the propagation of shape measurement biases was generalised to the convolutive case but the linear terms, bispectrum and E/B mode mixing terms were ignored.
In Massey et al. (2013) and Cropper et al. (2013) requirements were set on weak lensing experiments using an approximation for position-dependent shape measurement biases. In that case a form of propagation was determined for the constant case C E E ℓ = (1 + 2m 0 + m 2 0 )C E E ℓ + c 2 c , which was then replaced with a 'positiondependent' formulation proposed by Amara & Réfrégier (2008), We find that such an expression is similar to the full case when only linear terms in biases are assumed, whereas the relationship to the underlying position-dependent bias fields is much more complex. Furthermore in Massey et al. (2013) a worst-case scenario in sensitivity was assumed for M ℓ in which multiplicative biases mimicked the scaledependent behaviour of dark energy. These worst-case assumptions are conservative when designing an experiment and lead to requirements that will guarantee performance, but when assessing the actual performance of a survey they are not adequate.
We note that if the mean multiplicative bias is zero -as may be expected if pre-experiment simulations can determine any mean effect -then only second and third order convolutive terms remain. The power spectrum residuals caused by these terms, and the impact on cosmology, are expected to be much lower than than the mean terms for two reasons. First, because the terms are second order and so for m(Ω) ≪ 1 these are smaller. Second, because they are convolutions it is unlikely that a functional form will result that matches the cosmic shear power spectrum; therefore the impact on cosmological parameter inference is expected to be lower (Taylor & Kitching 2016).
for the three cases considered. The left plots show the comparison between the numerical case computed by transforming the modified shear field and performing a spherical harmonic transform (i.e. a forward model), and the analytic case computed using equation (7). The right plots show the contribution to the analytic case from each of the components in equation (7). In all cases the legends label the coloured lines.

Scaling
In Figure 3 we show how the analytic and forward modelling cases scale as a function of L, the maximum ℓ-mode . To compute the full case requires evaluations of terms that scale like L 4 × L 2 log L, where L is the maximum multipole; this is because of the W ℓℓ ′ mm ′ terms that have ∼ L 4 summations, and the spherical harmonic transforms that scale like L 2 log L. We find a slightly better scaling due to pre-computation of the spin-weighted spherical harmonic functions, but nonetheless the analytic calcu-lation scales like ∝ L 5 compared to the forward modelling that scales like ∝ L 2 .
We evaluate simple examples for L = 32, but scaling to a reasonable value of L > 1000 would result in prohibitively long calculations. On the other hand the forward modelling of systematic effects i.e. the evaluation of equation (1) in real space and a direct spherical harmonic transform to produce a power spectrum is tractable for L > 1000 and we therefore advocate this approach in Taylor et al. (2019) and Paykari et al., (in prep).

CONCLUSIONS
In this paper we derive a complete expression for the impact of constant and spatially varying multiplicative and additive shape measurement biases on the cosmic shear power spectrum. In doing so we find several terms that have thus far been overlooked, in particular terms relating to cross-correlations between biases and shear.
In performing the full calculation we find that to linear order spatially varying biases are well approximated by the sum of the product of the the power spectrum and the mean of the multiplicative bias field, and the crosscorrelation term between the additive bias field and the shear field. We note that the cosmic shear power spectrum is not sensitive to constant or mean additive biases.
We compare the computation of the full analytic expression with that obtained using a forward modelling approach using simplified simulations and find good agreement. Furthermore we use these simplified simulations to demonstrate how each term in the full expression contributes to the total. However, in performing the full calculation we also find that its computation scales as the maximum multipole ∝ L 5 . This means that its evaluation for large L > 1000 is unfeasible. Therefore we recommend that any assessment of the impact of biases on cosmic shear power spectrum must be performed using a forward modelling approach.
We thank the developers of SSHT and massmappy, Jason McEwen, Chris Wallis and Boris Leistedt for making their code publicly available. We thank Dipak Munshi for useful discussions. We thank Peter Schneider for comments on an early draft. PP is supported by the UK Science and Technology Facilities Council. TK is supported by a Royal Society University Research Fellowship. HH acknowledges support from Vici grant 639.043.512 financed by the Netherlands Organization for Scientific Research (NWO).
where we have expressed γ(Ω) = |γ(Ω)|e i2Φ(Ω) where Φ(Ω) is the angle between the orientation of the elliptical shape induced by the shear and the x-axis of local Cartesian frame in which measurement has been made. To this expression we apply an amplitude change m(Ω) corresponding to an incorrect measurement of the ratio of the semi-major axes, and a small rotation φ m (Ω). This preserves the spin-2 nature of the measured field. When we express the product of two complex numbers m(Ω)γ(Ω), it should be understood that the multiplicative bias fields take the form where m R (Ω) and m I (Ω) in the first equation are the real and imaginary parts respectively which are coupled as expressed in the subsequent equations. We note that we do not label these m 1 (Ω) and m 2 (Ω) since they do not map solely to the γ 1 (Ω) and γ 2 (Ω) components of the shear field where γ(Ω) = γ 1 (Ω) + iγ 2 (Ω). We refer to Pujol et al. (2019); Huff & Mandelbaum (2017b) for further discussion of the propagation of more complex multiplicative biases. We note that if the amplitude of the (residual) biases are small, and that the rotation angle is random, then it is reasonable to assume that m R (Ω) = m I (Ω) ; if one assumes instead small residual biases and applies the small angle approximation then m R (Ω) ≈ m(Ω) and m I (Ω) ≈ m(Ω)φ m (Ω) ≈ 0. In the case that m R (Ω) = m(Ω) and m I (Ω) = 0 then this would result in m 1 (Ω) = m(Ω) and m 1 (Ω) = m 2 (Ω) if expressed in this way.
We note that if one applies a bias of the form m 1 (Ω)γ 1 (Ω) + im 2 (Ω)γ 2 (Ω) (i.e. a different independent scalar multiplicative biases applied to each of the shear components) then this cannot in general be expressed as an amplitude change with a rotation, and therefore can result in a change in the spin properties of the measured field. We note that one could create such an expression by m ′ (Ω)γ(Ω) + δm(Ω)γ * (Ω) where, m ′ (Ω) = [m 1 (Ω) + m 2 Ω)]/2 and δm(Ω) = [m 1 (Ω) − m 2 (Ω)]/2, however the second term would represent a parity change/mislabelling in the γ 2 component which would be a very large systematic effect. We explored the expected size of δm/m using image simulations that resemble Euclid based on Hoekstra et al. (2017) and found that δm/m ∼ 0.1 for both PSF anisotropy and for a simple of charge trailing between pixels. Hence in practice it appears from these initial studies that on can typically ignore δm.

A.2. Additive Bias
In the additive bias case one can add a field with spin-2 properties such that where |c(Ω)| and φ c (Ω) are systematic changes in the amplitude and rotation angle of the measurements respectively. In this case the two additive components can be expressed like c(Ω) = c 1 (Ω) + ic 2 (Ω) c 1 (Ω) = |c(Ω)|| cos(2Φ(Ω) + φ c (Ω)) c 2 (Ω) = |c(Ω)| sin(2Φ(Ω) + φ c (Ω)), where in the additive case the real and imaginary parts will add to the respective γ 1 and γ 2 parts of the shear and hence we label them as such (which is not the case for the multiplicative biases). In the constant case one can write c 0 = c 1,0 + ic 2,0 where c 1,0 and c 2,0 are constants.

B. THE TWO DIMENSIONAL CASE
The full expanded expression for the two dimensional case can be written in matrix form as The M, N and B terms are defined in equation (8) in the main body of the text.

C. THE TOMOGRAPHIC CASE
In equation (7) we consider the case of a single population of galaxies, however in reality one typically will define several populations labelled as tomographic bins normally delineated in redshift. In this case equation (1) is labelled with a tomographic bin i such that and the power spectra are defined as where i and j label tomographic bins. We use notation where C XY ℓ,ij means that field X is associated with i, and field Y is associated with j. We note that in this case there is a difference between the E B power spectrum and the BE power spectrum for tomographic bins i j. The full expression in expanded form is then: The various matrices in the above expression are where X = (+, −), Y = (+, −), G = (E, B) and H = (E, B). In this case the linearised expressions, assuming no underlying B-modes, or EB power, are This paper was built using the Open Journal of Astrophysics L A T E X template. The OJA is a journal which provides fast and easy peer review for new papers in the astro-ph section of the arXiv, making the reviewing process simpler for authors and referees alike. Learn more at http://beta.theoj.org.