Analysis of the Amplitude of the Sunyaev–Zel'dovich Effect out to Redshift z = 0.8

We use the sample of clusters obtained by Wen et al. (2012), from the SDSS (Sloan Digital Sky Survey)-III survey. The spatial coverage of that survey is ∼14,000 square degrees and the catalog contains 132,684 clusters in the redshift range of $0.05\leqslant z\leqslant 0.80$. According to those authors, the catalog is more than 95% complete for clusters with a mass of ${M}_{200}\gt {10}^{14}$ ${M}_{\odot }$ in the range of $0.05\leqslant z\leqslant 0.42$ and contains a false detection rate of less than 6%. The catalog presents photometric redshifts for all the clusters, while only ∼30% have determination of redshift based on the spectroscopy of their brighter cluster galaxy (an updated version of the catalog with 52,683 spectroscopic redshifts was presented by Wen & Han 2013). From these cases, it was estimated for the photometric redshifts, a systematic offset $\lt 0.004$ and a standard deviation $0.025\mbox{--}0.030$ for $z\lt 0.45$ and $0.030\mbox{--}0.060$ for $z\gt 0.45$. Throughout the paper, we use the masses M₂₀₀ as estimated from optical richness, as the mass within the radius r₂₀₀ (the radius within which the mean density of a cluster is 200 times the critical density of the universe), derived through (Wen et al. 2012)

$$\begin{eqnarray}&&{M}_{200}=3.2\times {10}^{12}\times {\mathrm{Richness}}^{1.17}\ {M}_{\odot }.\end{eqnarray} \tag{ 1 }$$

The "Richness" is defined by $\tfrac{{L}_{200}}{{L}_{* }}$, where L_* is the characteristic luminosity of galaxies in the r band, and L₂₀₀ is the total r-band luminosity within r₂₀₀ in the SDSS r-band by summing luminosities of member galaxy candidates brighter than an evolution-corrected absolute magnitude ${M}_{r}^{e}(z)={M}_{r}(z)+Q\,z\leqslant -20.5$, with a passive evolution of Q = 1.62.

In Equation (1), the exponent has an error bar of ${\rm{\Delta }}e=0.03$ (Wen et al. 2012), which means that the correct masses ${M}_{200}^{\mathrm{corr}.}$ are related to the calculated M₂₀₀ by ${M}_{200}^{\mathrm{corr}.}={M}_{200}\,\times {\left(\tfrac{{M}_{200}}{3.2\times {10}^{12}{M}_{\odot }}\right)}^{{\rm{\Delta }}e/e}$. The error in the normalization of the mass might be important, but the dependence with the mass has only an exponent of $\lesssim 0.025$, which is unimportant for the present calculations. However,  Ford et al. (2014) give an exponent 1.4 ± 0.1 instead of 1.17 ± 0.03. If there was an error of ${\rm{\Delta }}e\approx 0.2$, we would have an extra dependence proportional to ${M}_{200}^{0.17}$, which is more important, but it will not wreck the main results obtained in this paper: adopting a Ford et al. (2014) exponent would imply that we would need to substitute M₂₀₀ for ${M}_{200}^{0.85}$ in this paper.

To extend the studies to lower redshifts, we also use the sample of 1059 Abell clusters with available spectroscopic redshifts.³

Our analysis uses a cluster catalog that covers masses much lower than those of the PLANCK catalog. Wen et al. clusters have also been used by the Planck Collaboration (2016d), but in a more restrictive way, without paying attention to the CMB temperature or the dependence on the redshift and without the more careful dust substraction that we have carried out here.

The Planck satellite (Planck Collaboration 2016a) was an ESA mission launched in 2009 that surveyed the full sky in nine independent frequency bands provided by the LFI (30, 44, and 70 GHz) and by the HFI (100, 143, 217, 353, 545, and 857 GHz) instruments, and with better angular resolutions as compared to previous surveys in this range of frequencies between 33 arcmin (at 30 GHz) and 5 arcmin (at 217–857 GHz). These characteristics make it an ideal observatory not only to study the CMB anisotropies, but also the SZ effect in galaxy clusters and many other topics. Its frequency coverage allows us to trace the negative ($\nu \lt 217$ GHz) and the positive ($\nu \gt 217$ GHz) parts, as well as the null ($\nu =217$ GHz) of the SZ spectrum. Furthermore, its fine angular resolution allows us to explore high-redshift clusters. This implied a great improvement over WMAP (Bennett et al. 2013), the previous NASA mission which, with a highest frequency of 94 GHz and finest angular resolution of 13.2 arcmin only allowed us to study the SZ effect in nearby and rich clusters. Another important asset of Planck wide frequency coverage is that it permits a much better removal of foreground contaminants.

