Cosmic Shears Should Not Be Measured In Conventional Ways

Cosmic Shears Should Not Be Measured In Conventional Ways

Jun Zhang, Eiichiro Komatsu

Texas Cosmology Center, University of Texas at Austin, Austin, TX 78712, USA
Department of Astronomy, University of California, Berkeley, CA 94720, USA
E-mail:jzhang@astro.as.utexas.edu
Abstract

A long standing problem in weak lensing is about how to construct cosmic shear estimators from galaxy images. Conventional methods average over a single quantity per galaxy to estimate each shear component. We show that any such shear estimators must reduce to a highly nonlinear form when the galaxy image is described by three parameters (pure ellipse), even in the absence of the point spread function (PSF). In the presence of the PSF, we argue that this class of shear estimators do not likely exist. Alternatively, we propose a new way of measuring the cosmic shear: instead of averaging over a single value from each galaxy, we average over two numbers, and then take the ratio to estimate the shear component. In particular, the two numbers correspond to the numerator and denominators which generate the quadrupole moments of the galaxy image in Fourier space, as proposed in Zhang (2008). This yields a statistically unbiased estimate of the shear component. Consequently, measurements of the n-point spatial correlations of the shear fields should also be modified: one needs to take the ratio of two correlation functions to get the desired, unbiased shear correlation.

keywords:
cosmology: gravitational lensing - methods: data analysis - techniques: image processing: large scale structure
pagerange: Cosmic Shears Should Not Be Measured In Conventional WaysAppendix D – Derivation of the Multiplicative Biases Resulting From a Misuse of Z08pubyear: 2010

1 Introduction

Weak gravitational lensing refers to the weak and systematic shape distortions of background source images (galaxies, CMB, etc.) by the foreground inhomogeneous density distributions on cosmological scales. Since this effect only involves gravity, it has been widely used as a direct probe of matter density fluctuations of our Universe (see, e.g., Hoekstra & Jain 2008 for a recent review).

The weak lensing effect can only be probed statistically due to the fact that the intrinsic projected shape of each galaxy is always somewhat random and anisotropic. A central theme in the study of weak lensing is to find unbiased cosmic shear estimators on galaxy images. This is indeed very challenging because the shape distortion due to weak lensing is generally much weaker than the intrinsic variations of the galaxy shapes.

In the early stage of this field, most of the work focused on issues regarding the use of the quadrupole moments of a galaxy image as a shear estimator (Tyson et al. 1990; Bonnet & Mellier 1995; Kaiser et al. 1995; Luppino & Kaiser 1997; Hoekstra et al. 1998; Rhodes et al. 2000; Kaiser 2000). Ever since then, a number of other shear estimators have been considered in the literature, including moments defined by a certain set of orthogonal functions (Bridle et al. 2001; Bernstein & Jarvis 2002; Refregier & Bacon 2003; Massey & Refregier 2005; Nakajima & Bernstein 2007), the spatial derivatives of the galaxy surface brightness field (Zhang 2008, 2010a), etc..

Conventionally, for each shear component, the shear estimator is simply one number derived from a galaxy image, whose statistical mean is supposed to be equal to the true shear value, provided that the intrinsic galaxy image is statistically isotropic. Unfortunately, even in the absence of the PSF, we show that such shear estimators at least do not exist in simple forms, making them hard to use in practice for precise shear measurements (e.g., in the presence of noise). In the presence of the PSF, we argue that such shear estimators do not likely exist. We give reasons for the above statements in §2. (The readers who are just interested in our new way of measuring shears may skip this section.)

In §3, we present a new form of shear estimators: instead of having only one number from each galaxy image for each shear component, one can keep two numbers, and use the ratio of their averages over many galaxies to accurately measure the cosmic shear. We find that this new way of measuring the shear can be easily implemented by using the method of Zhang (2008) (Z08 hereafter). The new type of shear estimator requires weak lensing statistics such as the n-point correlation functions of the shear field to be carried out in a slightly unusual way, but with little additional cost. This is discussed in §4. In §5, we give numerical examples. Finally, we summarize in §6.

2 Conventional Shear Estimators

Conventional shear estimators are defined as a class of shear estimators, which average over a single quantity per galaxy to estimate each shear component. Most of the existing shear measurement methods belong to this class. For example, the method by Kaiser et al. (1995) and its extensions basically use the quadrupole moments of each Gaussian-Profile-Weighted galaxy image to measure the shear components; in the shapelets method (Refregier 2003), the value of each shear component is estimated from best-fitting a shapelets model to each observed galaxy image; the method of Bernstein and Jarvis (2002) evaluate the shear components by fitting each galaxy shape (also the PSF) with a series of orthogonal 2D Gaussian-based functions (see, e.g., Massey et al. 2007 for more examples). A common feature of these methods is that they all generate one quantity per galaxy for each shear component.

