The ellipticity of galaxy cluster halos from satellite galaxies and weak lensing

TL;DR

This study measures the ellipticity of galaxy cluster halos using satellite galaxy distributions and weak lensing, finding consistent shapes and exploring their implications for dark matter and baryonic physics.

Contribution

It introduces a methodology combining simulations and observations to accurately measure cluster halo ellipticity and compares satellite galaxy shapes with lensing data.

Findings

Mean ellipticity of 0.271 with systematic uncertainties

Satellite galaxy shapes match dark matter halo ellipticity

Misalignment angle between galaxy and halo is about 30 degrees

Abstract

We study the ellipticity of galaxy cluster halos as characterized by the distribution of cluster galaxies and as measured with weak lensing. We use monte-carlo simulations of elliptical cluster density profiles to estimate and correct for Poisson noise bias, edge bias and projection effects. We apply our methodology to 10,428 SDSS clusters identified by the redMaPPer algorithm with richness above 20. We find a mean ellipticity $= 0.271 \pm 0.002$ (stat) $\pm 0.031$ (sys) corresponding to an axis ratio $= 0.573 \pm 0.002$ (stat) $\pm 0.039$ (sys). We compare this ellipticity of the satellites to the halo shape, through a stacked lensing measurement using optimal estimators of the lensing quadrupole based on Clampitt and Jain (2016). We find a best-fit axis ratio of $0.56 \pm 0.09$ (stat) $\pm 0.03$ (sys), consistent with the ellipticity of the satellite distribution. Thus cluster galaxies trace the shape of the dark matter halo to within our estimated uncertainties. Finally, we restack the satellite and lensing ellipticity measurements along the major axis of the cluster central galaxy's light distribution. From the lensing measurements we infer a misalignment angle with an RMS of ${30^\circ \pm 10}^\circ$ when stacking on the central galaxy. We discuss applications of halo shape measurements to test the effects of the baryonic gas and AGN feedback, as well as dark matter and gravity. The major improvements in signal-to-noise expected with the ongoing Dark Energy Survey and future surveys from LSST, Euclid and WFIRST will make halo shapes a useful probe of these effects.

Equations

$$\lambda =\sum _i{p}_{\mathrm{mem},i}=\sum _i\frac{\lambda u({\text{x}}_i|\lambda )}{\lambda u({\text{x}}_i|\lambda )+b({\text{x}}_i)}.$$

$$p_{\mathrm{cen}}({\text{y}}_i)=p_{\mathrm{free}}({\text{y}}_i)\frac{u_{\mathrm{cen}}({\text{y}}_i)}{u_{\mathrm{cen}}({\text{y}}_i)+(\lambda -1)u_{\mathrm{sat}}({\text{y}}_i)+u_{\mathrm{bg}}({\text{y}}_i)}$$

$$P_{\mathrm{cen},i}\propto p_{\mathrm{cen}}({\text{y}}_i)\prod _{j\ne i}(1-p_{\mathrm{cen}}({\text{y}}_i)),$$

$$1=\sum _i{P}_{\mathrm{cen},i}.$$

$$e_1=\frac{I_{xx}-I_{yy}}{I_{xx}+I_{yy}}$$

$$e_2=\frac{2I_{xy}}{I_{xx}+I_{yy}},$$

$$I_{xx}=\frac{{\sum }_i{x}_i^2{w}_i}{{\sum }_i{w}_i},I_{yy}=\frac{{\sum }_i{y}_i^2{w}_i}{{\sum }_i{w}_i},I_{xy}=\frac{{\sum }_i{x}_i{y}_i{w}_i}{{\sum }_i{w}_i},$$

$$e=\sqrt{e_1^2+e_2^2},$$

$$\tan 2\varphi =\frac{2I_{xy}}{I_{xx}-I_{yy}}.$$

$$e=\frac{1-q}{1+q}.$$

$${\Sigma }_{\mathrm{NFW}}\left(\sqrt{x^2+y^2/q^2}\right),$$

$${\Sigma }_{\mathrm{NFW}}(r)={\int }_{-\infty }^{\infty }{\rho }_{\mathrm{NFW}}(\sqrt{r^2+z^2})\mathrm{d}z,$$

$$R_{\lambda }=1\mathrm{Mpc}/h{\left(\frac{\lambda }{100}\right)}^{0.2}.$$

$$e_{\rm true}\pm\delta e_{\rm true}=f(e_{\rm measured}\pm\delta e_{\rm true})\$$

$$\Sigma (R,\theta )$$

$$\equiv$$

$$ϵ\approx -\frac{2{\Sigma }_2(R)}{{\eta }_0(R){\Sigma }_0(R)},$$

$${\gamma }_1$$

$${\gamma }_2$$

$${\Sigma }_{\mathrm{crit}}{\gamma }_1(R,\theta )$$

$${\Sigma }_{\mathrm{crit}}{\gamma }_2(R,\theta )$$

$$I_1(R)$$

$$I_2(R)$$

$${\Sigma }_{\mathrm{crit}}$$

$$\Delta {\Sigma }^{4\theta }$$

$$=$$

$$\Delta {\Sigma }^{\mathrm{const}}$$

$$=$$

$${ϵ}^{\mathrm{observed}}=D\times {ϵ}^{\mathrm{true}}.$$

$$D\equiv \int \mathrm{d}{\theta }_{\mathrm{off}}P({\theta }_{\mathrm{off}})\cos (2{\theta }_{\mathrm{off}}),$$

$$⟨\Delta {\Sigma }^{4\theta }⟩(R_j)=B(R_j)\frac{{\sum }_i{\Sigma }_{\mathrm{crit},ij}\left(\frac{w_{1,ij}{\gamma }_{1,ij}}{\cos 4{\theta }_{ij}}+\frac{w_{2,ij}{\gamma }_{2,ij}}{\sin 4{\theta }_{ij}}\right)}{{\sum }_i{w}_{1,ij}+w_{2,ij}},$$

$$w_{1,ij}={\cos }^2(4{\theta }_{ij}){\Sigma }_{\mathrm{crit},ij}^{-2}/({\sigma }_{\mathrm{shape}}^2+{\sigma }_{\mathrm{meas},ij}^2)$$

$$w_{2,ij}={\sin }^2(4{\theta }_{ij}){\Sigma }_{\mathrm{crit},ij}^{-2}/({\sigma }_{\mathrm{shape}}^2+{\sigma }_{\mathrm{meas},ij}^2),$$

$$⟨\Delta {\Sigma }^{\mathrm{const}}⟩(R_j)=B(R_j)\frac{{\sum }_i(w_{ij}{\Sigma }_{\mathrm{crit},ij}{\gamma }_{1,ij})}{{\sum }_i{w}_{ij}},$$

Full text

††pubyear: 2017††pagerange: The ellipticity of galaxy cluster halos from satellite galaxies and weak lensing–C

1 Introduction

Since dark matter contributes more mass than gas or stars to the halos of galaxies and galaxy clusters, it has a strong influence on the formation of visible galaxies. This means the correlation between dark matter and visible galaxy properties – such as mass, shape, and alignment – can provide constraints on theories of galaxy formation. In practice, though, such correlations are challenging.

However, galaxy clusters, composed of many satellite galaxies orbiting within a common dark matter halo, provide a unique solution. While this cluster halo is again not directly visible, cluster satellites are expected to trace the halo’s shape. Thus, correlations between the brightest cluster galaxy (BCG), or central galaxy (CG), its major axis and the locations of cluster satellites provide valuable information. Early attempts to detect this correlation were made by Kim et al. (2002) with a small cluster sample. The first detection of alignment between the BCG major axis and cluster satellites was made by Brainerd (2005). Previous results were noisy but showed hints of anti-alignment. Later work by Yang et al. (2006) confirmed the major axis alignment of Brainerd (2005). Faltenbacher et al. (2007), Azzaro et al. (2007), and Wang et al. (2008) found further confirmation using the larger SDSS DR4 data set. These studies were based on rather small groups with typically a few members per group. This was followed by Faltenbacher et al. (2008) which used N-body simulations to explain the same measurements.

Agustsson & Brainerd (2006) and Allgood et al. (2006) also used N-body simulations to explain the alignments found in data. Kang et al. (2007) improved on this work with more detailed simulations and mock galaxy catalogs, including using a semi-analytic model of galaxy formation to place galaxies in N-body halos. Kang et al. (2007) found that aligning galaxies perfectly with dark matter halos over-predicted the observed alignment signal. Thus some misalignment would be necessary to explain the data.

Wang et al. (2008) compared measured alignments to Monte Carlo simulations of both triaxial (using the Jing & Suto 2002 model) and isotropic halos. These simulations allowed Wang et al. (2008) to take into account discrete sampling issues due to the small number of satellite galaxies. The data strongly preferred triaxial models which are projected to 2D elliptical models and furthermore allowed Wang et al. (2008) to fit for the halo axis parameters. Niederste-Ostholt et al. (2010) confirmed the significant correlation between BCG and satellite alignment using SDSS clusters with more than 20 members. They also discerned an $\sim 4\sigma$ trend with BCG dominance: clusters with larger differences in brightness between the brightest and 2nd and 3rd brightest members had a stronger BCG-satellite alignment.

