Shear nulling after PSF Gaussianisation: Momentbased weak lensing measurements with subpercent noise bias
Key Words.:
cosmology: observations  gravitational lensing: weak  methods: statisticalAbstract
Context:Current optical imaging surveys for cosmology cover large areas of sky. Exploiting the statistical power of these surveys for weak lensing measurements requires shape measurement methods with subpercent systematic errors.
Aims:We introduce a new weak lensing shear measurement algorithm, shear nulling after PSF Gaussianisation (SNAPG), designed to avoid the noise biases that affect most other methods.
Methods:SNAPG operates on images that have been convolved with a kernel that renders the point spread function (PSF) a circular Gaussian, and uses weighted second moments of the sources. The response of such second moments to a shear of the preseeing galaxy image can be predicted analytically, allowing us to construct a shear nulling scheme that finds the shear parameters for which the observed galaxies are consistent with an unsheared, isotropically oriented population of sources. The inverse of this nulling shear is then an estimate of the gravitational lensing shear.
Results:We identify the uncertainty of the estimated centre of each galaxy as the source of noise bias, and incorporate an approximate estimate of the centroid covariance into the scheme. We test the method on extensive suites of simulated galaxies of increasing complexity, and find that it is capable of shear measurements with multiplicative bias below 0.5 percent.
Conclusions:
1 Introduction
The effect that masses can act as lenses and bend the path of light rays is called gravitational lensing. In the weak lensing regime first considered by Tyson et al. (1990) we statistically measure the slight distortion of the shapes of background galaxies by foreground lenses, called the shear. The subtle effects of weak gravitational lensing on galaxy shapes are an immensely powerful tool in observational astronomy. Amongst other applications, weak lensing has been an invaluable tool for cosmology through measurements of shearshear correlations, called cosmic shear, which are connected to the dark matter power spectrum. After its first detection 15 years ago (Bacon et al. 2000; Van Waerbeke et al. 2001; Wittman et al. 2000; Kaiser et al. 2000) cosmic shear has been extensively used in cosmological studies (e.g. Kilbinger 2015; Hildebrandt et al. 2016; Jarvis et al. 2016).
Currently, large (1000 deg) cosmic shear surveys are ongoing, such as the Kilo Degree Survey (de Jong et al. 2013), the Dark Energy Survey (The Dark Energy Survey Collaboration 2005), and Hyper SuprimeCam (Miyazaki et al. 2012); more hemispheresized missions are planned, such as LSST (Ivezic et al. 2008), WFIRST (Spergel et al. 2015), and Euclid (Laureijs et al. 2011). These surveys will observe unprecedented numbers of galaxies, pushing down statistical errors, and hence requiring percent (for ongoing missions) to subpercent level accuracies (for future missions) on the measured galaxy shapes.
In order to conduct weak lensing studies a crucial point is to measure the shapes of faint background galaxies with high accuracy as well as high precision in the face of inevitable noise, finite image resolution, and pixel effects. The first weak lensing techniques used the moments of the galaxy’s image to estimate its shape and are known as momentbased methods (e.g. Kaiser 1992; Kaiser et al. 1995 (hereafter KSB); Rhodes et al. 2000). These techniques need to use a weighting function with which to cut off the moment integrals so that the moments are not dominated by noise. Having to correct for the effect of the weight function and the PSF convolution are the main challenges for this class of techniques. The widely used KSB method employs an approximate deconvolution scheme, which assumes that the PSF is nearly Gaussian. Newer momentbased methods have improved upon the PSF correction (Melchior et al. 2011), and there have been methods that change the PSF to make the measurement more exact (as explained in Hirata & Seljak 2003 and used by Mandelbaum et al. 2013; Okura & Futamase 2015, 2016).
An alternative class of techniques relies on models of galaxies which are convolved with a PSF and then fit to the galaxy image and are hence known as modelfitting methods (e.g. Kuijken 1999; Miller et al. 2013; Zuntz et al. 2013). These techniques have the benefit of an accurate treatment of the PSF, but in return require realistic models of galaxies. The model of a galaxy is usually a parametric model (e.g. a linear combination of Sersic profiles) and if it does not resemble the intrinsic galaxy, the results can be biased (Bernstein 2010; Voigt & Bridle 2010). A similar class of techniques, known as shapelets methods (Bernstein & Jarvis 2002; Refregier & Bacon 2003; Kuijken 2006), use a set of basis functions which can, in theory, model any galaxy morphology by invoking ever higher order functions. However, in practice the order has to be truncated as the higher functions are dominated by noise, again leading to an unrepresentative galaxy model. In addition, noise in the galaxy image biases all shape measurement methods due to the nonlinear dependence of the galaxy’s ellipticity (the usual description of its shape) on the surface brightness (Refregier et al. 2012; Melchior & Viola 2012; Viola et al. 2014).
In order to quantify these uncertainties and to find ways of calibrating the different techniques, the weak lensing community started shape measurement challenges in which teams competed by using their methods to obtain the most unbiased shear estimate. This started with a general census and benchmark tests in the STEP challenges (Heymans et al. 2006; Massey et al. 2007) and continued with GREAT challenges (Bridle et al. 2010; Kitching et al. 2012; Mandelbaum et al. 2015), which focused on the understanding of different sources of bias. After the most recent GREAT3 challenge it appears that the development in shape measurement algorithms is slowly reaching the goals set by ongoing cosmic shear surveys.
The recent improvement in accuracy was mainly due to the advanced understanding of noise bias. Several authors have introduced correction schemes into their shape measurement methods which are able remove a large portion of the noise bias (Mandelbaum et al. 2015). An alternative route is to avoid biased shear estimators by using estimators with a linear response to the pixel data instead of traditional nonlinear variables, such as the ellipticity. Several authors have used the second moments of the galaxy’s image brightness to estimate the shear (Zhang & Komatsu 2011; Bernstein & Armstrong 2014; Viola et al. 2014). Recently, Bernstein et al. (2016) have reported that their Bayesian method based on moments is able to reach subpercent accuracy even with low signaltonoise () galaxies. However, the drawback of any Bayesian analysis is the requirement of accurate priors, for which external deep observations would be required. This requirement also means that there is no shear estimate for single galaxies, as then knowledge of the intrinsic galaxy profile would be needed, but only a shear estimate for an ensemble of galaxies.
In this paper we propose a novel shape measurement method which may help to reach the ambitious goals of future cosmic shear experiments. Shear nulling after PSF Gaussianisation (SNAPG) is a momentbased method based on a circular Gaussian PSF and weight function, and requires the images to be preprocessed with a PSF Gaussianisation routine. For such galaxies we have an analytic relation between the moments of the galaxy and the shear. Shearing a population of galaxies introduces anisotropy to their ellipticity distribution. Using the analytic expressions, SNAPG reintroduces isotropy to this population by applying a nulling shear to the weighted second moments. The inverse of the nulling shear is then the shear estimate. Such a nulling technique was first advocated by Bernstein & Jarvis (2002). We propose an analytic correction to mitigate the bias due to centroid errors (Bernstein & Jarvis 2002), which is directly computed from the galaxy image.
SNAPG is similar to the Bernstein et al. (2016) method, but instead of a Bayesian framework it uses a nulling technique extract the shear from the second moments of a population of galaxies. It does not require a prior on the intrinsic moments of the galaxy population, but instead relies on the more general requirement that galaxy ellipticities are isotropic. Our novel method thus only produces a shear value for an ensemble of galaxies, but has the benefit that no auxiliary data is needed.
In Sect. 2 we introduce the SNAPG concept and a correction for the bias due to centroid errors. Section 3 describes the image simulations we use to test SNAPG, and in Sect. 4 and Sect. 5 we present the results of the test runs. This is followed by a detailed discussion in Sect. 6, and a summary in Sect. 7.
2 Theory
2.1 Principles of SNAPG
Our novel method combines elements from a number of shear measurement methods. It follows KSB and Luppino & Kaiser (1997) in its use of Gaussianweighted second moments, and uses a nulling technique to estimate the shear (Bernstein & Jarvis 2002). We explain the basics of momentbased methods, such as KSB and Luppino & Kaiser (1997) in Sec. 2.2.
Because the ellipticities used in KSB are nonlinear functions of the pixel values, pixel noise makes them biased estimators. In SNAPG we work with the second moments of galaxies instead, which even in the presence of pixel noise are unbiased estimators as long as the pixel noise in the image is unbiased (as has been previously explored by Zhang & Komatsu 2011 and Viola et al. 2014). For mathematical tractability we require that the PSF in the images is Gaussian and circular; this allows us to work out analytically how the weighted moments respond to any preseeing shear (Sec. 2.3). A similar exercise was done by Rhodes et al. (2000), but here we do not make any simplifying assumptions and find an expression which is valid for all values of the gravitational shear. Note that we do not try to find the intrinsic unweighted moments as in Luppino & Kaiser (1997); instead, we are only interested in the response of the weighted moments to a shear, similar to Bernstein et al. (2016).
The centroid of a source needs to be chosen before the second moments can be calculated. As this position is determined from noisy data, there is a noise dependent shift in the centroid. We show how to incorporate the uncertainty on the centroid of the galaxy into the shear estimator in the approximation where this uncertainty is distributed as a bivariate Gaussian in Sect. 2.4.
The response of the weighted second moments to shear can be used to find the inverse shear which counteracts the gravitational lensing shear. Hence, given the true shear, we can use the inverse shear to compute the second moments of the galaxy before it was lensed. As the true gravitational lensing shear is unknown, we cannot use each galaxy as an independent shear estimator. Rather, we use an ensemble of sheared galaxies as a probe of systematic alignments and calculate the nulling shear that needs to be applied to this ensemble to render their intrinsic ellipticity distribution isotropic. The nulling shear is then the opposite to the true shear affecting these galaxies (Sec. 2.5). Our approach differs from the approach of Viola et al. (2014) and Zhang & Komatsu (2011), who average the numerator and denominator of the ellipticity separately to avoid introducing biases. In SNAPG only the numerator of the ellipticity is used as a measure of the isotropy of the ellipticity distribution, and its response to shear calculated in order to null the signal.
Because real PSFs are not circular Gaussians, SNAPG can only be applied to images that have been convolved with a suitable Gaussianisation kernel (Sec. 2.6). Such a convolution is a linear operation on the pixels, so does not introduce noise bias in the second moments. However, it does correlate the pixel noise, the effect of which can be tracked and corrected for.
2.2 Lensing basics
Here we introduce the basic expressions regarding general shear estimation via the ellipticity of a galaxy as we refer to them often throughout this section. For a more detailed weak lensing review see Bartelmann & Schneider (2001).
A gravitational potential changes the path of light rays moving through it, thereby changing the observed direction of incoming light rays. For extended luminous objects different light rays can be deflected differently and thus we will observe a distorted image of a distant object. To the first order this distortion consists of a stretch (shear) and a magnification (convergence). The deflection angle of light rays from the source depends on the gradient of a suitably defined lensing potential, . The relation between the position of the source and the position of the observed image is known as the lens equation
(1) 
Given that the deflection angles in weak lensing are small, the distortion can be expressed in terms of a Jacobian matrix
(2) 
The parameters
(3)  
are the gravitational lensing convergence and the two components of the shear , , respectively. Without information on the intrinsic size of the lensed source, only the reduced shears can be measured.
The shear affects a galaxy’s polarisation according to
(4) 
where and are the observed and the intrinsic, unlensed polarisation. As the intrinsic shape of a galaxy cannot be measured and the weak lensing shear is very small, the shear has to be statistically obtained from a large number of galaxies experiencing the same distortion. Assuming that galaxies are randomly oriented, the intrinsic polarisations should average out, .
Momentbased methods construct the polarisation of an object from the second moments of image brightness
(5) 
These moments are defined as the noiseless unweighted moments on the intrinsic galaxy image ,
(6) 
with the weight function . However, in practice a galaxy is observed convolved with a PSF ,
(7) 
and the weight function that goes to zero at large is required for the moments not to be dominated by the noise on the image. The aim of momentbased methods is then to estimate the intrinsic polarisations by correcting for the weight function and PSF.
2.3 Effect of preseeing shear on observed Gaussweighted moments
In the case when the PSF is Gaussian, we can reconstruct what the second moments would have been if the galaxy had been sheared.
The weighted second moments of the observed image are
(8)  
with a weight function that depends only on . The order of integration can be swapped and Eq. 8 rewritten as
(9) 
and so relate the weighted second moments directly to the intrinsic galaxy shape as an integral weighted by the expression in square brackets (which depends only on the weight function and PSF).
A gravitationally lensed source has a distorted image: the intrinsic image is transformed to
(10) 
by the distortion matrix .
In order to measure the gravitational shear, we need to know the weighted second moments that we would observe if the galaxy had been distorted by a distortion matrix before PSF convolution; these can be written as
by means of Eq. 10 and the transformation . We now show how the moments of the sheared source can be derived from the observed, PSFconvolved image , by constructing a new weight function which satisfies, for arbitrary ,
(12)  
It is easy to see from Eq. 2.3 that integrating the observed (PSFconvolved but unsheared) image times the weight function will give the moments . Equation 12 shows that can be constructed from the original weight function and the PSF by the following sequence of operations:

Convolving with the PSF;

Distorting the result of the previous step with distortion matrix and divide by ;

Deconvolving the result of the previous step by the PSF.
This recipe is valid as long as the deconvolution in the final step is well defined.
We do not attempt to solve the general problem, but concentrate on the simpler case where both the PSF and the weight function are round Gaussians:
(13) 
and
(14) 
In Appendix A we derive an expression (Eq. 34) for the convolution of with , where and are Gaussians of arbitrary covariance matrices and respectively. By substituting , , and in Eq. 34, we can write the lefthand side of Eq. 12 as
(15)  
The final step is now to deconvolve this expression by the PSF in order to obtain an expression for the weight function that satisfies Eq. 12. We first calculate the result of convolving with a general Gaussian PSF of covariance matrix ; the deconvolution we seek is then obtained by setting (note the sign). The first term (involving ) is straightforward: convolving two Gaussians results in a new Gaussian with covariance matrix equal to the sum. The second term can be calculated using the result of Eq. 34 in Appendix A by setting the matrix defined there to .
After some work we find
(16)  
with
(17) 
The weight function is only useful in practice if it tends to zero at large . This is the case as long as the distortions are small enough so that both eigenvalues of are positive, which is true when
(18) 
As long as the weight function is wider than the PSF (, a reasonable choice if one wants to avoid unnecessarily noisy measurements), this means a useful can be constructed for at least up to 0.3. We show the form of the weight function for a grid of in Fig. 1.
2.4 Bias as a consequence of centroiding errors
Applying the filter in Eq. 16 to an image yields unbiased estimates for the postseeing weighted second moments of the source about as long as the noise on each pixel is unbiased. However, in reality the true centre of the source is unknown, and must be estimated from the image itself. The associated scatter in the centroid biases the second moments. In this section we quantify that bias.
Suppose that is a noisy estimate of the centroid of the observed image . Then, using this centroid our estimate for is
(19) 
If the error distribution of the centroids is then the expectation value of is
(20)  
i.e. the weight function that determines is the original weight function convolved with the centroid error distribution . Hence, conversely, an unbiased estimate of is obtained by using a weight function obtained by deconvolving by the centroid error distribution . We assume that is Gaussian, of covariance . Remembering that expression 16 for was itself obtained by deconvolving Eq. 15 by the PSF, we see that is the deconvolution of Eq. 15 by , i.e. by a Gaussian of covariance matrix . As noted under Eq. 17, this deconvolution is simply accomplished by using Eq. 16 with a modified matrix
(21) 
It remains for us to quantify the covariance matrix of the centroid error for a given source. This will depend on the recipe used to determine the centre.
We centre each source by finding the peak of the correlation of its (noisy) image with a suitable centring kernel , equivalent to finding the optimum positional match between and . The centroid found this way satisfies
(22)  
where we used a Taylor expansion on about , assumed to be the true centre of our source. To derive the noise properties of we first separate into the noisefree observed image and a noise field , and obtain the firstorder relation between and :
(23) 
To calculate the covariance we define as the integral on the lefthand side of Eq. 23, and we assume the backgroundlimited case in which pixel noise is stationary, of constant covariance matrix across a source image. Squaring Eq. 23 and averaging over all possible noise realisations then yields
(24) 
and hence the covariance matrix needed in Eq. 21 is given by
(25) 
Here, depends only on the kernel function and the pixel noise properties , and can be estimated from the noisy image . If the kernel is circular and the noise covariance matrix isotropic , where is the root mean square of the noise background, then becomes a scalar.
A convenient choice is a Gaussian
(26) 
for which
(27) 
and
(28) 
2.5 SNAPG shear nulling estimator
In the previous sections we constructed the filter which, when applied to an observed GaussianPSF smeared image, yields the Gaussweighted second moments that would have been observed (with the same PSF) had the galaxy been distorted by distortion matrix . We have also quantified the noise bias on due to centroiding errors, and constructed a modified filter that compensates for it. In what follows we will drop the “hat” notation and assume that the centroid error correction is applied.
The weight function can be used to construct a shear estimator. If a galaxy is sheared by some known distortion matrix , then we can use the inverse of to find the intrinsic second moments of the galaxy. For a large ensemble of galaxies their combined intrinsic ellipticities (or equivalently their Stokes parameters) average out to zero. Then the search is for the distortion matrix which can null the Stokes parameters of a sheared population of galaxies. The inverse of that distortion matrix is a good estimator of the shear those galaxies experience. To efficiently search for the distortion matrix we use a nulling scheme similar to one already used in shape measurements (Bernstein & Jarvis 2002). In practice, a trial distortion matrix is chosen and the corresponding weight function is computed (see Fig. 1), with which the Stokes parameters for the ensemble of galaxies are calculated. Based on the (an)isotropy of the Stokes parameters a new distortion matrix is chosen, and the previous steps are repeated to reassess the isotropy. This procedure converges in roughly four trials, after which the inverse of the distortion matrix is taken as the shear estimate.
Because galaxies have a wide range of brightness, the Stokes parameters of a galaxy population have a large variance. This translates into a large variance in the nulling shear and increasing precision would require large numbers of galaxies. Alternatively, the moments could be weighted by flux or size to reduce the variance, but this would introduce a bias in the shear. In our current tests such a weight is not required, but we discuss possible solutions for future work in Sec. 6.5.
2.6 PSF Gaussianisation
As indicated in the beginning of Sec. 2.3, SNAPG relies on the assumption of a circular Gaussian PSF. Such a PSF is never present in observational data and thus we need to transform the actual observed PSF into the required PSF. We employ a Gaussianisation process which creates a circular Gaussian PSF by convolving the observed PSF with an appropriate kernel.
Gaussianisation starts by creating a shapelets model of the PSF. Shapelets are a set of basis functions of GaussHermite polynomials, which can be linearly combined to model astronomical objects (Refregier 2003). Convolution in shapelet space is a straightforward procedure, making shapelets an ideal basis for the Gaussianisation process. We use the shapelet implementation of Kuijken (2006) to create a shapelet model of the PSF. In practice, bright stars can be used to obtain a model of the PSF. A best fit circular Gaussian of the shapelets model of the PSF is determined. Then a convolution kernel is found that convolves the PSF into the best fit Gaussian. The resulting kernel is applied to the whole image to create galaxies with circular Gaussian PSFs. See Kuijken et al. (2015) for more detail on the process of PSF Gaussianisation.
It is worth noting that this procedure is different from the one presented by Hirata & Seljak (2003). They assume a Gaussian form for the intrinsic shape of the galaxy when calculating the corrections for PSF nonGaussianity, whereas our procedure is valid for any galaxy morphology. However, it does rely on wellsampled data and was designed with only groundbased PSFs in mind. It is unclear how the procedure would perform for diffractionlimited space telescopes.
The convolution mixes information from neighbouring pixels and hence introduces a correlation between the noise on different pixel values. The resulting noise covariance matrix is given by the original image’s pixel variance, multiplied by the autocorrelation function of the convolution kernel. It is important to propagate this noise covariance into the centroid error estimate (Eq. 24).
3 Image simulations
Set  PSF  Galaxy type  PSF HLR  Galaxy HLR  