In §2.1, we start our discussion with the well-known examples of shear estimators consisting of quadrupole moments of galaxy images, and show what the issues are. In §2.2, we show that, even in the absence of the PSF, any conventional shear estimator has to reduce to a highly nonlinear form, making it hard to use in practice. In §2.3, we provide arguments as to why we think that conventional shear estimators do not likely exist when a PSF is present.

2.1 A Review of the Problem

The use of galaxy quadrupole moments as shear estimators has been a central topic in weak lensing for many years. It is therefore easier to start our discussion with the quadrupole moments. To present the issues clearly, let us first consider the case without the PSF or any photon noise. For convenience, we use or instead of for coordinates in 2D in this paper.

Suppose that the surface brightness field of the lensed galaxy image is on the image plane, and that of the original (pre-lensing) galaxy image is on the source plane, where and are the position angles on the image and source planes, respectively. We have the following relations:

 fL(→xL)=fS(→xS) (1) →xL=A→xS

where , and , which are the spatial derivatives of the lensing deflection angle. can also be written as , where is sometimes called the lensing potential. Matrix can be alternatively written in terms of the convergence and the two shear components and .

The quadrupole moments of the lensed galaxy image are defined as follows:

 Qij=∫d2→xxixjfL(→x) (2)

where the origin of the coordinates has been chosen to be the center of the light, i.e.,

 ∫d2→x→xfL(→x)=0 (3)

Let us also define the ellipticities of the image as:

 ϵ1 = Q11−Q22Q11+Q22 (4) ϵ2 = 2Q12Q11+Q22

In the absence of the PSF, the quantities and are often thought to be good estimators for and up to the first order in the shear. Let us find out if they are indeed unbiased shear estimators. The observed quadrupole can be rewritten from using eq.(1) in eq.(2):

 Qij = ∫d2→xxixjfS(A−1→x) = |det(A)|∫d2→x(A→x)i(A→x)jfS(→x)

Note that the last step of the above equation is achieved by redefining as . Keeping up to first order in , , and in eq.(2.1), we get:

 Q11−Q22 = (1+4κ)(QS11−QS22)+2γ1(QS11+QS22) (6) Q12 = (1+4κ)QS12+γ2(QS11+QS22) Q11+Q22 = (1+4κ)(QS11+QS22)+2γ1(QS11−QS22) + 4γ2QS12

where are the quadrupole moments of the original galaxy image defined as:

 QSij=∫d2→xxixjfS(→x) (7)

Note that the two light centers defined in the image and source planes coincide. Based on eq.(6), we find:

 ϵ1 = ϵS1+2γ11+2γ1ϵS1+2γ2ϵS2 = ϵS1+2γ1[1−(ϵS1)2]−2γ2ϵS1ϵS2 ϵ2 = ϵS2+2γ21+2γ1ϵS1+2γ2ϵS2 = ϵS2+2γ2[1−(ϵS2)2]−2γ1ϵS1ϵS2

where

 ϵS1 = QS11−QS22QS11+QS22 (9) ϵS2 = 2QS12QS11+QS22

Given that the surface brightness distribution of the original galaxy image is statistically isotropic, we have and . Therefore, we find

 ⟨ϵ1⟩ = 2γ1[1−⟨(ϵS1)2⟩] (10) ⟨ϵ2⟩ = 2γ2[1−⟨(ϵS2)2⟩]

This result, Eq.(10), clearly shows that and are not unbiased shear estimators, as and in the multiplicative factors depend on the galaxy morphology distribution, and cannot be reduced to constant factors. (Also see Eq.(3.29) of Bernstein & Jarvis (2002), or Eq.(9.5.26) of Weinberg (2008).)

One can construct an unbiased estimator of the shear, if one keeps three quantities from each lensed galaxy image: , , and , and use the ratios of their averages. Assuming statistical isotropy of intrinsic galaxy shapes in eq.(6), and keeping up to first order in shear/convergence, we have (also see Eq.(9.5.30) of Weinberg (2008)):

 12⟨Q11−Q22⟩⟨Q11+Q22⟩ = γ1 (11) ⟨Q12⟩⟨Q11+Q22⟩ = γ2

This form of shear estimators is not conventional, as one has to keep more than one quantities from each galaxy image for each shear component. It is this class of estimators we shall discuss in this paper in detail.