We have used the DR2 data release (Planck Collaboration 2016a) in the form of full-sky maps projected in Healpix pixelization (Górski et al. 2005), downloaded from the Planck Legacy Archive,⁴  and with pixel resolutions of ${N}_{\mathrm{nside}}=1024$ (pixel size of 3.4 arcmin) for the LFI and ${N}_{\mathrm{nside}}=2048$ (pixel size of 1.7 arcmin) for the HFI. In the LFI, we have used only the 70 GHz frequency band. The lower frequencies were dropped from the analysis owing to their coarse angular resolution. On the other hand, we use the six HFI frequency bands, the highest ones being very important for subtraction of the thermal dust component. We select the data within the region outside of mask 2 of PLANCK, 40% of the sky, which avoids the low Galactic latitude regions and other lines of sight with high Galactic dust columns (Planck Collaboration 2014b). By removing the lowest latitude regions, apart from excluding the highest Galactic dust regions, we also exclude the fields with higher foreground of CO. The foreground and primary CMB anisotropies are not correlated with the positions of clusters, they only introduce uncorrelated noise in the analysis (we carry out some experiments using the Planck PR1 CMB map to check it); CMB anisotropies would be important if we analyzed individual clusters, but for a stacking of clusters this noise is quite low, so we do not remove them.

To extend the spectral range analyzed and to better constrain the spectral shape of the dust function emission, we also used the IRAS map at 100 μm. In particular, we used the new generation of IRAS processed images, which incorporates improvement on the zodiacal light subtraction, on destripping and have a zero level compatible with COBE-DIRBE. The generation of the map is described by Miville-Deschenes & Lagache (2005; IRIS maps). For this map, we also used the maps in HEALPIX (Górski et al. 2005) format with resolution ${N}_{\mathrm{nside}}=2048$.

Due to the low emission of a single cluster, we used a statistical approach in which the emission of all the clusters within given ranges of redshift and masses (as estimated from optical richness by Wen et al. 2012) are stacked around their centers and averaged (see Gutiérrez & López-Corredoira 2014 for a full description of the method).

Possible asymmetries or irregularities in each cluster are randomly oriented and independent between clusters, so stacking the emission from many clusters will largely smooth those asymmetries, and then in the stacked map we can ignore any angular dependence and consider only the radial profile. Two different binning methods (denoted through this paper as modes A and B) were used. In A, the ranges in redshift and masses were selected in order to have roughly the same number of clusters in each bin. In B, we defined bins with the same width in redshift (0.12) and mass (0.3 dex). The number of low-mass or low-redshift clusters is comparatively higher than those of high-mass or high-redshift clusters (see Tables 4, 5). In the catalog by Wen et al. (2012), this procedure produces bins with a very different number of objects.

A sky subtraction is done by subtracting for each stacked map of clusters the average background. Galactic dust was not previously subtracted. The subtraction of a Galactic dust foreground would also partially or totally subtract the dust emission of the cluster (without K-correction) and we prefer to subtract separately later, taking into account the redshift of the clusters.

We proceed to apply an aperture photometry of the stacked maps for each frequency ν (see Figure 1). For all clusters, we measure the flux in the central ${\theta }_{\max }=10\buildrel{\,\prime}\over{.} 5$, which gathers most of the flux of the clusters. The procedure is as follows.

Zoom In Zoom Out Reset image size

First, for each of the bins, we calculate the surface brightness ${f}_{\nu }(\theta )$ in the stacked maps with rings between angular distances θ and $\theta +d\theta$ (we use $d\theta =1^{\prime}$), where R is the angular distance to the cluster center. Then, we calculate the average sky emission (${S}_{\mathrm{sky},\nu }$) as the weighted average of surface brightness in rings between $\theta =20^{\prime}$ and $\theta =36^{\prime}$, where the weight is proportional to the number of pixels in each ring. Cluster emissions are within $\theta \lt 20^{\prime}$, and a sky ring close to the cluster is preferable to minimize the effect of the average gradients of sky emission; that is the reason for the radii selection. With that, we get a first profile of the cluster: ${f}_{1,\nu }(\theta )={f}_{\nu }(\theta )-{S}_{\mathrm{sky},\nu }$. An example is given in Figure 1 for the stacked bin A1 in the frequency 545 GHz.

Several clusters may be observed in a given line of sight. This is corrected by the method explained in Appendix A.1. Also, the loss of flux due to beam dilution is corrected with the method explained in Appendix A.2.

The fluxes F_ν obtained by this method were compared for some single clusters detected also by Planck to the values given by Luzzi et al. (2015), and we get very similar results. For instance, for the cluster Abell 85, we get the fluxes of −7, −6, −4, 0, +10 μK arcmin² respectively at 70, 100, 143, 217, 353 GHz, similar to those obtained by Luzzi et al. (2015, Figure 4 bottom) within the error bars (typically less than 20%–30% for each point, except for the point of 70 GHz in which the error is larger).

A minor correction is carried out by recalculating the effective frequencies of each of the filters, rather than taking the nominal value.

For each of the filters, we calculate the frequency averaged on the detected flux convolved with the filter transmission as