Recently, Huang et al. (2016) carried out a methodical study of SDSS redMaPPer clusters. They measured a 35∘ average offset between central galaxy (CG) major axis and major axis of satellite distribution. In addition, they found additional correlations between several cluster properties and alignment. Although these correlations are weak – roughly $0.1$ for central dominance, absolute magnitude, central size, and centering probability – they were detected at $10\sigma$ significance. This is made possible by the large size of the SDSS DR8 data.

The Huang et al. (2016) measurements are state-of-the-art for cluster alignment studies. However, the physical interpretations presented in Huang et al. (2016) still need to assume that cluster satellite galaxies are an unbiased probe of the underlying dark matter. While this may be true qualitatively, effects such as Poisson sampling bias (which is explained in this paper later on) can cause an offset between the dark matter major axis and the major axis of the satellite distribution. This is especially limiting for lower richness clusters since smaller number of satellite galaxies causes the larger offset. Even for richer clusters, where Poisson noise is less problematic, a direct probe of the dark matter would complement the information from visible satellites.

Gravitational lensing provides such a direct probe of the dark matter halos. Early attempts by Hoekstra et al. (2004); Parker et al. (2007) using the first quadrupole estimators were not able to detect the effect. Mandelbaum et al. (2006), with the quadrupole estimators of Natarajan & Refregier (2000), were somewhat more successful. Mandelbaum et al. (2006) also used the larger SDSS data set, probing alignments of galaxies and dark matter with several colour and luminosity bins. They found the strongest alignment for the brightest, red lenses, with a significance $\sim 2-3\sigma$. Other work by van Uitert et al. (2012) and Schrabback et al. (2015) with the same estimator and other surveys found little alignment.

Lensing ellipticity measurements on cluster and group scales have been more successful. Oguri et al. (2010) studied the ellipticity of $\sim 20$ massive clusters using an elliptical NFW fits to individual cluster shear profiles. They reported a $\sim 7\sigma$ detection of cluster ellipticity for a subset of the clusters that were well fit by the elliptical NFW profile. With a larger sample of 4300 clusters, Evans & Bridle (2009) obtained a $\sim 3.5\sigma$ detection of stacked cluster ellipticity. Note that Oguri et al. (2010) has a much higher number density for the source galaxies than Evans & Bridle (2009) does. More recently van Uitert et al. (2016) reported a $3-4\sigma$ detections of halo ellipticity of GAMA groups with KiDS shear catalogs. van Uitert et al. (2016) used several methods for stacking the groups, including aligning with the cluster BCG and the major axis of group members.

Previously, in Clampitt & Jain (2016) we measured the halo ellipticity of Luminous Red Galaxies (LRGs, $4\sigma$ significance) and redMaPPer clusters ($3\sigma$) using new estimators. These estimators involved shears, $\gamma_{1}$ and $\gamma_{2}$, in a Cartesian coordinate system defined by aligning the major axes of all lenses (LRGs) in the stacked measurement. We take the estimator of Clampitt & Jain (2016) one step further and derive the optimal halo ellipticity estimators. The quadrupole naturally has two optimal estimators: one measures a $4\theta$ angular variation in background shears, and the other measures a constant-angle signal. These estimators are subject to different systematics, providing natural checks of our measurement.

In addition, we stack the lensing measurement using the major axis of the CG light profile, as well as the major axis of the satellite distribution. These two ways of stacking are also subject to different systematics, providing further cross-checks. Our results are robust: all four methods (2 estimators and 2 ways of aligning lenses) of measuring halo ellipticity for these redMaPPer clusters are in agreement.

In section 2 we describe the clusters and lensing shear data used in this work. Section 3 presents our methodology and results for SDSS cluster ellipticity of satellite distribution. Section 4 defines the optimal quadrupole lensing estimators and presents the results of the measurement on SDSS clusters. Section 5 compares our results to other work. In section 6 we discuss the results and their implications.

2 Data

2.1 redMaPPer Clusters in SDSS

We use the public SDSS DR8 redMaPPer cluster catalog presented in Rozo et al. (2015) (See also Rykoff et al. 2014 for more details of the algorithm).

Using red galaxies having luminosity $L(z)\geq L_{*}(z)$, where $L_{*}$ is the characteristic luminosity defined in Rykoff et al. (2014) (see their Sec. 4), redMaPPer assigns each cluster a richness, $\lambda$, which quantifies the probable number of cluster satellites. $\lambda$ is given as

[TABLE]

Note that $\lambda$ is the only unknown in the equation, which can be numerically solved. Here, $u(\textbf{x}_{i}|\lambda)$ is the normalized density profile of the cluster and $b(\textbf{x}_{i})$ is the density of the uniform background, while $\textbf{x}_{i}$ is a vector that represents the projected radius, the luminosity and the colour of the $i$-th candidate member galaxy (See Rykoff et al. (2014) for details on how they model $u(\textbf{x}|\lambda)$ & $b(\textbf{x})$ and determine $\lambda$). Note that we round the richness value to the closest integer when performing MC simulations (Sec. 3).

We use all clusters in the range $20<\lambda<200$ and with redshifts between $0.1<z<0.41$. The upper redshift limit is necessary to ensure sufficiently many background sources for lensing measurements. For the $i$-th potential cluster central galaxy (CG), redMaPPer calculates a centering probability, $P_{{\rm cen},i}$ as follows:

[TABLE]

[TABLE]

with a normalization condition

[TABLE]

Note that $P_{{\rm cen}}$ and $p_{{\rm cen}}$ are different from each other. Here, $p_{\rm free}$ is the probability that a galaxy is partially masked by another cluster. $u_{\rm cen}$, $u_{\rm sat}$ and $u_{\rm bg}$ are the normalized distributions of central, satellite and background galaxies, respectively, where $\textbf{y}_{i}$ is a vector containing magnitude, photo-z and local galaxy density of the $i$-th galaxy (see Rykoff et al. (2014) for details on the filters they apply to model $u_{\rm cen}$, $u_{\rm sat}$ and $u_{\rm bg}$).

We only use redMaPPer clusters for which the CG is unambiguous, i.e., clusters with $P_{\rm cen}>0.9$. We have tested and confirmed that the cut on $P_{\rm cen}$ does not alter our fiducial results significantly.

2.2 Galaxy shapes

For lensing measurements we use an SDSS shear catalog with 34.5 million sources Sheldon et al. (2009a). For each source we also have a photometric redshift probability distribution, described in Sheldon et al. (2012), which is necessary to apply the optimal lensing weight to each lens-source pair. The total redshift distribution peaks at $z\sim 0.35$ with a tail towards higher redshifts.

We also use the shear catalog to find the major axes of the CG light distributions. This will be necessary to align our Cartesian coordinate systems for each lens prior to performing a stacked lensing measurement. The shear is estimated with the second moment method with an elliptical gaussian weight that is iteratively readjusted to the size and shape of the object, with Petrosian radius and photometric centroid as the initial guesses for the size and position of the source (see Sheldon et al. (2009a) for details). We match the shear catalog of Sheldon et al. (2009a) to all CG in the redMaPPer catalog. This results in 6,681 matching CGs (64% of CGs have a successful match) in our lens sample. We have confirmed that these 6,681 clusters are an unbiased subsample of the full sample (see Sec. 4.5 for the details of this test). Essentially all of these matched CGs have quality shear measurements, based on the small shape measurement error of each galaxy. We calculate the orientation of CGs and its uncertainty using: $\tan{2\phi}=e_{2}/e_{1}$, where $\phi$ and $e_{i}$ represent the orientation angle and the shear components of the object. The distribution of uncertainty of the orientation angles peaks at zero and 68% (95%) of them are smaller than $1.65^{\circ}$ ($5.75^{\circ}$). See Sheldon et al. (2009a) for more details on the shear catalog, and Sheldon et al. (2004) for more detailed descriptions of the shear measurement method.

Throughout the paper we use a spatially flat $\Lambda$CDM cosmology with $\Omega_{m}=0.3$ and show distances in physical (not comoving) Mpc/$h$.

3 Measurement of Ellipticity of Satellite Distribution

3.1 Raw ellipticity measurement

The contours of redMaPPer satellite density are shown in Fig. 1, with the position re-scaled with the cluster size, $R_{\lambda}$ (see Sec. 3.2.2). The left panel shows the stacked contours with random orientations. The middle panel shows the result when aligning the stack along CG major axes (estimated as in Sec. 2.2). The right panel shows the result of aligning along the satellite major axes (estimated as in Sec. 3.1). The randomly oriented stack has circular contours. Even by eye the nonzero ellipticity of the other two panels is clear. Furthermore, the stack along satellite major axes is visibly more elongated than the stack along CG major axes. Later we will show this is a systematic elongation due to Poisson noise in the satellite distribution (see Sec. 3.2). Such systematics make interpretation of the contour plots difficult.

Given a set of $\lambda$ Cartesian satellite coordinates, e.g., $(x_{i},y_{i})$ for satellite $i$, we can estimate the ellipticity of the cluster’s surface density. We define the two components of ellipticity as

[TABLE]

[TABLE]

where the second moments are given by