One may wonder whether unbiased shear estimators in the conventional form ever exist. The answer is yes, at least when the PSF is absent. For example, we find the following unbiased shear estimators:

 14⟨ln(1+ϵ11−ϵ1)⟩ = γ1 (12) 14⟨ln(1+ϵ21−ϵ2)⟩ = γ2

Eq.(12) can be checked by applying Taylor expansion of () to the first order in shear/convergence using eq.(2.1). Eq.(12) defines a special type of conventional shear estimators that are accidentally found by us. It is now immediately interesting to ask if there exist other types of unbiased shear estimators in the conventional form. We study this issue specifically in the next two sections. If the readers wish to go directly to the relevant sections on the new estimator, read on from §3.

For notational convenience, we shall abbreviate “conventional shear estimator” as “CSE” in the rest of the paper. Once again, by CSE we mean the shear estimators that are made of just one number measured from a galaxy image for each shear component.

2.2 CSE in the Absence of the PSF

In preparation for our main theme of this section, we discuss the spin properties of cosmic shears and their estimators in §2.2.1. We then study the forms of CSEs in the absence of the PSF in §2.2.2.

2.2.1 The Spin of the CSE

To study the forms of the CSEs, it is useful to first consider their properties under coordinate rotations. Suppose we rotate the coordinates clockwise by an angle . The new coordinates are related to the old one via the following relation:

 xθ=xcosθ−ysinθ (13) yθ=xsinθ+ycosθ

If we write the position vector as a complex number of the form , where is the complex unit, the coordinate transformation under rotation can then be written as:

 xθ+iyθ=(x+iy)exp(iθ) (14)

For notational brevity, we shall generally use to denote the value of any quantity in the new coordinates that are rotated clockwise by an angle with respect to the original coordinates.

Let us now discuss how cosmic shears and their CSEs transform under coordinate rotation. The definitions of the shear components (shown in the beginning of §2.1) involve spatial derivatives; thus, their transformation rules under coordinate rotation can be found from the chain rule:

 ∂∂x=∂xθ∂x∂∂xθ+∂yθ∂x∂∂yθ (15) ∂∂y=∂xθ∂y∂∂xθ+∂yθ∂y∂∂yθ

From eq.(13), we get:

 ∂xθ∂x=∂yθ∂y=cosθ (16) ∂yθ∂x=−∂xθ∂y=sinθ

Therefore, we have:

 ∂∂xθ+i∂∂yθ=(∂∂x+i∂∂y)exp(iθ) (17)

Taking the square of eq.(17), we find:

 (∂2∂(xθ)2−∂2∂(yθ)2)+i(2∂2∂xθ∂yθ) = [(∂2∂x2−∂2∂y2)+i(2∂2∂x∂y)]exp(i2θ)

Therefore, the shear components, which are 2nd order derivatives of the lensing potential , also transform under coordinate rotation as:

 γθ1+iγθ2=(γ1+iγ2)exp(i2θ) (19)

Because of this property, we usually say that cosmic shears form a spin-2 quantity. In general, a complex quantity, say , is called a spin- quantity if it transforms as under a clockwise coordinate rotation of angle .

Now let us discuss shear estimators. It is straightforward to see that the CSEs defined in eq.(12) do not form a spin-2 quantity. More generally, assuming that and are the CSEs for and , respectively, then, unlike , is not necessarily a spin-2 quantity. However, it turns out that we can regularize any CSEs by turning them into components of a spin-2 quantity. We have the following lemma:

Lemma 2.1

Based on any pair of CSEs , one can build a new pair of CSEs to form a spin-2 quantity through the following procedure:

 Γ′1+iΓ′2=12π∫2π0dθexp(−i2θ)(Γθ1+iΓθ2) (20)
Proof

Firstly, under a clockwise coordinate rotation by angle , we have:

 Γ′θ01+iΓ′θ02 = 12π∫2π0dθexp(−i2θ)(Γθ0+θ1+iΓθ0+θ2) = 12π∫2π0dθ′exp[−i2(θ′−θ0)](Γθ′1+iΓθ′2) = (Γ′1+iΓ′2)exp(i2θ0)

Therefore, form a spin-2 quantity. To show that are a pair of CSEs, let us take the ensemble average on both sides of eq. (2.2.1):

 ⟨Γ′1+iΓ′2⟩ = 12π∫2π0dθexp(−i2θ)⟨Γθ1+iΓθ2⟩

Since by definition, and measure the shear values in the rotated coordinates, we have,

 ⟨Γθ1⟩=γθ1=γ1cos2θ−γ2sin2θ (23) ⟨Γθ2⟩=γθ2=γ1sin2θ+γ2cos2θ

In a more compact form, we can write eq.(23) as:

 ⟨Γθ1+iΓθ2⟩=(γ1+iγ2)exp(i2θ) (24)