$$\begin{eqnarray}&&{\nu }_{\mathrm{eff}.}(\nu )=\displaystyle \frac{{\int }_{0}^{\infty }d\nu ^{\prime} \,\nu ^{\prime} \,{T}_{\nu }(\nu ^{\prime} ){F}_{\mathrm{dust}}(\nu ^{\prime} )}{{\int }_{0}^{\infty }d\nu ^{\prime} \,{T}_{\nu }(\nu ^{\prime} ){F}_{\mathrm{dust}}(\nu ^{\prime} )},\end{eqnarray} \tag{ 3 }$$

where ${T}_{\nu }(\nu ^{\prime} )$ is the transmission function of the filter with nominal frequency ν. Once we have these new frequencies ${\nu }_{\mathrm{eff}.}(\nu )$, we set the values of the fluxes ${F}_{{\nu }_{\mathrm{eff}.}(\nu )}$ as equal to the previously measured F_ν and we redo the fit given by Equation (2), obtaining new parameters ${A}_{\mathrm{dust}}$, ${\beta }_{\mathrm{dust}}$, ${T}_{\mathrm{dust}}$ that take into account this correction of frequencies. Those are the values that we use.

The analysis of the CMB temperature is not the main aim of this paper. Nonetheless, with the present numbers of Tables 4 and 5, we could see that our data are compatible with the standard model prediction ${T}_{\mathrm{CMB}}=2.726(1+z)$ K. In particular, we fitted our data (including Abell clusters) to ${T}_{\mathrm{CMB}}={T}_{0}{(1+z)}^{\beta }$, with ${T}_{0}=2.78\pm 0.14$, $\beta =-0.14\pm 0.19$ for Mode A, and ${T}_{0}=2.88\pm 0.11$, $\beta =-0.25\pm 0.22$ for Mode B (see Figure 3). This confirms previous results obtained from the Sunyaev–Zel'dovich effect on CMBR (Luzzi et al. 2009, 2015; Avgoustidis et al. 2012; Hurier et al. 2014; Saro et al. 2014; de Martino et al. 2015) or with absorptions in the line of sight of quasars (Molaro et al. 2002; Noterdaeme et al. 2011; Muller et al. 2013). Therefore, our analysis is a corroboration of previous results using clusters with much lower masses than in previous papers, though our errors are worse since we have not optimized our data selection for the analysis of T_CMB.

Zoom In Zoom Out Reset image size

We bin the data with the same redshift and different masses (as estimated from optical richness), or with the (approximately) same mass and different redshifts, doing weighted averages in all cases (the weight of each point is inversely proportional to the square of the error bar). We use different functional shapes to explore the best fits:

$$\begin{eqnarray}&&\displaystyle \frac{{L}_{\mathrm{SZ}}}{{(1+z)}^{4}}={L}_{1}{\left(\displaystyle \frac{{M}_{200}}{{10}^{14}{M}_{\odot }}\right)}^{{\alpha }_{1}},\end{eqnarray} \tag{ 7 }$$

$$\begin{eqnarray}&&\displaystyle \frac{{L}_{\mathrm{SZ}}}{{(1+z)}^{4}}={L}_{2}{(1+z)}^{{\alpha }_{2}}.\end{eqnarray} \tag{ 8 }$$

The dependence on mass and redshift of both variables is shown in Figures 4 and 5. We include in the calculation $\mathrm{Error}(z)=\tfrac{\mathrm{rms}(z)}{\sqrt{{N}_{\mathrm{cl}}}}$. Also, we have tried a fit with the 25 points of Wen et al. clusters without binning, with weights, and a double dependence of the shape:

$$\begin{eqnarray}&&\displaystyle \frac{{L}_{\mathrm{SZ}}}{{(1+z)}^{4}}={L}_{3}{\left(\displaystyle \frac{{M}_{200}}{{10}^{14}{M}_{\odot }}\right)}^{{\alpha }_{3}}{(1+z)}^{{\alpha }_{4}},\end{eqnarray} \tag{ 9 }$$

The obtained numbers of these fits are in Table 1.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Table 1.  Parameters of the Dependence on Mass M₂₀₀ (as Estimated from Optical Richness by Wen et al. 2012) and Redshift of ${L}_{\mathrm{SZ}}$ (Units: 10¹⁸ W K⁻¹ GHz⁻²) for Bins of Mode A and Mode B

Parameter	Mode A	Mode B
L1	3.65 ± 0.81	3.93 ± 0.77
${\alpha }_{1}$	2.57 ± 0.30	1.98 ± 0.12
L2	1.78 ± 2.64	12.0 ± 10.2
${\alpha }_{2}$	7.24 ± 3.91	1.83 ± 2.12
L3	2.81 ± 1.13	5.14 ± 2.41
${\alpha }_{3}$	1.66 ± 0.35	1.78 ± 0.22
${\alpha }_{4}$	3.90 ± 0.79	1.49 ± 1.35

Download table as:  ASCIITypeset image

There are two dependences with detection $\gt 2\sigma$:

• 1.  
  There is a dependence of the amplitude luminosity with ${M}_{200}$. We obtain power laws with exponents ${\alpha }_{1}$ or ${\alpha }_{3}$ between 1.66 and 2.57, being the best measurement with lower error 1.98 ± 0.12 (${\alpha }_{1}$ for Mode B), which is not far (though still at 2–3σ) from the self-similar value (5/3; Kaiser 1986), from what is expected for this mass range according to simulations (Bonaldi et al. 2007; Sembolini et al. 2013, 2014), and from what is obtained with observations (e.g., Andersson et al. 2011; Saliwanchik et al. 2015; Sereno et al. 2015). This value of 5/3 stems from ${L}_{\mathrm{SZ}}\propto {N}_{{\rm{e}}}\langle {T}_{{\rm{e}}}\rangle$ (Equation (6)), and ${N}_{{\rm{e}}}\propto {M}_{\mathrm{gas}}\propto {M}_{200}$, $\langle {T}_{{\rm{e}}}\rangle \propto {M}_{200}^{2/3}$ (Sembolini et al. 2013). With the same data and different analysis by the Planck Collaboration (2016d, Section 4.5), an exponent of 1.92 ± 0.42 is obtained, which is a similar value to ours, though our error bars are much smaller. Indeed, the Planck Collaboration (2016d, Section 4.5) used a slightly different mass proxy based on the same data: ${M}_{200}^{\mathrm{Wen}-\mathrm{Planck}}=(1.1\times {10}^{14}\ {M}_{\odot }){\left(\tfrac{{N}_{200}}{20}\right)}^{1.2}$, where N₂₀₀ is the number of galaxies brighter than a given magnitude within radius R₂₀₀, instead of our Equation (1). With a fit for all of the clusters of Wen et al. (2012), we can see that both amounts are related, on average, as ${M}_{200}^{\mathrm{Wen}-\mathrm{Planck}}=0.77\,{M}_{200}^{0.918}{(1+z)}^{-0.833};$ if we apply this relationship, we get for Mode B ${\alpha }_{1}=1.94\pm 0.20$ (instead of 1.98 ± 0.12).
• 2.  
  There is a dependence of $\tfrac{{L}_{\mathrm{SZ}}}{{(1+z)}^{4}}$ with redshift, once we discount the dependence on mass (${\alpha }_{4}$) for mode A; in mode B, there is also some signal but it is less significant. This is also reflected in the exponent ${\alpha }_{2}$, but here it is not significant and part of the dependence might be due to a different relative weight of the high-mass cluster at higher redshifts. We plot $\tfrac{{L}_{\mathrm{SZ}}}{{(1+z)}^{4}{\left(\tfrac{{M}_{200}}{{10}^{14}{M}_{\odot }}\right)}^{{\alpha }_{3}}}$ versus redshift, assuming the measured value of ${\alpha }_{3}$, in Figure 6. There is clearly a trend to increase $\tfrac{{L}_{\mathrm{SZ}}}{{(1+z)}^{4}}$ for a given cluster mass with redshift, at the $5\sigma$ level for mode A. The differences of the points with the fit are most likely to irregularities in the mass dependence at low mass.

Zoom In Zoom Out Reset image size

The dependence of ${L}_{\mathrm{SZ}}/{(1+z)}^{4}$ on mass is expected, but its dependence with the redshift is not. According to Equation (6), we see that L_SZ is proportional to ${N}_{{\rm{e}}}\langle {T}_{{\rm{e}}}\rangle$. Theoretical analyses of these dependences, assuming hydrostatic equilibrium and isothermal distribution for dark matter and gas particles within standard cosmology, give ${N}_{{\rm{e}}}\propto {M}_{\mathrm{gas}}$ and ${T}_{{\rm{e}}}\propto {[{M}_{200}E(z)]}^{2/3}$, where $E(z)=\sqrt{{{\rm{\Omega }}}_{m}{(1+z)}^{3}+{{\rm{\Omega }}}_{{\rm{\Lambda }}}};$ hence $\tfrac{{L}_{\mathrm{SZ}}}{{(1+z)}^{4}}\propto {f}_{\mathrm{gas}}{M}_{200}^{5/3}E{(z)}^{2/3}$ (Sembolini et al. 2013). The factor $E{(z)}^{2/3}$ would be between 1.0 and 1.3 for redshifts between 0 and 0.7, which is much lower than the average factor ${(1+z)}^{{\alpha }_{4}}$ with ${\alpha }_{4}\approx 4$, which is between 1.0 and 8.3 for redshifts between 0 and 0.7. The simulations using Marenostrum-MultiDark SImulations of galaxy Clusters (MUSIC; Sembolini et al. 2013) do not predict our observed evolution. Nonetheless, Sembolini et al. (2013) claim to observe in their simulations a more conspicuous evolution for lower masses, which is in line with what we observe; the difference is that they see very little evolution in $\tfrac{{L}_{\mathrm{SZ}}E{(z)}^{-2/3}}{{(1+z)}^{4}}$ and we see a factor of ∼6 for that quantity between z = 0 and z = 0.7.