[TABLE]

with weights $w_{i}=1/(x_{i}^{2}+y_{i}^{2})$. We define the centre of the coordinates as the centre of each redMaPPer cluster in the catalog. Besides, values of the second moments above are weight-dependent. Therefore, the weight can bias the measured ellipticity from the true ellipticity***The true ellipticity refers to the average ellipticity of the underlying distribution according to which satellites are generated inside $R_{\lambda}$. Note that we assume that the ellipticity is constant with radius inside a cluster in our MC simulation. To confirm that our assumption is accurate enough for the purpose of measuring the average ellipticity, we verified in the appendix B.3 that our MC simulation retrieves the distribution of the measured ellipticity very well.. However, we apply the same weight in our Monte Carlo (MC) simulation (see Sec. 3.2) to properly account for the effect. Therefore, the true retrieved ellipticity would not be biased. The two components of ellipticity can also be written in terms of the total magnitude

[TABLE]

and position angle of the major axis relative to the x-axis

[TABLE]

It will sometimes be useful to convert the magnitude of the ellipticity to an axis ratio $q$, using

[TABLE]

Next we divide the redMaPPer clusters into richness bins and measure the average ellipticity according to Eqs. (5) – (8). In order to apply these equations, we must choose which potential cluster members to include in the measurement. For each cluster and each red galaxy, redMaPPer estimates $p_{\rm mem}$, the probability that the red galaxy is a member of that cluster. This membership probability condenses galaxy colour, luminosity, and projected distance from the cluster centre into a number between 0 and 1 for each red galaxy. In order to minimize the number of possible interlopers (see Sec. 3.3), we need to remove galaxies with small $p_{\rm mem}$ values. However, $p_{\rm mem}$ tends to decrease quickly with increasing projected distance from the cluster centre. Therefore, if we apply a $p_{\rm mem}$ cut too aggressively, most of the member galaxies on the outskirts would be filtered out and our ellipticity estimates would only apply to the cluster centre. To balance these competing effects, we choose an intermediate cut of $p_{\rm mem}>0.5$ for our fiducial results. We have tested other cuts and find that the average ellipticity across all $\lambda$ bins is not sensitive to the precise cut. See Appendix A and Sec. 3.2.2 for a fuller discussion of the effects of $p_{\rm mem}$ cuts.

On the other hand, in redMaPPer algorithm, the projection effect due to nearby halos is likely to be aspherical. For instance, if a projected halo is located to the right of the parent halo, we get an ellipticity enhancement in that direction. However, this ellipticity enhancement is randomly oriented, so when stacking a large number of clusters, we obtain an effectively circularized distribution of the “projected” halos. It does not change the calculation we make in this paper.

In Fig. 2 the black points show the result for the average ellipticity. The error bars on the observed ellipticity are the standard deviation of the mean within each richness bin. Looking at these black points, the ellipticity seems to decrease with richness. However, these raw measurements are subject to several effects that bias the ellipticity of satellite distribution relative to the underlying dark matter halo ellipticity. In the following sections, we describe the use of simulations to correct for three of these potential biases:

• •

“noise bias” from measuring ellipticity with a finite number of satellites,

• •

“edge bias” due to redMaPPer satellites being restricted to lie within a circular aperture, and

• •

bias due to interlopers: red galaxies in the foreground or background that nonetheless make it into the cluster catalog.

3.2 Noise and Edge Bias Corrections

Note that edge and noise bias are connected: noise bias depends on the number of members, and applying a circular edge cut removes some members. For clarity we first describe the effects of noise bias (Sec. 3.2.1) and edge bias (Sec. 3.2.2) in isolation. Then in Sec. 3.2.3 we show the combined effect and the results of correcting both biases in the SDSS data.

3.2.1 Noise Bias

The redMaPPer clusters in our fiducial sample have between 20 and $\sim$ 200 satellites per cluster. Especially on the low-richness end, the small number of satellites results in a bias towards higher ellipticity values: we call this effect “noise bias.” This bias due to measuring ellipticity with a small number of tracers has been noticed in prior work. For example, Niederste-Ostholt et al. (2010) noted their Fig. 2 shows that smaller clusters have a larger ellipticity than larger clusters, a result they interpret as coming from this bias. In addition, Wang et al. (2008) studied groups with as few as four satellites by using MC simulations to reproduce the full ellipticity distribution to data. Although their $\chi^{2}$ values for elliptical halos did not indicate the triaxial model was accurate, these were clearly a much better match to the data than MC isotropic halos.

Like Wang et al. (2008) we use MC simulations of elliptical halos to deal with this Poisson noise bias. We use the simulations to correct the observed mean ellipticity. (Note we ignore the distribution of ellipticity values and assuming each richness bin has just one true ellipticity, but this is checked in Appendix B.3.) We generate simulated clusters with a range of known ellipticity values, then apply our ellipticity estimator. The simulations use an elliptical surface density profile,

[TABLE]

[TABLE]

where $\Sigma_{\rm NFW}(r)$ is the 2D NFW profile at projected radius r (Navarro et al., 1997), $\rho_{\rm NFW}$ is the 3D NFW profile and z is line-of-sight distance. The shape of the profile depends only on mass and concentration. We use the mass-richness relation of Simet et al. (2016). Concentration is a weak function of cluster mass and therefore richness (Mandelbaum et al., 2008a), ranging from $\sim 4.2$ to 5.5 for our clusters (Duffy et al., 2008).

In our MC simulations, we use the same richness bins as in Fig. 2. For a given richness bin, we randomly sample richness values from the SDSS redMaPPer catalog to properly account for fractional abundances of clusters with different richness values. Then for each such cluster, we treat the normalized surface density of Eq. (11) as a probability distribution, and randomly draw $N_{\rm p>0.5}$†††The number of satellites after the $p_{\rm mem}>0.5$ cut. sets of $(x,y)$ coordinates from that distribution. We repeat until the total number of simulated members reaches $\sim 5\times 10^{5}$, which is more than enough for convergence.

The magnitude of measured ellipticity $e$ is estimated from each cluster following Eqs. (5) – (8), then the mean and standard deviation of the mean are calculated. The results are shown in the left panel of Fig. 3 for a range of different richness bins. Comparing the input (true) and output (measured) values of ellipticity, it is clear that noise bias always increases the observed cluster ellipticity. The bias is less problematic for larger $\lambda$ values: for fixed true ellipticity, doubling the number of satellites decreases the bias by about a factor of 2. For all the richness bins, the observed ellipticity asymptotes to the true ellipticity as ellipticity increases. Also, as richness increases, the measured ellipticity is closer to the true, input ellipticity. Note that in the case where we include interlopers (Sec. 3.3 and the left panel of Fig. 4), the uniform interloper distribution biases the measured ellipticity to a lower value. We discuss the interloper correction in Sec. 3.3 in more detail.

Note that correlations between satellites could cause noise bias to be overestimated. Such correlations have been measured by Fang et al. (2016), so we have estimated the magnitude of the effect and confirmed it is negligible. Details of this estimate are shown in Appendix B.2.

3.2.2 Edge Bias

A second bias results from the requirement that redMaPPer galaxies fit within a circular aperture of size given by Rykoff et al. (2014),

[TABLE]

This will cut off potential satellites along the major axis, causing an underestimate of the true ellipticity.

As discussed in Sec. 3.1, a choice of $p_{\rm mem}$ cut can alter the shape of the edge of clusters. If we apply too aggressive a cut on $p_{\rm mem}$, it will remove almost all the member galaxies on the outskirts. In this case, it becomes a soft-edge cut which depends on the value of $p_{\rm mem}$ cut instead of the hard-edge cut at $R=R_{\lambda}$ of redMaPPer ’s radial filter used to select members. However, the soft-edge is more complicated to model and requires an algorithm-level investigation which is beyond our scope. The $p_{\rm mem}>0.5$ cut we apply here preserves enough member galaxies around the edge so that the hard edge-cut of redMaPPer algorithm be still valid.

We can correct for this bias by applying the same $R_{\lambda}$ cut when measuring the correction factor from simulated clusters. Applying this same $R_{\rm\lambda}$ cut to the simulations described in Sec. 3.2.1 leaves $\sim 60\%$ of simulated members. We repeat the same process with this edge cut until the total number of simulated members reaches $\sim 5\times 10^{5}$ per richness bin. This is enough points for convergence in the mean value and error.

For illustration, we first isolate the effect of edge bias by showing only the component of ellipticity ($e_{1}$) which is aligned with the true major axis of the simulated clusters. The result is shown in the centre panel of Fig. 3. As expected, edge bias always tends to decrease the measured ellipticity. There is no significant dependence on richness. The fractional error due to edge bias decreases slightly with input ellipticity.

3.2.3 Results after correcting noise and edge bias

In the right panel of Fig. 3 we show the combined effect of noise and edge bias in the MC simulations. Note that edge and noise bias are connected, since noise bias depends on the number of members and applying a circular edge cut removes some members. In the redMaPPer data, the richness value indicates the number of members after the circular edge cut. So to obtain the corresponding correction from simulations we generate members until the number within the circular aperture $R_{\lambda}$ reaches $N_{\rm p>0.5}$.