Using eq. (24) in eq. (2.2.1), we get:

 ⟨Γ′1+iΓ′2⟩=γ1+iγ2 (25)

which proves that are indeed also a pair of CSEs.

Due to the invariance under a coordinate rotation of angle , one can always decompose any shear estimator into components of integer spins. Eq.(20) essentially defines a procedure of isolating the spin-2 components of any CSEs using Fourier transformation. Since the cosmic shears form a spin-2 quantity, only the spin-2 components of any CSEs are the relevant/principle components of the estimators. This point is further supported by the fact that the ensemble averages of spin- () components of any CSEs are zero, i.e.,

 12π∫2π0dθexp(−inθ)⟨Γθ1+iΓθ2⟩ = 12π∫2π0dθexp(−inθ)(γ1+iγ2)exp(i2θ) = 0 ( if n≠2 )

Therefore, we only need to focus on spin-2 CSEs from here on.

2.2.2 Spin-2 CSEs

In the weak lensing limit, i.e., when the cosmic shear parameters are small, any spin-2 CSEs can be Taylor expanded to the first order in shear as follows:

 Γ1(γ1,γ2,κ) = (Γ1)0+γ1(∂Γ1∂γ1)0+γ2(∂Γ1∂γ2)0+κ(∂Γ1∂κ)0 Γ2(γ1,γ2,κ) = (Γ2)0+γ1(∂Γ2∂γ1)0+γ2(∂Γ2∂γ2)0+κ(∂Γ2∂κ)0

where ( is any quantity) denotes the value of at . Since , , , and are all spin-2 quantities, their ensemble average must vanish. On the other hand, the coefficients associated with and in eq.(2.2.2) can be decomposed into spin-0 and spin-4 components as follows:

 ⟨Γ1(γ1,γ2,κ)⟩=γ1⟨A+B1⟩+γ2⟨C+B2⟩ (28) ⟨Γ2(γ1,γ2,κ)⟩=γ1⟨B2−C⟩+γ2⟨A−B1⟩

where

 A=12(∂γ1Γ1+∂γ2Γ2) (29) C=12(∂γ2Γ1−∂γ1Γ2) B1=12(∂γ1Γ1−∂γ2Γ2) B2=12(∂γ2Γ1+∂γ1Γ2)

As shown in Appendix A, is a scalar, is a pseudo scalar, is a spin-4 quantity. The ensemble averages of and must vanish. The ensemble average of vanishes if galaxy images have parity symmetry along any direction in the plane of the sky statistically, which is assumed to be true in this paper. Consequently, for to be spin-2 CSEs, we only require . This actually implies that for any individual galaxy, . The reason is that of any single galaxy does not change under coordinate rotation of random angles, and is equal to because the galaxies generated by rotations of a single galaxy form a complete set of statistically isotropic samples (i.e., there are no special directions). As a result, any spin-2 CSEs must satisfy the following necessary condition:

 ∂Γ1∂γ1+∂Γ2∂γ2=2 (30)

This is also a sufficient condition, because directly implies that .

In general, the CSEs are functions of a certain number of shape parameters (e.g., the multipole moments of an image). The functions can be very complicated, and are certainly not fixed by the requirement given by eq.(30). However, for galaxies whose shapes are described by only three parameters (perfect ellipses), any CSE should reduce to a function of just three variables. In this case, we find that any pair of spin-2 CSEs must reduce to a unique form, which is sufficient for us to judge whether CSEs are convenient in practice: namely, if we find that the resulting form is highly non-linear even for such a simple case, then it is reasonable to conclude that CSEs are not so useful for accurate shear measurements from more realistic galaxy shapes as well. This is shown in the rest of this section. For clarity, we refer the readers to Appendix B for the mathematical details/proofs for some of the statements made hereafter in this section.

Let us consider a set of galaxies whose surface brightness profiles can be parametrized as with , where and are the coordinates, is a function of a fixed form, and are the three parameters determining galaxy shapes. For the images to be ellipses, we require the following three things: 1. decays sufficiently fast when becomes large; 2. ; 3. . For example, if ( is the step function) and satisfy the above conditions, the galaxy surface brightness is then distributed evenly inside the ellipse defined by . When such images are weakly lensed, the three conditions are not violated, and becomes without changing the form of . In other words, weak lensing does not introduce additional degrees of freedom to the galaxy shapes. Note that otherwise, one has to consider using more than 3 parameters to construct shear estimators. Among the three parameters, there are indeed only two degrees of freedom useful for shear measurement: the ratios of the parameters. This is because the overall amplitudes of only change the galaxy size, not its shape. As shown in Appendix B, the ellipticities defined in §2.1 are directly equal to , therefore, we can write the shear estimators as functions of only and .