The photometric redshifts have some small errors for each individual cluster (statistical errors of ${\sigma }_{z}=0.025\mbox{--}0.030$ for $z\lt 0.45$ and ${\sigma }_{z}=0.030\mbox{--}0.060$ for $z\gt 0.45;$ and almost null systematic errors; Wen et al. 2012) but they are negligible and they cannot be responsible for the present result. The statistical errors almost cancel when we do the stacking (error of the redshift equal to ${\sigma }_{z}/\sqrt{N}$).

The incompleteness of the cluster sample is estimated by Wen et al. (2012, Section 2.5). The detection rate is nearly 100% for clusters with masses ${M}_{200}\gt 2\times {10}^{14}$ M_⊙ up to redshifts of $z\approx 0.5$, and ∼75% for clusters of ${M}_{200}\gt 6\,\times {10}^{13}$ M_⊙ up to redshifts of $z\approx 0.42$. In fact, more than 95% of clusters of ${M}_{200}\gt {10}^{14}$ M_⊙ are detected for $z\lt 0.42$. This incompleteness should not affect our results because our results do not depend on the counts of clusters and the undetected clusters are not expected to be different from the detected clusters. The only effect may be a variation of the average mass with redshift, but this is already taken into account in the fits because we calculate the average for each bin. Indeed, the variation of the average mass within the five bins An, Bn with $n=m,m+5,m+10,m+15,m+20$, with fixed $m=1\mbox{--}5$ (see Tables 4 and 5), is very small, reflecting that the incompleteness is only slightly varying with redshift at least within $z\lt 0.4$ (see Figure 6 of Wen et al. 2012).

We check that the incompleteness effects are not important by repeating our fit of Equation (9), removing the bins with $\langle {M}_{200}\rangle \lt {10}^{14}$ ${M}_{\odot }$ or $\langle z\rangle \gt 0.42$, that is, we use only the bins A4, A5, A9, A10, A14, A15; and B2-5, B7-10, B12-15. The result is ${\alpha }_{3}=3.23\pm 0.96$, ${\alpha }_{4}=3.42\pm 1.61$ for Mode A, and ${\alpha }_{3}=1.96\pm 0.27$, ${\alpha }_{4}=-1.79\pm 2.16$ for Mode B. For other fits, there are not enough points. These numbers are compatible with those obtained for the whole sample (Table 1) within 1.5σ, though the error bars here are much larger due to the use of a lower number of points for the fits. Anyway, we still see a self-similar value in the dependence with the mass (${\alpha }_{3}\,=1.96\pm 0.27$ in Mode B) and we still see a significant evolution with redshift (${\alpha }_{4}=3.42\pm 1.61$ in Mode A), so we can say that the selection effects are not responsible for these trends.

Another possibility to explain the present result is that the masses ${M}_{200}$ were not correctly calculated and that an average variation of $\langle {M}_{200}\rangle$ with the redshift exists. We cannot attribute this to the evolution of galaxies, since Wen et al. (2012) have taken it into account. However, their approach, which is only an extrapolation of the evolution at low z, may possibly introduce some bias into their estimation of mass. In particular, at $z\gtrsim 0.42$, some galaxies are not observed and the catalog is incomplete, so the calculated richness will contain some bias.

In order to explore the hypothesis of this miscalculation of masses in Wen et al., we have cross-correlated the ${M}_{200}$ masses given by those authors with the masses of 414 common objects (coincident position within 1 arcmin) of the catalog MCXC of clusters detected by Planck and with masses ${M}_{500,\mathrm{MCXC}}$ derived from X-ray (Piffaretti et al. 2011). Figure 7 (top) gives the dependence of the ratio of $\tfrac{{M}_{200}}{{M}_{500,\mathrm{MCXC}}}$ with redshift. There is also some dependence on mass. On average, using again a fit of the type $\tfrac{{M}_{200}}{{M}_{500,\mathrm{MCXC}}}={R}_{5}{\left(\tfrac{{M}_{200}}{{10}^{14}{M}_{\odot }}\right)}^{{\alpha }_{5}}{(1+z)}^{{\alpha }_{6}}$, we get ${R}_{5}=1.39\pm 0.11$, ${\alpha }_{5}=0.34\pm 0.04$, and ${\alpha }_{6}=-1.72\pm 0.27$. These values of ${\alpha }_{5}$ and ${\alpha }_{6}$ may change slightly if we select clusters within some limited range of redshift or masses instead of the whole sample of 414 clusters, but the results are similar: there is no significant variation of these exponents with redshift and the variation with mass is more conspicuous, but compatible with the global analysis within the errors. It is a very significant departure from an almost constant $\tfrac{{M}_{200}}{{M}_{500,\mathrm{MCXC}}}$ that would be expected theoretically if the estimations of both masses were correct: $\approx 1.40{(1+z)}^{0.088}$ for Navarro–Frenk–White (NFW) density distribution, or 1.58 for a singular isothermal distribution (see Appendix C). Taking NFW, we see that the average correction to apply to M₂₀₀ given by Wen et al. in order to be compatible with the values of M₅₀₀ of MCXC is