This is the correction we apply to the data, resulting in the blue points in Fig. 2. The result is well fit by a constant $0.240\pm 0.002$ (axis ratio 0.613) over an entire decade in richness, from 20 to 200. This range corresponds to a factor of 20 in halo mass, from $0.88\times 10^{14}$ to $1.88\times 10^{15}M_{\odot}/h$. These “corrected” ellipticity errors are propagated via the mapping between observed and true ellipticity in simulations,

[TABLE]

where $f$ represents the mapping from measured ellipticity to true ellipticity from the simulation and $\delta e_{\rm true}\,(\delta e_{\rm measured})$ is the uncertainty in $e_{\rm true}\,(e_{\rm measured})$. Since we assume a constant model for our fit, if there is any relatively strong tendency on ellipticity as a function of richness (see blue dots in Fig. 2), we would underestimate the fitting error by ignoring the non-zero slope of ellipticity as a function of richness. In this sense, we show the constant fitted line just as an index how elliptical the clusters are in average, but it doesn’t mean that clusters in different richness bins would have same ellipticity.

In this section we have studied noise and edge bias using MC simulations. These biases are straightforward to accurately simulate and correct. The next section describes a more difficult problem, the correction of biases due to interlopers or “fake members.”

3.3 Interloper correction

To test the effects of interlopers, we repeat our MC simulations including a population of interlopers along with the cluster members. For a cluster of richness $\lambda$, we take the members having $p_{\rm mem}$ values greater than 0.5 and compute $N_{\rm real}=\sum_{i}p_{{\rm mem},i}$, where the sum is only over those chosen members after $p_{\rm mem}>0.5$ cut.

As before, these real members are drawn from an elliptical NFW distribution. Then we add $N_{\rm fake}=\sum_{i}(1-p_{{\rm mem},i})$ interlopers. The interlopers are drawn from a uniform circular distribution, i.e., a disk with radius $R_{\lambda}$. We use the same richness bins as before, and sample richness and $p_{\rm mem}$ values from the actual redMaPPer catalog for consistency between simulation and measurement.

The resulting correction is shown in Fig. 4. Compared to edge and noise bias alone (Fig. 3) adding bias due to interlopers causes an underestimate of the true ellipticity. We then apply this correction to the raw data (black points) in Fig. 2. The result (green points) has a higher ellipticity, than the raw data or the previous correction without interlopers (blue points). This is expected, since the presence of interlopers dilutes the ellipticity further, requiring an even larger correction than in the case with edge and noise bias alone. The new ellipticity estimate is $0.303\pm 0.003$, about 26% higher than the estimate of 0.240 for which interlopers were neglected.

The difference between cases with and without interlopers is $\Delta e=0.063$. Since we do not know the actual effect of interlopers, we simply take the average of these extreme cases. The true value of ellipticity should lie between these two extreme values, assuming no other source of systematics. The result is an average ellipticity $e=0.271\pm 0.031$. Here we take half the difference as an estimate of the systematic error.

Note that these corrections for the mean ellipticity do not take into account the full distributions of ellipticity within each richness bin, but only the average. However, after obtaining the best-fit ellipticity values with the simulations, we verified that the simulated ellipticity distributions used for the correction accurately match the real distributions.

3.4 Summary

We have corrected raw ellipticity measurement of satellites for noise and edge bias. We also have attempted to correct for interlopers, but this remains our most important systematic. Our best estimate of the ellipticity of the satellite distribution is 0.271. Our estimated systematic error is 0.031, from the difference between extreme cases with no interlopers and $N_{\rm fake}/(N_{\rm fake}+N_{\rm real})$ interloper fraction. The statistical error is only 0.002 (less than 1%) assuming constant ellipticity with richness (within richness bins it is in the 2-5% range). Fig.  2 suggests that the corrected ellipticity is indeed approximately constant.

4 Weak lensing ellipticity

4.1 Method

We use an estimator for halo ellipticity which builds on that developed in Clampitt & Jain (2016). We first summarize the model of Clampitt & Jain (2016) before describing more optimal estimators of ellipticity.

This model for the surface density of elliptical halos uses a multipole expansion:

[TABLE]

and we assume the coeffecient of the quadrupole $-\epsilon\eta_{0}/2\ll 1$, justifying the neglect of higher orders in the expansion. Here, $\eta_{0}$ is logarithmic slope of the monopole, ${\rm d}\ln{\Sigma_{0}}/{\rm d}\ln{R}$. Also, we assume the ellipticity, $\epsilon$, is constant over the whole range of $R$. Here $R$ is the projected distance from the centre of the halo, and $\theta$ is the angle relative to the halo’s major axis. We set

[TABLE]

thus allowing the quadrupole $\Sigma_{2}$ to be completely determined by the monopole $\Sigma_{0}$, up to a proportionality factor $\epsilon$, the magnitude of the ellipticity. In order to avoid confusion with various ellipticity definitions, we actually fit the more universal axis ratio, $q$, given by $\epsilon=(1-q^{2})/(1+q^{2})$.

Recall that the Cartesian components of shear can be derived from the tangential and cross components as

[TABLE]

Clampitt & Jain (2016) builds on Adhikari et al. (2015) to derive the following equations for the quadrupole shear in Cartesian coordinates:

[TABLE]

where

[TABLE]

[TABLE]

where $D_{A}(z_{s}),D_{A}(z_{l})$ and $D_{A}(z_{l},z_{s})$ represent angular diameter distance from observer to source, from observer to lens and from lens to source, respectively. We go a step further than Clampitt & Jain (2016) and define optimally-weighted‡‡‡By “optimal”, we mean that the full (internal plus external) quadrupole signal is captured with a weight that maximizes S/N and further that we have separated the internal and external quadrupole signal. halo ellipticity estimators that include all the information from both Cartesian components:

[TABLE]

and

[TABLE]

which Bernstein & Nakajima (2009) show come from internal and external quadrupoles, respectively. This useful distinction between quadrupoles internal and external to $R$ will be described in greater detail when we present our results in Sec. 4.4.

In the presence of misalignment between the light and dark matter major axes the measured ellipticity is actually an underestimate of the halo ellipticity:

[TABLE]

The misalignment factor is given by

[TABLE]

where $P\,(\theta_{\rm off})$ is the distribution of angle differences $\theta_{\rm off}$ between the major axes of dark matter and the visible proxy (measured major axis of each cluster’s satellite distribution according to Sec. 3.1). See Clampitt & Jain (2016) for the derivation of this equation.

For clarity, we note that we have used three different definitions of major axis in this paper. These are:

• •

Major axis of the underlying distribution according to which member galaxies are scattered (A),

• •

Major axis of distribution of member galaxies for a given cluster, calculated as in Sec. 3.1 (B)

• •

Major axis of dark matter halo (C).

B is subject to a Poisson noise bias with respect to A (see Sec. 3.2.1). Assuming A and C are perfectly aligned, we calculate the misalignment between B and C. For the latter in Appendix C we give the relation between misalignment and richness, $D(\lambda)$, as well as misalignment distributions for several values of ellipticity and $\lambda$. The MC simulations described in Sec. 3 were used to obtain these relationships.

Since we have 10428 different clusters in the catalog, the expected signal would be the average of those clusters. Using Eqs. (25-28) and the misalignment factor described above, we calculate the expected average signal for a range of ellipticity to be fitted to our measurements, adopting the misalignment factors calculated with the ellipticity of the satellite distribution. We discuss our model in more detail in Sec. 4.3.

4.2 Ellipticity estimator

Equation (25) implies the optimal estimator for the varying component of the quadrupole lensing signal is

[TABLE]

where

[TABLE]

where i runs over each lens-source pair inside each j-th radial bin, $R_{j}$ and $\Sigma_{{\rm crit},ij}$ follows the definition of Eq. (24). $B(R)$ is the boost factor accounting for the correlated galaxies that are accidentally included in the lens (see the following paragraph for the definition and details). We have suppressed the $R_{j}$ arguments of the shear components for simplicity. The noise terms of the shear estimator $\gamma$ are the usual shape noise due to intrinsic galaxy ellipticities, $\sigma_{\rm shape}=0.32$, and the shape measurement error, $\sigma_{{\rm meas},ij}\sim 0.05$. Likewise, Eq. (27) implies the optimal estimator for the constant (with angle) quadrupole lensing signal is

[TABLE]

where

[TABLE]

note that $w_{ij}$ here is different from $w_{i}$ in Sec. 3.1.

The boost factor $B(R)>1$ is used to correct the dilution of signal at small scales from sources that are actually physically correlated with the lens. We make use of the random point catalog accompanied by the redMaPPer cluster catalog in Sec. 2.1. Then the boost factor is calculated according to (Sheldon et al. 2004; Mandelbaum et al. 2008a)

[TABLE]

where $i$ indicates lens-source pairs, $k$ indicates random-source pairs and $N_{{\rm lens},j}$ ($N_{{\rm rand},j}$) is the number of lens clusters (random points) inside the j-th radial bin. Note that we see no evidence for azimuthal variation in the boost factor, in agreement with the results of van Uitert et al. (2016). The boost factors for satellite stacking and CG major axis stacking is the same, since the richness and redshift distributions of both samples are nearly identical.