We can further show that must take the following form:

 Γ1+iΓ2=(ϵ1+iϵ2)g(u) (31)

where , and is a one-variable complex function, whose form is to be determined later in this section. To see why eq.(31) is true, one can use the Taylor expansion to write and as power series of and 111We do not consider shear estimators that cannot be Taylor expanded as power series of the galaxy shape parameters in this paper. Such shear estimators likely suffer numerical instabilities in practice.. Consequently, one can write as power series of and , whose spins are and respectively. Since is a spin-2 quantity, each term in the power series must also be a spin-2 quantity. Therefore, in each term of the power series, the power on must be larger than that on by exactly one, i.e., each term must take the form of , where is any non-negative integer, and is a coefficient which can be any complex number at this point. As a result, the shear estimators must have the form defined in eq.(31). To find out what is, let us write it as . Eq. (31) then becomes:

 Γ1=ϵ1g1(u)−ϵ2g2(u) (32) Γ2=ϵ1g2(u)+ϵ2g1(u)

Using the constraint in eq.(30), we find:

 2 = ∂Γ1∂γ1+∂Γ2∂γ2 = ∂Γ1∂ϵ1∂ϵ1∂γ1+∂Γ1∂ϵ2∂ϵ2∂γ1+∂Γ2∂ϵ1∂ϵ1∂γ2+∂Γ2∂ϵ2∂ϵ2∂γ2 = 2(2−u)g1(u)+4u(1−u)dg1du

It is interesting to note that eq.(2.2.2) does not place any constraints on , i.e., it can be any real function. This is because simply adds unnecessary even-parity terms into odd parity ones, and vice versa, without affecting the ensemble averages and the spin of the shear estimators. For convenience, we set hereafter.

Eq.(2.2.2) is a typical first-order ordinary differential equation. It can be solved by introducing an integrating factor which satisfies:

 k(u)(2−u)=ddu[2k(u)u(1−u)] (34)

Multiplying both sides of eq.(2.2.2) with , we get:

 ddu[2k(u)u(1−u)g1(u)]=k(u) (35)

It is now straightforward to solve both eq.(34) and eq.(35). The results are:

 k(u)∝(1−u)−3/2 (36) g1(u)=1u(1+C√1−u)

where is a real number constant. To guarantee that and do not diverge when and approach zero, we need . Finally, we find the unique form for the spin-2 CSEs:

 Γ1+iΓ2=(ϵ1+iϵ2)1−√1−ϵ21−ϵ22ϵ21+ϵ22 (37)

Regarding the uniqueness, it is useful to note that if we transform the CSEs defined in eq.(12) into spin-2 shear estimators using the procedure given in eq.(20), we achieve the same shear estimators as those shown in eq.(37).

We have shown that the principle components (spin-2) of any pair of CSEs have to take specific and highly nonlinear forms for galaxies of elliptical shapes. This feature makes CSEs not convenient in practice (e.g., in the presence of noise).

2.3 CSEs in the Presence of the PSF

Any CSEs which correct for the PSF effect also have to reduce to the forms given in eq.(37) in the limit of zero PSF size when the galaxy images have pure elliptical shapes. For this reason, the conclusion in the previous section is already sufficient to argue against the usefulness of CSEs in practice. For academic interests, we provide the following arguments for why CSEs may not even exist in the presence of the PSF:

In the presence of the point spread function, structural details of galaxy images on scales smaller than the size of the PSF are smeared out. This implies that there are only a finite number of shape parameters (e.g., multipole moments up to some order) available for constructing shear estimators. On the other hand, the derivatives of the lower order shape parameters (e.g., lower order multipole moments) with respect to the cosmic shears depend on the higher order shape parameters because of the PSF, suggesting the requirement for an infinite number of shape parameters to form the shear estimators. Combining the above two reasons, we find it unlikely to form CSEs when a PSF is present. The mathematical details of the above statements are given in Appendix C.

3 A New Way of Estimating Shears

As searching for optimal shear estimators is actively ongoing nowadays (Heymans et al. 2006; Massey et al. 2007; Bridle et al. 2009, 2010), it is important to realize that CSEs (“conventional shear estimators,” by which we mean the shear estimators that are made of just one number from a galaxy image for each shear component) are hard to use in practice due to their unavoidable complex forms even in the absence of the PSF (simpler forms, such as the quadrupole moments, are biased estimators, as shown in § 2.1).