$$\begin{eqnarray}{M}_{200}^{\mathrm{corr}.} & = & (1.01\pm 0.08){M}_{200}{\left(\displaystyle \frac{{M}_{200}}{{10}^{14}{M}_{\odot }}\right)}^{-0.34\pm 0.04}\\ & & \times {(1+z)}^{1.81\pm 0.27}.\end{eqnarray} \tag{ 10 }$$

Zoom In Zoom Out Reset image size

With this correction, we can calculate again the parameters ${\alpha }_{i}$ of the best fit, substituting M₂₀₀ for ${M}_{200}^{\mathrm{corr}.}$. The new results are given in Table 2. Here we have not grouped the bins of the same mass because, with the correction, the dispersion of values is higher. The systematic errors of Equation (10) are, however, considered in the error bars of the derived parameters; the statistical errors for each bin of clusters are considered negligible given the large number of clusters per bin. The result is again problematic but in a different sense: now we see that dependence of the luminosity with mass has an exponent larger than 5/3, between 2.2 and 2.7 (maximum significance of the difference: 2.8σ in ${\alpha }_{3}$ for Mode B); this is illustrated in Figure 8. Furthermore, either there is no evolution of luminosity with redshift (for mode A) or there is an evolution with redshift opposite to what was observed previously (${\alpha }_{4}\lt 0$ for mode B). We solve with it a result that was not understood with the masses of Wen et al. (2012) but we get a new problematic result: the self-similar value of the dependence of the luminosity mass with the mass is not obtained, which might be due to the fact that the clusters are not virialized and in the conditions of the simulations that obtain the exponent of 5/3. This is something that was already questioned by Czakon et al. (2015), but they relate it to possible biases. However, again, we may doubt the validity of X-ray masses of the MXCX catalog. The method to derive these masses may introduce some bias depending on the mass (von der Linden et al. 2014; Hoekstra et al. 2015).

Zoom In Zoom Out Reset image size

Table 2.  Parameters of ${L}_{\mathrm{SZ}}$ (Units: 10¹⁸ W K⁻¹ GHz⁻²) in the Dependence on the Corrected Mass ${M}_{200}^{\mathrm{corr}.}$ (See Equation (10)) and Redshift for Bins of Mode A and Mode B

Parameter	Mode A	Mode B
L1	3.22 ± 1.29	3.34 ± 1.69
${\alpha }_{1}$	2.20 ± 0.45	2.55 ± 0.36
L2	1.78 ± 2.64	12.0 ± 10.2
${\alpha }_{2}$	7.24 ± 3.91	1.83 ± 2.12
L3	2.87 ± 1.31	5.01 ± 2.61
${\alpha }_{3}$	2.51 ± 0.55	2.70 ± 0.37
${\alpha }_{4}$	−0.78 ± 1.12	−3.39 ± 1.52

Note. Errors include the errors of the fits and the uncertainties in the parameters of Equation (10).

Download table as:  ASCIITypeset image

At least one of the two mass estimators must be wrong, and it is possible that both are wrong. It cannot be a geometric or morphological question in the cluster given by the different radii up to which the mass or luminosity is integrated: we show in Appendix C that the difference between M₂₀₀ and M₅₀₀ should be a constant factor, almost independent of redshift. Our guess is that the evolution correction measured by Wen et al. (2012) in order to convert richness into mass is not accurate enough, and it may be possibly affected by some bias of incompleteness for clusters with $z\gt 0.42$, which do not allow us to count all the galaxies embedded in a cluster and thus produce an error in the estimation of the richness. In a previous work, Wen et al. (2009) claimed that their mass determinations were compatible with X-ray measurements, but these comparisons were limited to redshifts $z\lt 0.26$, and, moreover, they used a different definition of richness. It is at higher redshift where we find significant differences between the masses of Wen et al. (2012) and MCXC masses. The analysis done by Wen & Han (2015) has also found that there is a difference between Wen et al. (2012) masses and the X-ray masses, which may be corrected with a better estimation of the evolution factor in the cluster richness–luminosity relationship. On the other hand, for X-ray measurements of MCXC masses, we may suspect that the method of Piffaretti et al. (2011) assuming a particular mass–luminosity ratio may be biased or that the extrapolation of the cluster mass profile to its outer parts may contain some inaccuracy. Piffaretti et al. (2011) also warn that the ${M}_{500,\mathrm{MCXC}}$ values provided by the MCXC rely on the assumption that, on average, the Malmquist bias for the samples used to construct the MCXC is the same as that of the REXCESS sample. The X-ray luminosities are also affected by two main sources of uncertainty (Piffaretti et al. 2011): uncertainties in the underlying model used in the homogenization procedure and the measurement uncertainty in the original input catalog.

Weak-lensing masses of clusters (Applegate et al. 2014; von der Linden et al. 2014; Hoekstra et al. 2015; Sereno 2015; Melchior et al. 2016) are less affected by the bias of the present optical or X-ray derived masses, and they will be useful to solve this question.