The lensing cross-components are given by

[TABLE]

and

[TABLE]

Note that the quadrupole cross-components do not vanish identically for all lensing like they do for monopole. They only vanish if the mass distribution is symmetric under reflections about the $x$- and $y$-axes. Fig. 1 shows the contours of the satellite distribution are symmetric under reflections. Furthermore, since we are stacking $\sim$ 10,000 clusters (with major axes aligned along the $x$-axis) it is a safe assumption that the mass contours are also symmetric after stacking. The asymmetry of individual cluster halos will be smoothed out in the stack so we expect these cross-components to qualify as null tests.

Note that we have not weighted each lens-source pair by lens ellipticity. It is difficult to obtain a reliable ellipticity for individual cluster halos, due to the biases described in Sec. 3.2. And when stacking using the CG light distribution (as in Sec. 4.5) we see no evidence that weighting by lens ellipticity improves the signal-to-noise ratio (S/N). A similar test was carried out by van Uitert et al. (2016) in their study of group ellipticities from weak lensing. In Section 4.1 of van Uitert et al. (2016) find slightly weaker constraints when weighting by CG ellipticity. They claim this is not surprising since their KSB shape measurement upweights the bulge at the centre of the light. The bulge is rounder and also may not correlate as well as the outer part of the light profile with the dark matter.

For all lensing measurements we obtain the covariance using 200 jackknife patches (Norberg et al., 2009) generated from a kmeans algorithm§§§https://github.com/esheldon/kmeans_radec. This algorithm divides the lens clusters into 200 approximately equal area patches over the entire SDSS survey. We use the full resulting covariance matrix for all fits and quoted S/N values.

4.3 Model

Next we describe the model in more detail. To fit the halo ellipticity we use an NFW model for the monopole and, for a given true ellipticity $\epsilon$, we calculate the quadrupole according to Eqs. (26) and (28). Since an NFW profile depends on the mass of the halo, we need information on the mass of the clusters. Instead of using one effective mass representing the whole cluster sample, we calculate a more accurate model as follows.

For the mass-richness calibration of this cluster sample, we use the results of Simet et al. (2016),

[TABLE]

There are uncertainties in the parameters of the mass-richness relation. However, we only use the mean values of the parameters as shown in the above equation when calculating the mass of each cluster. We have tested that our result is insensitive to the scatter in the mass-richness relation. The corresponding NFW profiles are the monopoles of the clusters. Then using Eqs. (26) and (28) we obtain the quadrupole for internal and external contributions, $\Delta\Sigma^{4\theta}$ and $\Delta\Sigma^{const}$ respectively. Finally, we apply misalignment factor D (Sec. 4.1 and below) to each cluster and average these model quadrupole signals ($D\Delta\Sigma^{4\theta}$ & D$\Delta\Sigma^{const}$) to get the final quadrupole. We generate models with ellipticity varying from 0 to 0.9, assuming the ellipticity is constant with radius. Using these models, we fit the axis ratio instead of ellipticity to the measured quadrupole signal and obtain the final constraint on axis ratio.

The misalignment factor D is calculated as follows. We estimate in Sec. 3 that the true underlying ellipticity of satellite distribution is 0.271. Using the same MC simulation as Sec. 3 assuming $e_{\rm true}=0.271$, we can calculate the misalignment factor D as a function of number of satellite galaxies used in major axis calculation ($N_{sat}$). With an assumption that $D(N_{sat})=D(\lambda)$ we then obtain the misalignment factor $D$ to be applied to every individual cluster. We plot $D(N_{sat})$ in Fig. 12. We propagate the systematic uncertainty $\Delta e_{\rm true}=0.031$ into the model through the calculation of $D$.

Note that our model does not include the 2-halo term which can alter the outer profiles, especially as we approach the maximum radius we measure, $5\;{\rm Mpc}/h$. If the 2-halo term impacts the monopole and quadrupole differently than the 1-halo term, it would impact the ellipticity.¶¶¶The 2-halo term refers to the effect of the surrounding matter around the mass profile of the halo of interest.

4.4 Results

The quadrupole lensing signal is shown in Fig. 5. The measurements of both components $\Delta\Sigma^{4\theta}$ and $\Delta\Sigma^{\rm const}$ are significant over scales $0.1\lesssim R\lesssim 2\;{\rm Mpc}/h$. Both components also have the right sign (negative) for lensing from halo ellipticity. The detection significance,

[TABLE]

is $5\sigma$, with more signal coming from the $\Delta\Sigma^{4\theta}$ component. The corresponding cross-components $\Delta\Sigma^{4\theta}_{\times}$ and $\Delta\Sigma^{\rm const}_{\times}$ are both consistent with zero, with $\chi^{2}$ per degree of freedom being 1.03 and 0.29 respectively. We have also checked that the lensing ellipticity constraints are not very sensitive to the centering probability: using 80% or 95% (instead of our fiducial value of 90%) shifts the best-fit value by well under $1\sigma$.

As expected from Eqs.(26) and (28) the constant term falls off faster than the $4\theta$ term. The $\Delta\Sigma^{\rm const}$ signal is consistent with zero at $\sim 1\;{\rm Mpc}/h$ while the $\Delta\Sigma^{4\theta}$ signal continues out to $\sim 3\;{\rm Mpc}/h$. Physically, this happens because of the nonlocal nature of both components of shear. Recall that for the monopole signal $\Delta\Sigma(R)$ or $\gamma_{+}$ only mass inside $R$ contributes to the signal. Circularly-symmetric (i.e., monopole) mass external to $R$ causes no lensing when averaged over the entire annulus. In contrast, mass arranged in a quadrupole produces two distinct signal components (Bernstein & Nakajima, 2009) depending on whether the quadrupole mass is internal ($\Delta\Sigma^{4\theta}$) or external ($\Delta\Sigma^{\rm const}$) to $R$. Thus, if the shear is produced by an elliptical halo confined to $R\lesssim 1\;{\rm Mpc}/h$, we would expect $\Delta\Sigma^{4\theta}$ to fall off rapidly beyond $R\sim 1\;{\rm Mpc}/h$ since at that point none of the halo is external to $R$. The measurements in Fig. 5 match the predictions for an NFW halo out to $3\;{\rm Mpc}/h$.

On the other hand, $\Delta\Sigma^{\rm const}$ is sensitive to quadrupoles external to the halo as well, perhaps produced by filaments or large-scale structure. In Fig. 5 we also show the best-fit model of the 1-halo term quadrupole. At intermediate scales $\sim 0.5\;{\rm Mpc}/h$ the $\Delta\Sigma^{\rm const}$ signal is slightly stronger than the best-fit model. This excess could include contributions from filaments or other external quadrupoles. Higher S/N data would be needed to conclusively show whether or not our measurement has a contribution from external quadrupoles. However, we can get some intuition for the contribution from Fig. 6 which compares the significance of the joint fit ($5\sigma$) to the significance of $\Delta\Sigma^{4\theta}$ alone ($4.3\sigma$). The two fits are not in tension.

When fitting the axis-ratio we calculate the likelihood $\mathcal{L}$ as

[TABLE]

where $d_{\rm d}$ and $d_{\rm m}$ are vectors of measured data points and model expectation respectively and $Cov^{-1}_{ij}$ is the (i, j) component of the inverse covariance matrix. Note that we combine $\Delta\Sigma^{4\theta}$ and $\Delta\Sigma^{const}$ signal when calculating $\mathcal{L}$. Then we apply the maximum-likelihood method to get the best-fit axis ratio q and its uncertainty, with a non-informative prior on q over [0,1.2]. The full posterior for the axis-ratio fit is shown in Fig. 6. The best-fit value is $q=0.558\pm 0.086$(stat) $\pm 0.026$(sys), or $e=0.284\pm 0.072$(stat) $\pm 0.021$(sys), where the statistical error is dominated by lensing shape noise∥∥∥Jackknife error estimation reflects both shape noise and cluster-to-cluster variance. To confirm that the statistical error is dominated by shape noise, we perform an additional test by randomly rotating the source shapes to wash out cluster-to-cluster variance. The errors on $\Delta\Sigma^{4\theta}$ and $\Delta\Sigma^{\rm const}$ then change by $\sim 5\%$ at most, which indicates our error is dominated by shape noise.and the systematic uncertainty originates from the uncertainty in ellipticity of satellite distribution, which is used to calculate the misalignment factor $D$. This is in an agreement within $1\sigma$ with the axis ratio of the satellite distribution, $0.58$, whose 68% confidence interval is also shown on Fig. 6 for comparison. Note that for the satellite result, we converted the error on ellipticity (0.031) estimated in Sec. 3 to an error on the axis ratio (0.039) via Eq. (10). Note that if we do not correct for Poisson sampling effect of satellites (the misalignment factor D in Eq. 29), the measured best fit axis ratio from lensing is $q=0.649\pm 0.067$, or $e=0.213\pm 0.049$. Comparing the two results with and without Poisson sampling correction, we deduce $D=0.213/0.284=0.750$ with an uncertainty of 0.059, where the uncertainty is evaluated using covariances from the same jackknife patches as in Sec 4.2. This translates into an RMS of $18^{\circ}\pm 2.5^{\circ}$ between the true major axes and the major axes of the observed satellite distribution (the stacking orientation).