Therefore, existing shear estimators of the conventional type must quantify the bias factor when estimating the shear, which can be achieved numerically (see, e.g., Erben et al. 2001, Bacon et al. 2001, or most recently, Heymans et al. 2006; Massey et al. 2007; Bridle et al. 2010) or estimated analytically (e.g., shear susceptibility in KSB [Kaiser et al. 1995] and derived methods, or responsivity factor in Bernstein & Jarvis 2002 and similar methods), although most people have been mainly focusing on the systematic errors caused by the photon noise and the PSF. However, to achieve percent or even sub-percent level accuracy in cosmic shear measurements, it does not seem enough to completely rely on numerical tests using computer-generated galaxies of limited morphology richness, or approximate analytical methods. Unfortunately, in the presence of PSF, most of the existing shear measurement methods are too complicated or too model-dependent (Voigt & Bridle 2010; Bernstein 2010) to allow for an accurate analytic analysis of the systematic errors in their shear estimators.

The method of Z08 (see also Zhang 2010a for the treatment of photon noise and the pixelation effect) is easily amenable to the corrections described in eq.(10), and can also account for the PSF correction. Not only is it simple, but also well defined regardless of the morphologies of galaxies and the PSF. We show here how to properly use this method (instead of using it as CSEs) to recover the cosmic shear in an unbiased way.

3.1 The Idea

The basic idea of Z08 is to use the spatial derivatives of the galaxy surface brightness field to measure the cosmic shears. It relies on the fact that gravitational lensing does not only distort the overall shape of the object, but also locally modifies the anisotropy of the gradient field of the surface brightness. As it allows for using the shape information from galaxy substructures, the method of Z08 can potentially improve on the signal-to-noise ratio of the shear measurements.

It is shown in Z08 that the shear measurement should be carried out in the Fourier space, in which any PSF can be transformed into the desired isotropic Gaussian form through multiplications, and the spatial derivatives of the surface brightness field can be easily measured. The cosmic shear can be estimated using the following relations:

 12⟨(∂1fO)2−(∂2fO)2⟩⟨(∂1fO)2+(∂2fO)2+Δ⟩=−γ1 (38) ⟨∂1fO∂2fO⟩⟨(∂1fO)2+(∂2fO)2+Δ⟩=−γ2

where

 Δ=β22→∇fO⋅→∇(∇2fO) (39)

is the scale radius of the isotropic Gaussian PSF , which is defined as:

 Wβ(→θ)=12πβ2exp(−|→θ|22β2) (40)

is the surface brightness field. denotes . As shown in Appendix D, the method of Z08 effectively utilizes the quadrupole moments in the Fourier space to measure the cosmic shears.

3.2 A New Unbiased Estimator

Now, here is an important point: in order to implement this method, we must make it clear what we mean by the angular brackets in eq.(38). First, we need to measure the derivatives of the surface brightness and average them within a single galaxy. Let us denote this averaging by , and write:

 12⟨(∂1fO)2−(∂2fO)2⟩g⟨(∂1fO)2+(∂2fO)2+Δ⟩g ⟨∂1fO∂2fO⟩g⟨(∂1fO)2+(∂2fO)2+Δ⟩g

Of course, these are still extremely noisy as they use only one galaxy. The question is then, “how do we average these quantities over many galaxies to obtain an unbiased estimator of the shears?”

If one uses these quantities as if they were the CSEs, then one would simply average them over many galaxies. However, this will produce a biased estimator: 222This bias was not detected in the numerical calculations presented in Z08, as the number of galaxies () used for the average was still too small.

 ⟨12⟨(∂1fO)2−(∂2fO)2⟩g⟨(∂1fO)2+(∂2fO)2+Δ⟩g⟩en=−γ1(1−δ1) (41) ⟨⟨∂1fO∂2fO⟩g⟨(∂1fO)2+(∂2fO)2+Δ⟩g⟩en=−γ2(1−δ2)

where and are the ensemble averages of functions of multipole moments of the galaxy images in Fourier space, and denotes the ensemble average over many galaxies. The derivation of the forms of and is given in Appendix D. and are generally nonzero and dependent on the galaxy morphology.

Instead, we need to take the ensemble averages of the numerator and the denominator separately first, and then divide them to obtain an unbiased estimator:

 12⟨⟨(∂1fO)2−(∂2fO)2⟩g⟩en⟨⟨(∂1fO)2+(∂2fO)2+Δ⟩g⟩en=−γ1 (42) ⟨⟨∂1fO∂2fO⟩g⟩en⟨⟨(∂1fO)2+(∂2fO)2+Δ⟩g⟩en=−γ2

This is the main result of this paper, and the form of the unbiased estimator that we propose to use for the actual analysis of the weak lensing data.