We cross-correlate the ${M}_{200}$ masses given by Wen et al. (2012) with the weak-lensing masses of clusters in the compilation by Sereno (2015): there are 194 objects with data in both samples. On average, the fit of type $\tfrac{{M}_{200}}{{M}_{200,\mathrm{lensing}}}\,={R}_{5}{\left(\tfrac{{M}_{200}}{{10}^{14}{M}_{\odot }}\right)}^{{\alpha }_{5}}{(1+z)}^{{\alpha }_{6}}$ gives: ${R}_{5}=0.68\pm 0.13$, ${\alpha }_{5}=0.29\,\pm 0.06$, and ${\alpha }_{6}=-1.84\pm 0.51$. Figure 7 (bottom) gives the dependence of the ratio of $\tfrac{{M}_{200}}{{M}_{200,\mathrm{lensing}}}$ with redshift. Therefore, the average correction to apply to M₂₀₀ given by Wen et al. in order to be compatible with these values of weak lensing is

$$\begin{eqnarray}{M}_{200}^{\mathrm{corr}.2} & = & (1.45\pm 0.27){M}_{200}{\left(\displaystyle \frac{{M}_{200}}{{10}^{14}{M}_{\odot }}\right)}^{-0.29\pm 0.06}\\ & & \times {(1+z)}^{1.84\pm 0.51}.\end{eqnarray} \tag{ 11 }$$

This correction is very similar to that of Equation (10) using X-ray, except only for the amplitude, which lead us to think that X-ray derived masses are not so biased. With this new correction, we can calculate again the parameters ${\alpha }_{i}$ of the best fit substituting M₂₀₀ for ${M}_{200}^{\mathrm{corr}.2}$. The new results are given in Table 3. Again, we see no significant evolution but some departure from the self-similar exponent of 5/3 in the mass dependence, though less significant in this case (maximum significance of the difference: 2.2σ in ${\alpha }_{3}$ for Mode B): a value of ${\alpha }_{1}$ or ${\alpha }_{3}$ between 2.1 and 2.5.

Table 3.  Parameters of ${L}_{\mathrm{SZ}}$ (Units: 10¹⁸ W K⁻¹ GHz⁻²) in the Dependence on the Corrected Mass ${M}_{200}^{\mathrm{corr}.2}$ (See Equation (11)) and Redshift for Bins of Mode A and Mode B

Parameter	Mode A	Mode B
L1	1.42 ± 1.02	1.42 ± 1.04
${\alpha }_{1}$	2.14 ± 0.53	2.40 ± 0.36
L2	1.78 ± 2.64	12.0 ± 10.2
${\alpha }_{2}$	7.24 ± 3.91	1.83 ± 2.12
L3	1.23 ± 0.90	2.02 ± 1.57
${\alpha }_{3}$	2.34 ± 0.54	2.51 ± 0.38
${\alpha }_{4}$	−0.53 ± 2.04	−3.12 ± 1.89

Note. Errors include the errors of the fits and the uncertainties in the parameters of Equation (11).

Download table as:  ASCIITypeset image

In this case, we cannot think about any possible systematic bias, and the fact of the coincidence with the result of X-ray measurements of masses corroborates this scenario as the most likely one: no (significant) evolution and ${L}_{\mathrm{SZ}}\propto \sim {M}^{3}$.

Using the value of L_i with higher signal-to-noise, ${L}_{1}=(3.93\pm 0.77)\times {10}^{18}$ W K⁻¹ GHz⁻² for mode B analysis, Equation (6) allows us to estimate the average redshift of the sample and that of cluster masses of ${M}_{200}={10}^{14}$ ${M}_{\odot }$ for which ${N}_{{\rm{e}}}\langle {T}_{{\rm{e}}}(K)\rangle =(9.07\pm 1.78)\times {10}^{76}$ per cluster. An electronic temperature of ${T}_{{\rm{e}}}\sim 2\times {10}^{7}$ K (Vikhlinin et al. 2006) gives around $5\times {10}^{69}$ electrons per cluster. With ${N}_{{\rm{e}}}={f}_{{\rm{e}}}{N}_{{\rm{H}}}$, where f_e is the number of electrons per hydrogen atom of the cluster gas and N_H is the number of hydrogen atoms, and estimating the mass of the gas as ${M}_{\mathrm{gas}}={N}_{{\rm{H}}}\,{m}_{\mathrm{proton}}\left(1+4\tfrac{Y}{1-Y}\right)$ (Y is the helium abundance; we assume it is equal to 0.25; we neglect the weight of the atoms heavier than helium), we get ${M}_{\mathrm{gas}}\sim {10}^{13}/{f}_{{\rm{e}}}$ ${M}_{\odot }$, that is, $\tfrac{{M}_{\mathrm{gas}}}{{M}_{200}}\sim 0.10\,{f}_{{\rm{e}}}^{-1}$. Roughly, ${f}_{{\rm{e}}}\sim 1$, because most atoms are hydrogen and are ionized, so this means that the gas mass is 10% of the total mass, as previously observed (Vikhlinin et al. 2006).

Several clusters may be observed in a given line of sight. One way to take into account this and calculate the emission per cluster is