We also fit the axis ratio to $\Delta\Sigma^{4\theta}$ and $\Delta\Sigma^{\rm const}$ separately as a consistency test. The result is shown by the blue and the red dashed curves respectively in Fig. 6. The best fit axis ratios of $\Delta\Sigma^{4\theta}$ and $\Delta\Sigma^{\rm const}$ are $q=0.585\pm 0.097$(stat) $\pm 0.024$(sys) and $q=0.365\pm 0.207$(stat) $\pm 0.044$(sys), consistent with our fiducial result.

Since both lensing and satellite axis ratios have been corrected for Poisson sampling, the posteriors in Fig. 6 represent our constraints on the true satellite distribution and the dark matter halo axis ratios, respectively, of redMaPPer  clusters. The agreement of the lensing with the satellite distribution is striking and rules out any significant misalignment between the satellite galaxies and the cluster halo. Likewise, the agreement between our measurements and CDM N-body simulations is excellent. In Sec. 5 we describe those comparisons in more detail. In Table 1 we summarize our main results for lensing and satellite axis ratios.

4.5 Stacking on CG major axis

In Fig. 7 we show the ellipticity of the cluster satellite distribution when stacking along the CG major axis. As in Fig. 2 the result is shown in bins of richness both before and after correction of edge bias. Unlike the case when stacking along the major axis of the satellite distribution itself, in this case the ellipticity of the satellite distribution is not subject to noise bias. Noise bias only arises when the same observable is used to orient the stack and to calculate the ellipticity.

To check that the 6,681 clusters we have with good CG shape measurements are an unbiased subsample, we measure ellipticity of satellite distribution (Sec. 3) using only those 6,681 clusters. We have confirmed that this selection alters the ellipticity of satellite distribution (0.271) at a sub-percent level, much smaller than our statistical error of $\sim 1\%$.

In Fig. 8 we show the lensing signal when stacking on the CG major axis. Again the signal components are consistent with coming from halo ellipticity, with a total significance of $4.4\sigma$. This is slightly less S/N than the satellite-stacked result in Fig. 5. There are two reasons for this: the number of clusters is smaller, increasing the noise, while we have not included the misalignment factor, decreasing the signal. The full posterior of the lensing result is shown in Fig. 9. The best-fit lensing axis ratio (jointly fitting $\Delta\Sigma^{4\theta}$ and $\Delta\Sigma^{\rm const}$) is $q=0.746\pm 0.058$, or $e=0.145\pm 0.038$ (black solid line). Note that the posterior of $\Delta\Sigma^{4\theta}$-only fitting ($0.762\pm 0.064$, blue dashed line) and that of $\Delta\Sigma^{\rm const}$-only fitting ($0.714\pm 0.144$, red dashed line) are consistent to each other and with the joint-fit. Again, in this case the satellite and lensing axis ratios agree well within $1\sigma$ with each other.

Note that these numbers are not corrected for misalignment with the dark matter halo. Our simple MC simulations cannot be used to correct this misalignment – more complex hydrodynamic simulations which include the formation of the stars of the CG would be required. However, we can estimate the misalignment of the CG by calculating how much dilution would be required to match the CG-stacked ellipticity to the satellite-stacked ellipticity. This corresponds to a root-mean-squre (RMS) misalignment angle of $30^{\circ}\pm 10^{\circ}$ between the central galaxy and dark matter halo. This is a slightly larger misalignment than that of satellite galaxies, which have an RMS misalignment of $18^{\circ}\pm 2.5^{\circ}$. Moreover the CG case is a ‘true’ misalignment as the measurement error in its orientation is small, while for satellite galaxies the misalignment is due to Poisson sampling.

5 Comparison with other work

5.1 Simulations

Despali et al. (2017) studies halo shapes in the Le SBARBINE N-body simulations (Despali et al., 2016). They use an Ellipsoidal Overdensity halo finder which is similar to Spherical Overdensity finders except without imposing a spherical shape. Their Figure 11 shows 2D, projected halo axis ratios. Our SDSS sample covers a range of halo masses between $10^{14}$ and $10^{15}M_{\odot}/h$, at $z\sim 0.3$. According to their Figure 11, halos with these masses have axis ratios of 0.65 and 0.55, respectively. This gives an average of approximately 0.6 for our sample. This value is within $1\sigma$ of the lensing result shown in Table 1, $0.558\pm 0.086$. It is also within $1\sigma$ of the satellite result, $0.573\pm 0.039$.

The results of Despali et al. (2017) are consistent with previous N-body simulations, for example, the 2D projected ellipticity plots in Ho et al. (2006). This work found axis ratios of 0.645 (0.621) for halos above $10^{14}M_{\odot}/h$ at $z=0.25$ for cosmologies with $\Omega_{m}=0.3$ and $\sigma_{8}=0.9\,(0.7)$. This axis ratio also agrees with our result at the $\sim 1\sigma$ level. (Note that Ho et al. (2006) used a different ellipticity definition, $q=(1-e)^{2}$. We have obtained the axis ratio numbers via their Figure 1 which shows ellipticities $\sim 0.197\,(.212)$ for these two cases.)

Velliscig et al. (2015) show results for halo axis ratios in hydrodynamical simulations. Results were obtained over a wide range in halo mass, $10^{10}-10^{15}M_{\odot}/h$, by using multiple simulations with different resolutions and sizes. (Small, high-resolution simulations are used to study small halos while large, lower-resolution simulations study rarer, more massive halos.) Reading off their Fig. 3, they find average 3D (not projected) axis ratio $\sim 0.55$ for clusters with mass of $\sim 10^{14}M_{\odot}/h$. We can approximately convert this 3D axis ratio into a 2D axis ratio $\sim$ 0.7 using Fig. 9 of Despali et al. (2017). The inferred 2D axis ratio in Velliscig et al. (2015) is larger at about 2$\sigma$ than our lensing axis ratio $=0.558\pm 0.086$ – further detailed analysis is needed to identify possible sources of this discrepancy.

5.2 Observations

Observational studies of cluster shapes include both the satellite distribution and lensing measurements. A detailed and recent study is that of van Uitert et al. (2016) who measure the ellipticities of 2600 groups ($\lambda>5$) over 180 square degrees from the GAMA group catalog of Robotham et al. (2011). (Note that the van Uitert et al. (2016) groups are smaller than our clusters: only $\sim$ 100-200 of the objects in their sample have richness $\lambda>20$.) They use the KiDS weak lensing catalogs of Hildebrandt et al. (2017) to measure the lensing quadrupole. Stacking on the CG major axis, they find a $3.2\sigma$ detection of lensing ellipticity assuming a diagonal covariance matrix. (Using an alternative estimator that cancels the spurious alignment between lenses and sources, the signal strength is reduced to $2.5\sigma$; see Sec. 4.1 of that work.) Stacking the lensing measurement on the BCG major axis, their best-fit value is $e=0.38\pm 0.12$, corresponding to an axis ratio $=0.45^{+0.14}_{-0.12}$ over 40 kpc $<R<250$ kpc. While this result is consistent with our result, it should not be directly compared with ours because van Uitert et al. (2016) assumes a perfect alignment between the visible light (CG or satellites) and dark matter major axes.

Stacking along the major axis of the satellite galaxy distribution, van Uitert et al. (2016) find $e=0.49\pm 0.13$. In their Eq. (37) they estimate the dilution due to misalignment as $\sim 20-25\%$ when stacking on the satellite major axis. Revising their result to include this misalignment results in $e=0.65\pm 0.17$ or axis ratio $=0.21^{+0.14}_{-0.11}$. This is a $\sim$2$\sigma$ tension with our result, $0.558\pm 0.086{\rm(stat)}\pm 0.026{\rm(sys)}$. It is also a nearly $3\sigma$ tension with the Despali et al. (2017) N-body results, which obtain a 2D axis ratio of $\sim 0.65$ for group-size halos with mass $M\sim 10^{14}M_{\odot}/h$. Halo mass dependence does not resolve this tension as the more massive halos in the N-body results should be more elliptical, not less. Note that their assumption of a diagonal covariance also cannot explain this discrepancy since we do not see any strong off-diagonal component of our covariance.

One possible source for the higher ellipticity is due to the projected galaxies that are not in fact cluster members (for example see Farahi et al. 2016) and impact the lower richness groups more. These “accidentally” projected galaxies are more likely to lie on the outskirts of the lens halos and could bias the major axes toward their position. Then stacking the groups along the calculated major axes would cause a spurious quadrupole lensing detection at large radii. Table 1 in van Uitert et al. (2016) indeed shows that stacking along the major axis of satellite distribution leads to a $\sim 0$ ellipticity at small radii (40 kpc $<R<250$ kpc) but a large ellipticity of $e=0.49\pm 0.13$ for 250 kpc $<R<750$ kpc.