Of course, one could divide the left hand sides of eq.(41) by and to obtain an unbiased estimator. This is similar to correcting the measured shear for a multiplicative bias that is evaluated from the same ensemble of galaxies. In this sense, eq.(41) provides the exact definitions for the multiplicative biases for and . However, since and in eq.(41) involve many high order Fourier-space multipole moments of the surface brightness field, evaluation of these terms from simulations (which are incomplete anyway) can be highly uncertain. Even worse, the multiplicative bias mentioned here is not even a constant, but depends on the morphological distribution of the galaxies. This makes the conventional way of measuring shear correlation functions even more challenging, as one must take into account the correlations of the multiplicative biases, as will be shown in §4.

In summary, according to eq.(42), for each shear component, two quantities from each galaxy should be kept, and the ratios of their ensemble averages yield unbiased estimates for the corresponding shear components. Finally, it is important to note that, to efficiently use eq.(42), the surface brightness of each participating galaxy should be normalized to have roughly the same maximum value, so that faint galaxies are not much less weighted than their brighter counterparts. The details regarding the optimal weighting scheme as a function of the galaxy luminosity should also take into account the photon noise. This is a separate topic, and will be studied in a future work.

3.3 Comments on Errors due to Finite Number of Galaxies

Strictly speaking, Eq.(42) holds when we average over an infinite number of galaxies; however, as we shall show in this section, the error that we make by having a finite number of galaxies for averaging is much smaller than the statistical errors, and thus the estimator remains unbiased for practical applications.

For simplicity, we use eq.(11) rather than our main equation [eq.(42)] in the following discussion, but the conclusion will be the same for eq.(42).

Let us use to denote the average of the quantity over galaxies. From eq.(6), we get:

 ⟨Q11−Q22⟩N = (1+4κ)⟨QS11−QS22⟩N + 2γ1⟨QS11+QS22⟩N ⟨Q12⟩N = (1+4κ)⟨QS12⟩N+γ2⟨QS11+QS22⟩N ⟨Q11+Q22⟩N = (1+4κ)⟨QS11+QS22⟩N + 2γ1⟨QS11−QS22⟩N+4γ2⟨QS12⟩N

Consequently, we have:

 12⟨Q11−Q22⟩N⟨Q11+Q22⟩N = 12Δ1+γ1(1−Δ21)−γ2Δ1Δ2 (44) ⟨Q12⟩N⟨Q11+Q22⟩N = 12Δ2+γ2(1−Δ22)−γ1Δ1Δ2

where

 Δ1 = (45) Δ2 = 2⟨QS12⟩N⟨QS11+QS22⟩N

Here, the terms and contribute to random errors because their ensemble averages vanish, whereas the terms and lead to systematic biases because their ensemble averages do not vanish. Fortunately, as scales as and , the amplitudes of such systematic biases are always much smaller than the sizes of the statistical errors. Therefore, the results from this new type of shear estimators may be regarded as unbiased for practical applications. Numerical verifications will be given in §5.

4 Shear Statistics - n-Point Correlations

The cosmic shear field can only be probed statistically. This is mainly due to the intrinsic variations of the galaxy shapes and the spatial fluctuations of the shear components. As a result, the shear statistics is usually studied in the form of n-point spatial correlation functions of the shear field. The previous discussions and measurements in the literature are based on “conventional” shear estimators (CSEs), i.e., one often assumes that the following is true:

 ⟨Γ⟩en=γ (46)

where can be either or , and is a CSE for . For individual galaxies, eq.(46) implies:

 Γ=γ+Ψ (47)

where satisfies 333Note that for our purpose, it is not necessary to know the form of .. It is usually assumed that ’s of different galaxies do not correlate with each other444This is at least true if the relevant galaxies are separated by a large physical distance. Detailed discussions about the correlations of ’s belong to the topic of “Galaxy Intrinsic Alignment”, which is beyond the scope of this paper.. Therefore, the n-point correlation functions of the shear field can be directly measured by the correlations of ’s.

However, in §2.2, we have shown that such a at least does not exist in a convenient form. Instead, as proposed in §3.2, we can use the new form of shear estimators defined in eq.(42) to probe the cosmic shear in an unbiased way. Let us now find out how to measure the n-point shear correlation functions with the new form of shear estimators. Numerical examples are given in §5.

For notational convenience, the type of shear measurement in eq.(42) can be symbolized as follows:

 ⟨A⟩en⟨B⟩en=γ (48)

where can be either or , and and are properties of a galaxy, such as those defined in eq.(42). Similar to eq.(47), eq.(48) implies the following:

 A=γB+C (49)