$$\begin{eqnarray}&&{f}_{2,\nu }(\theta )={f}_{1,\nu }(\theta )-\displaystyle \frac{1}{2\pi }{\int }_{0}^{2\pi }d\phi ^{\prime} {\int }_{0}^{\infty }d\theta ^{\prime} \,\theta ^{\prime} \end{eqnarray} \tag{ 12 }$$

$$\begin{eqnarray*}&&\times \,C(\theta ^{\prime} )\times {f}_{1,\nu }(\sqrt{{\theta }^{2}+\theta {^{\prime} }^{2}-2\,\theta \,\theta ^{\prime} \,\cos (\phi ^{\prime} )}),\end{eqnarray*}$$

$$\begin{eqnarray*}&&C(\theta )=\sqrt{\langle \sigma \sigma \rangle (\theta )-\langle \sigma {\rangle }^{2}},\end{eqnarray*}$$

where σ is the surface density of clusters. These amounts are obtained from the analysis of our sample. An example of the correction is given in Figure 9 for the stacked bin A1 in the frequency of 545 GHz. For small correlations, as in our case, one iteration is enough; otherwise, we could get a better approximation by calculating in the same way ${f}_{3,\nu }(\theta )$, ${f}_{4,\nu }(\theta )$,...

Zoom In Zoom Out Reset image size

Finally, we derive the integrated flux of the cluster within ${\theta }_{\max }$ by means of

$$\begin{eqnarray}&&{F}_{\nu ,{\theta }_{\max }}=2\,\pi \,{\int }_{0}^{{\theta }_{\max }}d\theta \,\theta \,{f}_{2,\nu }(\theta ).\end{eqnarray} \tag{ 13 }$$

The error bars are calculated with the propagation of errors of all the processes.

There is some small amount of flux that we loose for angular distances larger than ${\theta }_{\max }:$ This is stronger for low frequencies, where the telescope beam profile is wider. Assuming a Gaussian profile, their corresponding sigmas are 567 for 70 GHz (Planck Collaboration 2016b); 412 for 100 GHz, 311 for 143 GHz, 214 for 217 GHz, 210 for 353 GHz, 206 for 545 GHz, 197 for 857 GHz (Planck Collaboration 2016c); and 085 for 2998 GHz (Miville-Deschenes & Lagache 2005). The closest and most massive clusters may have an angular size of R₂₀₀ comparable to ${\theta }_{\max }:$ for instance, at z = 0.12, 10.5 arcmin correspond to 1.4 Mpc (for Abell clusters, the angular size is larger, but they will not be used in the statistics of the amplitude of the flux). Planck beam is very close to Gaussian, especially for frequencies below 545 GHz. Moreover, stacking of clusters produce profiles with high circular symmetry, reducing fluctuations. However, there may be departures from Gaussianity due to the shape of the clusters and the sum of clusters with different angular size in the stacking, due to the dispersion of redshifts and masses in each bin; we neglect them since this correction of flux is of second order.

Numerical experiments on the integration of the flux showed us that increasing the value of the maximum radius of integration with ${\theta }_{\max }\gt 10\buildrel{\,\prime}\over{.} 5$ does not improve the result, because the tail of the Gaussian introduces mostly noise. The fit to a Gaussian would be a solution, but it reduces the accuracy with respect to our method of integration. Instead, we do the following approach of extrapolation: assuming Gaussianity in ${f}_{2,\nu }$ with rms σ, we get from Equation (13) ${F}_{\nu ,{\theta }_{\max }}\propto {\sigma }^{2}\left[1-\exp \left(-\tfrac{{\theta }_{\max }^{2}}{2{\sigma }^{2}}\right)\right]$ and consequently

$$\begin{eqnarray}&&{F}_{\nu }\equiv {F}_{\nu ,\infty }={F}_{\nu ,{\theta }_{\max }}\times \mathrm{factor},\end{eqnarray} \tag{ 14 }$$

$$\begin{eqnarray*}&&\mathrm{factor}=\displaystyle \frac{1}{1-y};\ \ y\equiv \exp \left(-\displaystyle \frac{{\theta }_{\max }^{2}}{2{\sigma }^{2}}\right).\end{eqnarray*}$$

Using another angle ${\theta }_{0}=\tfrac{{\theta }_{\max }}{\sqrt{2}}\approx 7\buildrel{\,\prime}\over{.} 5$, we can derive y from the ratio $r=\tfrac{{F}_{\nu ,{\theta }_{\max }}}{{F}_{\nu ,{\theta }_{0}}}$:

$$\begin{eqnarray}&&r=\displaystyle \frac{1-y}{1-\sqrt{y}},\end{eqnarray} \tag{ 15 }$$

thus,

$$\begin{eqnarray}&&{F}_{\nu }={F}_{\nu ,{\theta }_{\max }}\displaystyle \frac{1}{1-{\left(\tfrac{{F}_{\nu ,{\theta }_{\max }}}{{F}_{\nu ,{\theta }_{0}}}-1\right)}^{2}}.\end{eqnarray} \tag{ 16 }$$

This correction is relatively small, $\lesssim 10$%, being the most significant for 70 GHz.