The earlier study of Evans & Bridle (2009) extends to slightly larger groups and clusters in SDSS. Their sample goes down to $\lambda=10$ with a mean richness $\sim 10^{14}M_{\odot}/h$. They find a best-fit axis ratio of $0.48^{+0.14}_{-0.09}$, assuming an NFW model. For the galaxy distribution, they find an axis ratio of 0.6. These numbers are within $1\sigma$ of Despali et al. (2017) and our results, although again cannot be compared directly as they do not include a correction for the Poisson sampling of the satellites.

Huang et al. (2016) found an RMS misalignment of $35.07^{\circ}\pm 0.28^{\circ}$ between the major axis of the satellite distribution and of the CG light. Assuming the misalignment between DM and CG major axes is independent of the misalignment between DM and the satellite distribution, we find an RMS misalignment between the major axis from the satellite distribution and the CG light to be $(30^{2}+18^{2})^{0.5}=35\pm 10$ degrees. Our result is consistent with the result from Huang et al. (2016).

6 Discussion

We have measured the ellipticity of SDSS clusters stacked along the major axes of their satellite distribution and those of the CGs. We have used two different methods (ellipticity of the satellite distribution and ellipticity from quadrupole lensing), and found that the results from the two methods are consistent with each other.

To measure the ellipticity of the satellite distribution, we began by modeling the clusters with MC simulations. The model includes satellites which follow an elliptical NFW profile and interlopers which are uniformly distributed. This model allows us to correct for several biases in the raw ellipticity measurements. First, Poisson noise, due to the finite number of member galaxies in a cluster, causes the raw measured ellipticity to be larger than its true value. Second, the redMaPPer algorithm’s hard edge (at $R=R_{\rm\lambda}$) can decrease measured ellipticity by excluding galaxies along the cluster’s major axis with $R>R_{\rm\lambda}$. We believe we have modeled and corrected these biases sufficiently accurately so that they do not contribute to our systematic uncertainty, but interlopers are more difficult to model. To estimate the error due to interlopers, we considered two extreme cases: one with a large interloper fraction (typically 40%), and one with no interlopers at all. The difference between these cases is about a $\sim 10\%$ uncertainty in ellipticity, which is a conservative estimate of the effect of interlopers. The raw and corrected ellipticity measurements for both cases are shown in Fig. 2. Our final result for cluster ellipticity of the satellite distribution is $0.271\pm 0.002({\rm stat})\pm 0.031({\rm sys})$, which corresponds to an axis ratio of $0.573\pm 0.002({\rm stat})\pm 0.039({\rm sys})$.

In the context of the satellite distribution we have tested a number of other effects. These are generally less important than the biases mentioned above and are summarized in Appendix B. They include ellipticity of the $p_{\rm mem}$ filter (B.1), subhalo clustering due to infall as a group (B.2), ellipticity variation due to intrinsic differences between halos or projection effects (B.3), and richness errors in the redMaPPer algorithm (B.4). In addition, there might be contamination from projected foreground structures which get included in the richness estimation. In this case, it is possible that the fraction of interlopers is larger than assumed here. We try to minimize this effect by selecting member galaxies with $p_{\rm mem}>0.5$ (we checked higher values as well, as discussed in Appendix B). Farahi et al. (2016) estimates the “total” non-member fraction, including both the non-member fraction that is imprinted in $p_{\rm mem}$ values by $\sum{(1-p_{\rm mem})}$ and the non-member fraction that is included in the richness estimation, in SDSS redMaPPer clusters in their Table 2 to be $\sim 0.4$. This is smaller than our average non-member fraction $1-<p_{\rm mem}>=0.473$, therefore, adding more interlopers to account for this projection effect appears unnecessary.

So far, we have assumed that the lensing signal would be dominated by dark matter. However, rich clusters can have up to $\sim 15\%$ of the mass in hot gas, which can have a non-negligible effect on the lensing signal. If the gas is not well aligned with the shape of the dark matter halo, it would alter the ellipticity. Accounting for this effect is out of the scope of this work and we leave it to further studies using the results of simulations that include gas physics. We do compare with one study of AGN feedback below.

The second major aim of this paper is a measurement of the lensing quadrupole of cluster halos. We presented new optimal estimators which allowed a $5\sigma$ measurement of the cluster lensing quadrupole. Aligning the stacked lensing measurement with the satellite major axis, we found a mean axis ratio of $q=0.558\pm 0.086$(stat) $\pm 0.026$(sys), in good agreement with the axis ratio from the satellite distribution itself. This axis ratio, which has already been corrected for Poisson sampling of satellites, represents our constraints on the true dark matter halo axis ratios of redMaPPer clusters. Thus it can be directly compared to cluster-sized halos in N-body simulations – the agreement to the N-body results is within $1\sigma$ (see Sec. 5).

We repeated all satellite and weak lensing axis ratio measurements using the central galaxy (CG) light’s major axis to orient the stack. Again, agreement between the satellite ($q=0.772\pm 0.025$) and weak lensing ($q=0.746\pm 0.058$) axis ratio is within $1\sigma$. Here the lensing detection was slightly weaker at $4.4\sigma$ due to the smaller number of clusters available (about 1/3 of CGs had poor shape measurements) and greater observable-halo misalignment. We infer that the CG is misaligned with respect to the cluster satellite distribution by an RMS angle of $35^{\circ}\pm 10^{\circ}$, consistent with Huang et al. (2016). These two observable methods for orienting the stack are subject to different systematics. For example, the CG and nearby background sources can exhibit spurious alignments (Mandelbaum et al., 2006; Clampitt & Jain, 2016) due to PSF- and camera-based effects. The consistency between the ellipticity from satellite distribution and lensing using both observables to orient the clusters supports the robustness of our measurements.

One application of halo ellipticity measurements is to constrain models of baryonic physics. Our mean axis ratio of $q=0.558$ is in some tension with the no-AGN simulations of Suto et al. (2016). Based on their Figure 9 the integrated probability (DM curve, left panel) from 0 to $q=0.558$ is $\sim$ a few percent. In contrast, the simulations with AGN (right panel) show consistency within $\sim 1\sigma$ with our results. However, it is difficult to say with certainty due to the small number of clusters (40) in the Suto et al. (2016) simulations as well as their somewhat smaller mass ($>5\times 10^{13}M_{\odot}/h$) and the overcooling problem in their no-AGN simulations. Based on Despali et al. (2017) the shift in axis ratio due to the relatively lower mass cutoff is fairly small ($\sim 0.05$). More precise constraints on the models with AGN or other feedback will need to wait for larger simulations. Such simulations are underway by the Eagle (Barnes et al. (2017)) and Illustris (Nelson et al. (2015)) collaborations.

In future work it will be interesting to test for variation in halo ellipticity with radius. The effects of AGN feedback are expected to dominate the inner parts of the cluster, as are those of self-interacting dark matter which makes the inner density contours rounder (see Brinckmann et al. (2017) for a recent study from simulations). It will be a challenge to simultaneously test for these and other effects that go beyond the baseline CDM model, but the joint analysis of the monopole and quadrupole as a function of scale and other parameters should provide useful test. With current data the radial dependence is difficult to test because of the way edge bias interacts with $p_{\rm mem}$ (see Appendix A). Other work, for example Zu et al. (2016), using redMaPPer clusters have also noted difficulties from selection effects produced by cuts on $p_{\rm mem}$. Finally, studies of the alignment of the dark matter halo with large scale structure (LSS) can extend current work comparing light and LSS alignments. For example, Zhang et al. (2013) has shown, using SDSS DR7 data, that galaxies within filaments are aligned with the filament direction. The alignment is strongest for red, central galaxies of groups or clusters. It will be interesting to see whether the alignment between halos and filaments are even stronger. With ongoing surveys such as the Dark Energy Survey (DES Collaboration (2005)) we expect to have at least a factor of two improvement in the statistical signal-to-noise (which dominates for lensing). Future surveys from LSST (Abell et al. (2009); Ivezic el at. (2008)), Euclid (Laureijs et al. (2010)) and WFIRST (Spergel et al. (2013)) will provide huge gains in statistical accuracy and the opportunity to study trends with radius, richness and redshift.

Acknowledgments

We would like to thank Eric Baxter, Joe DeRose, Ian Dell’Anotonio, Daniel Gruen, Hung-Jin Huang, Benjamin Joachimi, Andrey Kravtsov, Rachel Mandelbaum, Surhud More, Yuedong Fang, Mike Jarvis, Masahiro Takada, Joop Schaye, Frank van den Bosch, Mark Vogelsberger, Anja von der Linden and Risa Wechsler for helpful discussions. We are very grateful to Erin Sheldon for allowing us to use his SDSS shear catalogs and source photometric redshift distributions. JC, BJ and GB are partially supported by the US Department of Energy grant DE-SC0007901 and funds from the University of Pennsylvania.

Appendix A Effect of Different $p_{\rm mem}$ Cuts