Well resolved  Gaussian  Exponential  1.76 pixels  2.5 pixels  5  100 
Barely resolved  Gaussian  Exponential  1.76 pixels  1.5 pixels  5  100 
GREAT08 RNK  Gaussianised Moffat  Exponential or de Vaucouleurs  1.72 pixels  2.1 or 10 pixels  200 
GREAT08 LNK  Gaussianised Moffat  Exponential or de Vaucouleurs  1.72 pixels  2.1 or 10 pixels  20 
To test the performance of SNAPG we create simulated images of galaxies with known applied shear. Following the image simulations of the GREAT challenges, we create a grid of isolated galaxies on postage stamps. This approach gives us a clean test of the performance of SNAPG without introducing errors related to blended galaxy isophotes (see Hoekstra et al. 2015 for a discussion on how blends affect shear measurements). The images of the GREAT challenges do not have circular Gaussian PSFs, so for a clean test of the SNAPG framework we use GalSim (Rowe et al. 2015) to create our own image simulations with perfect circular Gaussian PSFs.
Our simulated galaxy images are a grid of 100 x 100 galaxies, all with exponential profiles of the same size. The grid of postage stamps have a single galaxy randomly offset from the centre of the stamp. The postage stamp is large enough to avoid any bias due to truncation of the surface brightness profile. The half light radius (HLR) of the galaxies is 2.5 pixels. The flux of all galaxies is the same, so that when noise is added all galaxies will have the same . The modulus of the ellipticity of a galaxy is randomly drawn from a Rayleigh distribution of width 0.25, cut off at 0.6 to avoid artificial truncation by the edge of the postage stamp. The position angle of the galaxy is taken from a random uniform distribution between 0 and 180 degrees. The galaxy models are convolved with a Gaussian PSF with a half light radius of 1.76 pixels. The size of the galaxies is larger than the size of the PSF, so we call this set of images the well resolved sample. Each image has a constant shear applied to all 10000 galaxies, where the shear is taken from a grid of for each shear component separately, resulting in 289 different pairs.
We also create a similar set of images where the galaxy half light radius is set to 1.5 pixels. The half light radius of the Gaussian PSF is 1.76 pixels, so that the PSF is larger than the galaxy. We call this set of images the barely resolved sample. Fluxes are fixed and the ellipticity is sampled in the same way as described above. Here too, constant shears are applied to all galaxies on an image, and the shear is taken from the same grid.
These two suites of image simulations contain a total of galaxies. These galaxy images do not contain any noise, instead Gaussian noise is added as required for each test. Each set represents a different target for shape measurement. The well resolved images present our fiducial dataset as the galaxy shapes are not badly affected by the PSF and provide us with a benchmark test of the performance of SNAPG. The barely resolved images present a challenging sample, as galaxy shapes are heavily influenced by pixelisation and severely blurred by the PSF. These galaxies are a difficult target for most shape measurement methods and are sometimes cut from the sample owing to the uncertainty in the galaxy shapes. However, faint small galaxies are abundant in observations and their removal presents a serious loss of statistical power. Having a shear measurement technique able to reliably measure such objects will be a huge advantage for future weak lensing experiments.
As a visual aid to interpreting the different sets of simulated galaxy images we show some of our mock galaxy images in Fig. 2. The upper images show the fiducial well resolved galaxies and the lower images show galaxies from the challenging set of barely resolved galaxies. The images have Gaussian noise added so that the mean 100 (left panel) and 10 (right panel), where is defined as the FLUX/FLUXERR measured by SExtractor on default settings (Bertin & Arnouts 1996). A summary of the image properties can be found in Table 1.
The images of the GREAT08 challenge (Bridle et al. 2010) provide us with a test of the PSF Gaussianisation. In addition, we can compare the performance of SNAPG to other tested methods. We use the 15 LowNoise_Known (LNK) and 300 RealNoise_Known (RNK) sets of images from the challenge, where each image has 10000 isolated galaxies in postage stamps of 40 pixels across. All 10000 galaxies in an image have the same shear applied to provide ample statistics. The galaxies are either an exponential or a de Vaucouleurs profile with a fiducial for LNK and for RNK. The sizes of galaxies are set so that the PSF convolved galaxy size is 1.4 times larger than the PSF size. We use the PSF Gaussianisation algorithm explained in Sec. 2.6 to outfit the GREAT08 images with a circular Gaussian PSF. The Gaussianisation algorithm is applied to PSF set0001, which is a Moffat profile of full width half maximum 2.85 truncated at 5.7 pixels and ellipticity components and . The main properties of the GREAT08 images are summarised in Table 1.
4 Test runs
The new shear nulling method is coded in python and we apply it to the image simulations described in the previous section. The code returns the shear value that nulls the average distortion of each 10,000galaxy image. First we test the SNAPG formalism, then the centroid bias correction formalism, and finally a full implementation of SNAPG. Throughout this section we use =3 and =3 for the widths of the centroid and moment weight functions, respectively. As noted above, none of the images contains noise. Instead, we add noise for each test as storing each noise realisation presented storage problems. For every level of added noise we calculate the mean signaltonoise ratio () using SExtractor with default settings and defined as FLUX/FLUXERR. Each measurement we present was obtained by using the full set of 2.89 galaxies in each set of simulations described in the previous section.
The performance of SNAPG is measured by performing a linear fit using the functional form (Heymans et al. 2006) for each shear component . This procedure quantifies the shear bias as a multiplicative term (e.g. arising from method assumptions or noise) and an additive term (arising from imperfect corrections for the elliptical PSF). Because our simulations are ideal with a circular Gaussian PSF we do not expect any additive bias in these tests.
4.1 High signaltonoise tests of SNAPG
We start by quantifying the performance of SNAPG on the fiducial set of images with well resolved galaxies for =100 for the true centroid of the galaxy. In the left panel of Fig. 3 we show the measured residuals between the input shear and the measured one, and the true input versus the measured shear. We find and for the average of the two components of the shear. We also test the algorithm on the set of barely resolved galaxies and the result is plotted in the right panel. Again SNAPG retrieves the applied shear without detectable bias: and .
As expected SNAPG returns unbiased shear estimates for our ideal images with circular Gaussian PSF with high , showing that the pipeline works. These tests also show the potential of SNAPG as a shear measurement method, regardless of the size of the galaxy in relation to the size of the PSF. We note that at this high the results remain unchanged if the centroid is measured from the data.
4.2 Tests of centroid bias correction
We expect a bias in the shear estimate to originate from the random error on the measured centroid due to image noise. In order to test the effectiveness of the centroid bias correction proposed in Sec. 2.4, we perform a test with truly random centroid values. The well resolved images are reanalysed with SNAPG, but the centroids are artificially offset from the true centroid by a random Gaussian value. The error on the centroid is taken from a normal distribution with a standard deviation of 0.5 pixels. Such a distribution would occur for our simulated images with a for .
Besides introducing random centroid errors, we also add Gaussian noise to the images before measuring the shear. The addition of noise, which is uncorrelated to the centroid error, should not bias SNAPG as the moments are linear with respect to the noisy surface brightness. Hence, we analyse the images with SNAPG several times where each time Gaussian noise with a different root mean square is added to the images. The results of these tests are presented in Fig. 4 as the multiplicative bias in black and in red versus the mean of all galaxies. The dotted lines show the measured bias when no correction is applied. The dashed lines show the bias when the covariance matrix is set to the correct centroid covariance . As we expected, the application of the correction reduces the multiplicative bias from centroid errors to subpercent levels regardless of the noise on the galaxy images. Higher levels of noise increase the variance of the measurements but do not lead to a bias. The mean corrected multiplicative bias over all is and the measured additive bias is below for all .
4.3 Full test of SNAPG
We now add a centroid measurement algorithm to SNAPG and use it as input for SNAPG. The centroids of each galaxy are estimated by nulling the first moments of the galaxy and the centroid error covariance matrix is estimated using Eq. 25. Gaussian noise is added to the well resolved images and SNAPG is run to obtain a shear estimate. In contrast to the previous test, the noise is now directly related to the centroid error. Again we repeat this exercise for different noise levels and show the multiplicative bias as a function of the measured mean in Fig. 5. The dotted line shows the bias in in black and in red without applying the centroid bias correction and the solid line shows the corrected bias. Centroid errors lead to a bias of several percent for , which is decreased to % by the centroid bias correction. The measured additive bias is below 0.1% for all . For very low the estimate for the covariance matrix becomes dominated by noise and the formalism breaks down.
The centroid error correction removes a large part of the noise bias, but does not remove the bias completely. This can occur if the measured covariance matrix is not a good representation of the true centroid variance. To check this hypothesis, we compute the true centroid variance and compare its performance to our previous results. First we measure the centroid from a noisy image and compare this to the true centroid to compute the centroid error . We then estimate the centroid variance as the average over all 10000 galaxies in each image, square the centroid error, and compute the centroid variance as the squared centroid error averaged over all 10000 galaxies in an image . This true centroid error covariance is set into Eq. 21 as the covariance matrix for all 10000 galaxies in the image and the shear for the image is measured. We show these results as dashed lines in Fig. 5 and note that they are very similar to our previous results (solid lines). These findings indicate that the estimate of the centroid variance is good, but that there is unresolved noise bias in SNAPG.
We have also analysed the set of barely resolved images and plotted the multiplicative bias as a function of in Fig. 6. These results sketch a similar picture: noise bias can be reduced by roughly half, but not completely removed. However, even with residual noise bias, SNAPG has only percent level biases for very faint, very small objects. This achievement is remarkable and highlights the potential of SNAPG.
We have corrected for the nonlinearity due to the noisy estimates of the centroid and the second moments themselves have a linear relationship to the noise. Indeed it is shown in Fig. 4 that if the centroid bias is perfectly corrected for, the noise in the second moments does not introduce a bias. However, an additional nonlinearity in SNAPG remains: the correlation between the centroid error and the pixel noise. This correlation is not present in Fig. 4, but it is in Figs. 5 and 6. We now check whether this correlation is the origin of the residual bias after correction, by measuring the centroid and its variance from a different noise realisation than the second moments. We repeat this exercise again for different noise levels and show the results in Fig. 7. Even without the use of the correction (dotted lines) the bias is significantly decreased when compared to Fig. 5, highlighting the bias induced by the correlation. We see that the centroid bias can be corrected to subpercent levels even for galaxies with . Although the correction breaks down for galaxies that are too faint, the bias is consistent with zero down to the lowest .
5 Great08
In the previous section we have shown that the shear estimated using SNAPG is accurate up to the percent level for images with circular Gaussian PSFs. Noise bias can be removed using the centroid correction described in Sec. 2.4 if different noise realisations of the same galaxy are available. However, it remains to be shown how SNAPG would perform on real observations and specifically on images without a circular Gaussian PSF. Therefore, we now apply the SNAPG algorithm to the more realistic Gaussianised images of the GREAT08 competition (see Sec. 3). These tests will show how well the PSF Gaussianisation algorithm performs. As before we use and for the sizes of the weight functions. For the covariance matrix, is given by the autocorrelation function of the kernel which is used to Gaussianise the PSF in the GREAT08 images. This is convolved with the original covariance matrix and this convolution is used in Eq. 24 to compute the centroid covariance matrix.
We measure the shear with SNAPG on the 15 Gaussianised LNK images and the results are shown in the left panel of Fig. 8. There is a slight overestimation of the multiplicative bias of 12%, and there is a small additive bias inconsistent with zero and . The sign of the multiplicative bias and the nonzero additive bias point towards a PSF, which is not a circular Gaussian. The results of SNAPG measurements on the 300 RNK images are shown in the right panel of Fig. 8. Again there is positive residual multiplicative bias , but a slightly lower value than the one we found for the LNK images. This is probably the combination of the Gaussianisation process and the imperfect centroid bias correction we found in the previous section. As we do not possess different noise realisations of the GREAT08 images, we cannot remove the residual noise bias. There is also a small, but statistically significant discrepancy between the bias in and which is not seen in other tests.
We investigate the percent level bias found in the GREAT08 in more detail by looking at the shear bias for various PSF profiles. We simulated two images of galaxies of opposite shears () with a nonGaussian PSF. Six different Moffat profiles with 2, 3, 4, either circular or elliptical, with , were used as PSFs. The PSF half light radius was 1.76 pixels and the galaxy half light radii were 2.5 pixels, so that these images resembled the well resolved images. At the images underwent PSF Gaussianisation and afterwards the shear was estimated. We found that regardless of the original PSF, the bias in the shear is 2%, similar to the results from the GREAT08 images. At such a high this bias is not due to a centroid error and therefore we suspect an imperfect Gaussianisation of the PSF to be the cause. It is unclear which aspect of the PSF Gaussianisation routine causes the bias in the shear estimate, although it does seem to be robust against variations in the PSF profile.
6 Discussion
6.1 SNAPG formalism
We have introduced the SNAPG formalism and tested its performance as a shear measurement method. For galaxy images convolved with a round Gaussian PSF the effect of shear on weighted second moments of image brightness can be analytically calculated. This analytical treatment is used to create a pipeline which finds the gravitational lensing shear by nulling the polarisations for an ensemble of galaxies. This procedure thus finds an estimate for the shear experienced by the galaxies. On test images with high galaxies convolved with a circular Gaussian PSF, the method obtained shear estimates deviating from the input shears by only parts per thousand.
6.2 Noise bias
Like most shape measurement methods, SNAPG suffers a noise bias when applied to images of galaxies with low . However, by using only linear combinations of second moments instead of ratios of moments such as the polarisation or ellipticity, much of the noise bias can be avoided. This strategy allows SNAPG to obtain only a percent level bias in images with a . Noise in the data introduces errors in the centroid estimates, which in turn biases the shear estimates. We compute an analytic treatment to correct the centroids and show that it can significantly improve the performance of SNAPG for low galaxies. Remaining biases after correction for are in the range of less than one percent.
The residual biases increase with decreasing , which indicates that the centroid error correction does not account for the full effect of noise bias. We traced their origin to the correlation between the centroid errors and pixel noise in the second moments. By removing the correlation, we can greatly decrease the measured bias, and also correct for the remaining bias with our centroid bias correction to subpercent accuracy. For multiband surveys a possible solution is to use different filters for the estimates of the centroid and the measurement of the moments. In this way, the correlation between the centroid and the image is removed and without this correlation SNAPG can produce almost unbiased results. The impact on the bias of such a scheme will have to be investigated as galaxy colours and colour gradients may become an issue. In addition, this introduces a correlation with the photometric redshift estimate, which might pose a problem for cosmic shear measurements.
6.3 Galaxy resolution
The shape of a galaxy similar in size to the PSF is heavily distorted by the PSF, making it difficult to estimate the intrinsic shape. However, the analytic treatment of the PSF in the SNAPG formalism ensures that shear estimation is possible even for barely resolved galaxies. For galaxies 0.84 times smaller than the PSF, the shear was retrieved to similar accuracy, as were resolved galaxies for 710 galaxies. By being able to measure unresolved galaxies reliably, SNAPG is able to use the large population of faint small galaxies to boost statistical power.
6.4 PSF Gaussianisation
We have run SNAPG on the images of the ‘LowNoise_Known’ and ‘RealNoise_Known’ branches of the GREAT08 challenge. To make them suitable for SNAPG, the GREAT08 images were first passed through the PSF Gaussianisation. We find a slight overestimation of the shear for the LNK images with , of the order of 12%. The PSF Gaussianisation introduces a correlation in the noise, which is analytically corrected for. SNAPG can retrieve the shear from the RNK images with to an accuracy that is similar to that for the high images. Further tests revealed that this percent level bias is probably inherent to the PSF Gaussianisation routine that we have used. For a variety of PSF profiles the multiplicative shear bias remained constant around 2%. The PSF Gaussianisation appears to be the limiting factor for SNAPG to obtain subpercent shear bias and detailed investigation into this routine is necessary before SNAPG can be reliably applied to observations.
We can compare the performance of SNAPG to the performance of the other methods tested in the GREAT08 challenge. This will only provide an indication as we did not run our pipeline on all datasets in the challenge and shear measurement methods have evolved since. However, a comparison to figures C3 and C4 in Bridle et al. (2010) shows that the 12% bias SNAPG has obtained is at least competitive with other shear measurement methods. A more quantitative comparison to other (recent) shape measurement methods will require testing on image simulations which incorporate realistic observational features. However, optimistically the performance we find for SNAPG is sufficient to meet the requirements of the largest cosmic shear survey to date (Hildebrandt et al. 2016) without any calibration being required.
6.5 Shear precision
So far we have been concerned only with the accuracy of SNAPG, but an equally valid demand is high precision. To estimate the scatter in the shear estimate we use the simulated images of galaxies observed with the Hubble Space Telescope (HST) included in the GalSim software. These galaxies were observed as part of the COSMOS survey (Koekemoer et al. 2007) and we used galaxies between magnitudes 20 and 24.5, similar to the depths of the Kilo Degree Survey and the Dark Energy Survey. These galaxies were rescaled to a pixel size of 0.214 arcseconds and convolved with a circular Gaussian PSF. We find that the scatter in the shear estimate for this set of galaxies is roughly , where is the number of galaxies in the image. Thus the scatter in the SNAPG shear estimate for a fully realistic ensemble of galaxies is worse than an ellipticity based estimate; roughly 23 times more galaxies are needed by SNAPG for the same precision. This result is more optimistic than the increase by a factor of 10 found by Viola et al. (2014) in their analysis of a shear estimator based on Stokes parameters. Our use of a weight function reduces the variation of the moments, thereby shrinking the scatter in the Stokes parameter. In our tests we used identical weights for all sources, which naturally downweighs large, bright galaxies, which would otherwise dominate the ensemble average of second moments. Ideally, in order to optimise the of the individual moment measurement, the size of the weight function should match the observed size of the galaxy. However, fitting weight functions to individual galaxies is in itself a noisy process that may lead to a bias. We therefore advocate using the same weight function size for all galaxies (since most will be only partially resolved, it is not difficult to find a size that is nearly optimal for most of the galaxies by picking a small multiple of the PSF size; see also Eq. 18).
A possible improvement is to assign each galaxy a weight to reduce the variance in the shear estimate. We find that for our sample of HST galaxies weighting by the inverse of the true flux can reduce the required number of galaxies by a factor of 4. This would bring the precision of SNAPG close to the precision of shear estimates based on galaxy ellipticities. In practice, estimating this weight factor from the galaxy fluxes measured in other images (e.g. adjacent photometric bands in a multiband survey) uncorrelated with the lensing images will avoid introducing noise bias.
6.6 Variable shear
Observational weak lensing deals with varying shear fields, for instance in cosmic shear measurements or when measuring the mass of groups or clusters of galaxies. The traditional method is then to average the shear estimate for individual galaxies to obtain the lensing signal. This is not possible with SNAPG as it does not produce a shear estimate per source. In addition, SNAPG requires a large number of galaxies to obtain a precise shear estimate and satisfy the condition that the intrinsic ellipticities average to zero.
Instead of nulling a single shear value for an ensemble of galaxies, we therefore advocate nulling a parametrised model shear field for that ensemble. For example, to measure a galaxygalaxy lensing signal, the model would include parameters that describe the average shear profile of galaxies and their scaling with pertinent galaxy properties. The model parameters would then be varied until the average shear in a number of annular bins around the lensing galaxies is nulled, analogous to a traditional tangential shear stacking analysis. As another example, for cosmic shear measurements, the amplitudes of independent Fourier modes in the shear field could be nulled.
Developing this procedure will be left to the future.
7 Summary
We have presented a new momentbased method that attempts to combine the best aspects of earlier approaches to the problem of highaccuracy, precise shear measurement from galaxy images. Momentbased methods generally approximate the deconvolution of the PSF, but do not require any information beyond the data and generally run very fast. Model fitting methods perform exact forward modelling, including convolution with the PSF, but are expensive to run because they need to search through a large parameter space, and may suffer model bias. The shear nulling after PSF Gaussianisation or SNAPG technique deals analytically with the PSF deconvolution and as a momentbased method only requires a few measurements on the data. In addition, SNAPG incorporates a correction scheme to mitigate the effects of noise bias, a major hurdle to all shape measurement techniques.
Idealised test images show that SNAPG can retrieve shears to percent level accuracy for galaxies with low signaltonoise, even if they are smaller in size than the PSF. The main issue limiting this technique is the correlation between the noisy estimate of the centroid and the pixel noise, which may be mitigated by incorporating further data about the sources, such as images from neighbouring bands in a multiwavelength survey. In such a setup, SNAPG can obtain shear estimates to subpercent accuracy for galaxies with a Gaussian PSF.
Application to real data requires PSF Gaussianisation and if this routine is imperfect it can introduce percent level biases. This level of accuracy is comparable to what is required of the shape measurement algorithms used for ongoing surveys. As such, we expect SNAPG to be a useful asset for current and future weak lensing experiments.
7.1 Acknowledgements
We would like to thank Massimo Viola, Henk Hoekstra, and Peter Schneider for useful discussions. We also thank the anonymous referee, whose suggestions helped to improve this paper. RH acknowledges support from the European Research Council FP7 grant number 279396. AB was supported for this research partly through a stipend from the International Max Planck Research School (IMPRS) for Astronomy and Astrophysics at the Universities of Bonn and Cologne and through funding from the Transregional Collaborative Research Centre ‘The dark Universe’ (TR 33) of the DFG. KK acknowledges support from an Alexander von Humboldt Foundation research award.
References
 Bacon et al. (2000) Bacon, D. J., Refregier, A. R., & Ellis, R. S. 2000, MNRAS, 318, 625
 Bartelmann & Schneider (2001) Bartelmann, M. & Schneider, P. 2001, PhR, 340, 291
 Bernstein (2010) Bernstein, G. M. 2010, MNRAS, 406, 2793
 Bernstein & Armstrong (2014) Bernstein, G. M. & Armstrong, R. 2014, MNRAS, 438, 1880
 Bernstein et al. (2016) Bernstein, G. M., Armstrong, R., Krawiec, C., & March, M. C. 2016, MNRAS, 459, 4467
 Bernstein & Jarvis (2002) Bernstein, G. M. & Jarvis, M. 2002, AJ, 123, 583
 Bertin & Arnouts (1996) Bertin, E. & Arnouts, S. 1996, A&AS, 117, 393
 Bridle et al. (2010) Bridle, S., Balan, S. T., Bethge, M., et al. 2010, MNRAS, 405, 2044
 de Jong et al. (2013) de Jong, J. T. A., Verdoes Kleijn, G. A., Kuijken, K. H., & Valentijn, E. A. 2013, Experimental Astronomy, 35, 25
 Heymans et al. (2006) Heymans, C., Van Waerbeke, L., Bacon, D., et al. 2006, MNRAS, 368, 1323
 Hildebrandt et al. (2016) Hildebrandt, H., Viola, M., Heymans, C., et al. 2016, ArXiv eprints [\eprint[arXiv]1606.05338]
 Hirata & Seljak (2003) Hirata, C. & Seljak, U. 2003, MNRAS, 343, 459
 Hoekstra et al. (2015) Hoekstra, H., Herbonnet, R., Muzzin, A., et al. 2015, MNRAS, 449, 685
 Ivezic et al. (2008) Ivezic, Z., Tyson, J. A., Abel, B., et al. 2008, ArXiv eprints [\eprint[arXiv]0805.2366]
 Jarvis et al. (2016) Jarvis, M., Sheldon, E., Zuntz, J., et al. 2016, MNRAS[\eprint[arXiv]1507.05603]
 Kaiser (1992) Kaiser, N. 1992, ApJ, 388, 272
 Kaiser et al. (1995) Kaiser, N., Squires, G., & Broadhurst, T. 1995, ApJ, 449, 460
 Kaiser et al. (2000) Kaiser, N., Wilson, G., & Luppino, G. A. 2000, arXiv:astroph/0003338 [\eprintastroph/0003338]
 Kilbinger (2015) Kilbinger, M. 2015, Reports on Progress in Physics, 78, 086901
 Kitching et al. (2012) Kitching, T. D., Balan, S. T., Bridle, S., et al. 2012, MNRAS, 423, 3163
 Koekemoer et al. (2007) Koekemoer, A. M., Aussel, H., Calzetti, D., et al. 2007, ApJS, 172, 196
 Kuijken (1999) Kuijken, K. 1999, A&A, 352, 355
 Kuijken (2006) Kuijken, K. 2006, A&Ap, 456, 827
 Kuijken et al. (2015) Kuijken, K., Heymans, C., Hildebrandt, H., et al. 2015, MNRAS, 454, 3500
 Laureijs et al. (2011) Laureijs, R., Amiaux, J., Arduini, S., et al. 2011, ArXiv eprints [\eprint[arXiv]1110.3193]
 Luppino & Kaiser (1997) Luppino, G. A. & Kaiser, N. 1997, ApJ, 475, 20
 Mandelbaum et al. (2015) Mandelbaum, R., Rowe, B., Armstrong, R., et al. 2015, MNRAS, 450, 2963
 Mandelbaum et al. (2013) Mandelbaum, R., Slosar, A., Baldauf, T., et al. 2013, MNRAS, 432, 1544
 Massey et al. (2007) Massey, R., Heymans, C., Bergé, J., et al. 2007, MNRAS, 376, 13
 Melchior & Viola (2012) Melchior, P. & Viola, M. 2012, MNRAS, 424, 2757
 Melchior et al. (2011) Melchior, P., Viola, M., Schäfer, B. M., & Bartelmann, M. 2011, MNRAS, 412, 1552
 Miller et al. (2013) Miller, L., Heymans, C., Kitching, T. D., et al. 2013, MNRAS, 429, 2858
 Miyazaki et al. (2012) Miyazaki, S., Komiyama, Y., Nakaya, H., et al. 2012, in Proc. SPIE, Vol. 8446, Groundbased and Airborne Instrumentation for Astronomy IV, 84460Z
 Okura & Futamase (2015) Okura, Y. & Futamase, T. 2015, ArXiv eprints [\eprint[arXiv]1506.09156]
 Okura & Futamase (2016) Okura, Y. & Futamase, T. 2016, ArXiv eprints [\eprint[arXiv]1606.07837]
 Refregier (2003) Refregier, A. 2003, MNRAS, 338, 35
 Refregier & Bacon (2003) Refregier, A. & Bacon, D. 2003, MNRAS, 338, 48
 Refregier et al. (2012) Refregier, A., Kacprzak, T., Amara, A., Bridle, S., & Rowe, B. 2012, MNRAS, 425, 1951
 Rhodes et al. (2000) Rhodes, J., Refregier, A., & Groth, E. J. 2000, ApJ, 536, 79
 Rowe et al. (2015) Rowe, B. T. P., Jarvis, M., Mandelbaum, R., et al. 2015, Astronomy and Computing, 10, 121
 Spergel et al. (2015) Spergel, D., Gehrels, N., Baltay, C., et al. 2015, ArXiv eprints [\eprint[arXiv]1503.03757]
 The Dark Energy Survey Collaboration (2005) The Dark Energy Survey Collaboration. 2005, ArXiv Astrophysics eprints [\eprintastroph/0510346]
 Tyson et al. (1990) Tyson, J. A., Valdes, F., & Wenk, R. A. 1990, ApJ, 349, L1
 Van Waerbeke et al. (2001) Van Waerbeke, L., Mellier, Y., Radovich, M., et al. 2001, A&A, 374, 757
 Viola et al. (2014) Viola, M., Kitching, T. D., & Joachimi, B. 2014, MNRAS, 439, 1909
 Voigt & Bridle (2010) Voigt, L. M. & Bridle, S. L. 2010, MNRAS, 404, 458
 Wittman et al. (2000) Wittman, D. M., Tyson, J. A., Kirkman, D., Dell’Antonio, I., & Bernstein, G. 2000, Nature, 405, 143
 Zhang & Komatsu (2011) Zhang, J. & Komatsu, E. 2011, MNRAS, 414, 1047
 Zuntz et al. (2013) Zuntz, J., Kacprzak, T., Voigt, L., et al. 2013, MNRAS, 434, 1604
Appendix A Convolution calculations
In this Appendix we calculate the result of convolving with a Gaussian point spread function , where and are Gaussians of arbitrary covariance matrix.
First we consider the product of a noncircular Gaussian of covariance matrix with an offset one of covariance and centre :
(29) 
The sum can be rearranged to yield
(30) 
with
(31) 
and
(32)  
The terms in Eq. 30 not involving simplify to
(33)  
Using this result we can calculate the convolution of with a normalised Gaussian of covariance as
(34)  
where in the last line we have defined . We note that a deconvolution is simply accomplished by changing the sign of .