where satisfies . If we assume that the of any galaxy does not correlate with the ’s and ’s of other galaxies555Such correlations again belong to the topic of “Galaxy Intrinsic Alignment,” which is not considered in this paper., the n-point correlation functions of the shear field can be probed using the following relation:

 ⟨γ(→x1)γ(→x2)…γ(→xn)⟩en=⟨A(→x1)A(→x2)…A(→xn)⟩en⟨B(→x1)B(→x2)…B(→xn)⟩en (50)

The ensemble averages are taken over a large number of galaxies whose relative positions () are fixed. In practice, the n-point shear correlation functions can be measured using:

 ∑A(→xi1)A(→xi2)…A(→xin)∑B(→xi1)B(→xi2)…B(→xin) (51)

where the sum is taken over all the galaxy groups that satisfy the positional constraints. Note that the ratio is taken after the summations. The standard deviation () of the correlation function in such a measurement can be calculated as follows:

 σ2 ≐ ⟨[∑A(→xi1)A(→xi2)…A(→xin)∑B(→xi1)B(→xi2)…B(→xin)]2⟩en ≐ ⟨∑A2(→xi1)A2(→xi2)…A2(→xin)[∑B(→xi1)B(→xi2)…B(→xin)]2⟩en = N⟨A2(→x1)A2(→x2)…A2(→xn)⟩enN2⟨B(→x1)B(→x2)…B(→xn)⟩2en = 1N(⟨A2⟩en⟨B⟩2en)n

where is the total number of galaxy groups (e.g., the number of galaxy pairs for 2-point correlations) used.

To summarize, in the new type of shear measurement, the shear correlation function should be measured using the ratio of two ensemble averages, as shown in eq.(50). If in eq.(48) is viewed as a multiplicative bias, we need to measure the correlations of these multiplicative biases as well in order to get the correct shear correlation functions.

5 Numerical Tests

In this section, we show how accurately one can recover the cosmic shears and their 2-point correlation functions with the method proposed in §3 and §4. For a comparison, we also show the results using the method of Z08 in the “conventional” (but wrong) way, i.e., eq.(41), but without taking into account the biases, and .666Note that Z08 did not detect these biases because it did not use enough galaxies () in the numerical tests. Here, we use galaxies (or pairs of galaxies for the 2-point correlation function) for each test. Since our focus is to demonstrate principles, we do not include photon noise or the pixelation effect in this paper, and we approximate the PSF as an isotropic Gaussian. Including these effects is straightforward. (Note also that the conventional estimators yield biased results even in this idealized case.)

More comprehensive tests of the method of Z08 have been presented in Zhang (2010b), which further improves the accuracy of shear measurement by including the second order shear/convergence terms. As shown in that paper, the second order correction is proportional to the convergence . We simply set in the numerical simulations here. The residual systematic error on the recovered shear ( or ) therefore should be expected to have a magnitude comparable to the third order terms in shear (e.g., , , , ). The errors shown in the rest of this section are all at the confidence level.

5.1 Image Generation

The mock galaxy images we use in our numerical tests are generated by the algorithms introduced in Z08 and Zhang (2010a), i.e., each galaxy is generated as a collection of point sources. The reason is simple: one can accurately and easily mimic the lensing effect by displacing the points. It also allows us to generate galaxies of complex morphologies. There are two types of galaxies we use in this paper: 1. randomly oriented regular galaxies, each of which contains an exponential disk in the galactic plane (no bulge); 2. irregular galaxies being made of points generated by the trajectories of 2D random walks. For simplicity, the PSF is always an isotropic Gaussian function, whose scale radius is four times the grid size to avoid the pixelation problem. All the lengths in our simulations are in units of the grid size in the rest of this section. The dimension of the grid is .

5.2 1-Point Statistics

As our first example, we study how accurately a single input cosmic shear can be recovered by a large number of mock galaxies, i.e., the 1-point statistics. We use the regular type mock galaxies as introduced in §5.1. Each disk galaxy is composed of ten point sources which are randomly distributed within a radius of 7. The intensity of a point is an exponentially decaying function of its distance to the center of the disk with a decay length equal to 7. The galactic disk is then projected onto the source plane in a random direction. For each input shear value, we use mock galaxies to recover the shear.

To quantify the accuracy of shear recovery, we adopt the standard technique in the weak lensing community by using the “multiplicative bias” and the “additive bias” , which are defined as:

 γmeasured1=(1+m1)γinput1+c1 (53) γmeasured2=(1+m2)γinput2+c2

Our simulations use six sets of input shear values (, ). They are: (, ), (, ), (, ), (, ), (, ), (, ). The recovered shear values as well as the linear fitting results for and