For the edge bias correction (Sec.  3.2.2) to be accurate, we need sufficiently high galaxy number density around the edges, $R\sim R_{\rm\lambda}$, of redMaPPer clusters. With too aggressive a $p_{\rm mem}$ cut, the density profile contracts, falling off faster than NFW. It then fails to match our MC simulation with elliptical NFW out to $R=R_{\rm\lambda}$. In Fig. 10, we compare the stacked number density profiles with different $p_{\rm mem}$ cuts. For example, implementing $p_{\rm mem}>0.8$ leaves 10-50 times fewer member galaxies on the outskirts compared to the case with no $p_{\rm mem}$ cut at all. Considering the total number of member galaxies per cluster (order of $10-100$) and the decreasing tendency of $p_{\rm mem}$ with projected distance from the cluster centre, this suggests that there will be almost zero member galaxy left on the outskirts. On the other hand, a $p_{\rm mem}>0.5$ cut leaves about 25-60% of member galaxies in the outermost region. This better approximates a sharp boundary at $R=R_{\rm\lambda}$ as in our MC simulations. Therefore, we decide to use the $p_{\rm mem}>0.5$ cut as our fiducial case, since further relaxing the cut to a lower value of $p_{\rm mem}$ causes higher systematic error in ellipticity due to increased interlopers. For instance, if we do not apply any $p_{\rm mem}$ cut at all, the uncertainty on ellipticity increases by $\sim 150\%$, compared to the $p_{\rm mem}$>0.5 case, while the value of ellipticity stays at $\sim 0.27$.

Furthermore, the corrected value of ellipticity is not very sensitive to the choice of $p_{\rm mem}$ cut. We have checked that using different criteria for the $p_{\rm mem}$ cut makes less than a 5% difference on the measured ellipticity value.

These uncertainties due to interlopers and $p_{\rm mem}$ cuts are the most significant, so we have emphasized them here in the main text. In Appendix B we describe a few other potential sources of error that turned out to be negligible.

Appendix B Miscellaneous tests

B.1 Elliptical $p_{\rm mem}$ distribution

The value of $p_{\rm mem}$ strongly depends on the projected distance from the cluster centre. However, it could also have an angular dependence. In this case, the soft edge that a $p_{\rm mem}$ cut imposes becomes elliptical, engendering totally a different edge bias. Because the satellite distribution is typically misaligned with the underlying dark matter distribution, we cannot figure out how elliptical this additinoal soft edge (by a $p_{\rm mem}$ cut) of the stacked clusters would be, unless accurate an N-body simulation is provided. We calculate the monopole signal of this soft edge that our $p_{\rm mem}$ cut imposes by dividing the average of $p_{\rm mem}>0.5$ by that without any $p_{\rm mem}$ cut, as a function of projected distance from the cluster centre with an appropriate normalization. This monopole can be used to construct two extreme cases to estimate the possible systematic due to this uncertainty in the ellipticity of $p_{\rm mem}$ distribution : 1) perfectly circular edge and 2) elliptical edge with a specific ellipticity we want to apply. We can apply these two extreme cases of soft edges to our MC simulation and measure the true ellipticity of each case. For example, we measure the true ellipticity of 0.40 assuming a circular edge, whereas we measure the true ellipticity of 0.18 assuming an elliptical edge with ellipticity of 0.3.

B.2 Satellite infall in groups or filaments

We know that these redMaPPer members are clustered with each other (Fang et al., 2016). Our simulations don’t include this clustering of red galaxies, so they will slightly overestimate ellipticity by undercorrecting for noise bias. We estimate the size of this effect based on the measurements of Fang et al. (2016). Using the definition of the angular correlation function $w$, for an annulus at a projected distance $R$ with area $dA$, the number of satellites around a given satellite is ${\rm d}N(R)={\rm d}A\times\bar{n}_{nfw}(R)\times[1+w(R)]$, where $\bar{n}_{nfw}(R)$ denotes a smooth, unclumped NFW profile. Therefore, the excess due to correlations at a projected distance $R$ is just $w(R)$. If we integrate that over the cluster, centered on a typical satellite at $\sim 0.3\;{\rm Mpc}/h$ from the cluster centre, then we get the fractional excess in counts (and hence overestimate of noise bias) due to satellite clustering.

We can estimate this quantity using Figure 1 of Fang et al. (2016). The ratio of satellite correlations (red points) to the correlations associated with the host cluster (black points) ranges from 0.2 to 0 over $R=0.1-1\;{\rm Mpc}/h$. Since the average $\bar{w}(R)\sim\int R\,{\rm d}R\,w(R)$ is about 0.1, we would need to correct noise bias by a factor of 1.1. In other words, satellite correlations mean that there are 10 percent fewer independent points. Note that Fang et al. (2016) showed cross-correlations of cluster members with redMaGiC galaxies, while we are interested in the auto-correlations of cluster members. However, the results are the same: we have checked that the ratio of auto-correlations of cluster members to cluster-member correlations is the same as the analogous ratio shown in Figure 1 of Fang et al. (2016).

Next, we compare this estimate of 10% fewer independent points to the noise bias correction in Fig. 4. From Fig. 4 we estimate the effect of 10% fewer independent points is negligible. Noise bias varies most rapidly for smaller $\lambda$, so that’s where this effect will have the biggest impact. Focusing then on the difference between the $\lambda=26-34$ bin and the 20-26 bin we see the difference between these bins is about 5% of the true ellipticity (left panel, for the ellipticity values $\sim 0.25$). But this is an overestimate, since 10% of lambda $\sim 30$ is just a shift of 3. So we have, at worst, $\sim 3\%$ bias in ellipticity. The uncertainty in background galaxies gives us an $\sim 10\%$ systematic uncertainty (see Sec. 3.3). So the error of $\lesssim$ 3% from satellite clustering is much smaller and we neglect it.

B.3 Broad ellipticity distribution

In Sec. 3 we use MC simulations to correct for biases in ellipticity of satellite distribution. In doing so, we swept all variations of ellipticity away into an average correction factor. Such variations are expected from (i) intrinsic differences between halos (some are more elliptical than others) and (ii) differences in their orientiation with respect to the line-of-sight. Thus, our noise and edge bias corrections could be misestimated if the distribution of ellipticity is very broad or highly skewed. While the corrections in Figures 3 and 4 are fairly smooth functions of input ellipticity, we perform an additional test to verify the accuracy of our simulations.

If the ellipticity distributions of simulation and the data are similar to each other, using the average value will not result in any bias. Thus, we compare distribution of raw ellipticity of data to that of our MC simulation (Sec. 3) in each richness bin. We take ellipticity distributions from the same MC simulation described in Sec. 3, before average ellipticity values are calculated. The results are shown in Fig. 11. We compare three distributions, (a) raw ellipticity measurement of data (black distribution, the black points in Fig. 2 are the average values of this distribution), (b) measured ellipticity in our MC simulation with interloper galaxies considered (green distribution) and (c) measured ellipticity in the MC simulation without including any interloper (blue distribution). Note that for the comparison to be accurate, we must assume the true ellipticity retrieved in Sec. 3, 0.303 and 0.240 respectively for the case (b) and (c) (Fig. 2).

For the first four richness bins with $\lambda\lesssim 60$ the distributions match very well. The ellipticity distributions of some richer cluster bins are slightly broader than or skewed relative to the data. But these distributions are also noisier, and statistical errors on the ellipticity itself are higher. We conclude that any additional systematic error from averaging over the broad ellipticity distribution is well below the uncertainty due to interlopers, our most important systematic.

B.4 Impact of Richness errors on biases corrections

Edge bias (Sec. 3.2.2) could be inaccurate due to errors in cluster richness ($\lambda$) of redMaPPer . If $\lambda$ is underestimated, we consequently underestimate the size and the mass of the underlying dark matter halo. This would cause an undercorrection of the edge bias if the edge bias correction is a strong function of $\lambda$. However, fortunately, the edge bias is insensitive to $\lambda$. According to the middle panel of Fig. 3 and 4, one can readily conclude that the edge bias is very weak function of $\lambda$. Even in an extreme case assuming that $\lambda\sim 20-26$ of a cluster is mis-estimated to $\lambda\sim 100-200$, we would only have a few percent shift on the edge bias. It is well under our systematic uncertainty of $\sim 12\%$ from the uncertainty in the number of interlopers in clusters. Thus we conclude that richness error does not alter our conclusion.

On the other hand, noise bias is a function of number of observed member galaxies. If unobserved galaxies were included in the calculation, we would obtain a smaller value of the raw ellipticity measurement due to the increased number of galaxies (left panels of Fig. 3 and 4). However, larger number of member galaxies implies smaller noise bias correction. As a result, the same true ellipticity will be retrieved no matter how many unobserved member galaxies there are. Since our aim is to obtain one representative value of ellipticity over all the richness range, the error on richness does not affect our result.

Appendix C Misalignment

The dilution factor described in Sec. 4 is shown in Fig. 12 for ${\rm\epsilon_{true}}=0.271$, the true mean ellipticity of satellite distribution measured in Sec. 3. Since the major axes along which we align the clusters are calculated with satellite distributions, we must derive D assuming the true mean ellipticity to be that of the satellite distribution. Once the value of ellipticity is assumed, we use the same MC simulation as described in Sec. 3 to obtain the dispersion of misalignment angle of major axes of satellite distribution with respect to the true major axes of the halos. From this we calculate the misalignment factor D according to the Eq. (30). For illustration, we also show several alignment distributions for various $N_{\rm sat}$ and $\epsilon$ in Fig. 13.