Cluster Cosmology Constraints from the 2500 deg² SPT-SZ Survey: Inclusion of Weak Gravitational Lensing Data from Magellan and the Hubble Space Telescope

The South Pole Telescope (SPT) is a 10 m telescope located within 1 km of the geographical South Pole (Carlstrom et al. 2011). The ∼1' resolution and 1° field of view are well suited for a survey of rare, high-mass clusters from a redshift of z ≥ 0.2 out to the highest redshifts where they exist. From 2007 to 2011, the telescope was configured to observe with the SPT-SZ camera in three millimeter-wave bands (centered at 95, 150, and 220 GHz). The majority of this period was spent on the SPT-SZ survey, a contiguous 2500 deg² area within the boundaries 20^h ≤ R.A. ≤ 7^h and −65° ≤ decl. ≤ −40°. The survey achieved a fiducial depth of ≤18 μK arcmin in the 150 GHz band.

Galaxy clusters are detected via their thermal SZ signature in the 95 and 150 GHz maps. These maps are created using time-ordered data processing and map-making procedures equivalent to those described in Vanderlinde et al. (2010) and Reichardt et al. (2013). Galaxy clusters are extracted using a multiscale matched-filter approach (Melin et al. 2006) applied to the multiband data as described in Williamson et al. (2011) and Reichardt et al. (2013).

We use the same SPT-SZ cluster sample that was analyzed in dH16. Namely, this cosmological sample is a subset of the full SPT-SZ cluster sample presented in Bleem et al. (2015), restricted to redshifts z > 0.25 and detection significances ξ > 5. This cosmological sample has an expected and measured purity of 95% (Bleem et al. 2015). For clusters at redshifts below z = 0.25, confusion with primary CMB fluctuations changes the scaling of the ξ–mass relation.

We have improved the cluster redshift estimates from the original values provided in Bleem et al. (2015) to incorporate new spectroscopic measurements (Bayliss et al. 2016; Khullar et al. 2019; A. Mantz et al. 2019, in preparation), two updated high-redshift photo-z measurements with HST (Strazzullo et al. 2019), and improved photometric measurements. These improved photometric redshifts are enabled both via the recalibration of our Spitzer redshift models using the new spectroscopic data and by the use of optical data from the Parallel Imager for Southern Cosmology Observations (PISCO), a new imager installed on the Magellan/Clay telescope at Las Campanas Observatory (Stalder et al. 2014). PISCO—with a fast (∼20 s) readout, 9' field of view, and simultaneous four-band (griz) imaging capability—is optimized for efficient characterization of clusters and other systems identified from external surveys. As part of further efforts to characterize the SPT-SZ cluster sample, we have obtained approximately uniform PISCO imaging for the majority of the previously confirmed SPT-SZ clusters. Notably, these deeper optical data have allowed less constraining infrared-driven redshift estimates from Spitzer to be replaced by more robust estimates based on optical red sequence techniques for a significant number of clusters in the range 0.8 ≲ z ≲ 1. As a consequence, while the improved data and model calibration result in small changes in redshift estimates for systems at z ≲ 0.8 and z ≳ 1, at intermediate redshifts, replacing infrared-driven redshifts with more robust optical estimates leads to up to 1.5σ systematic shifts; see Figure 2. We will briefly come back to this issue in Section 4.3.

Zoom In Zoom Out Reset image size

We use X-ray measurements for a subsample of 89 clusters. Eighty-one of these were also used in our previous cosmological analysis (dH16). We decided not to use the X-ray data available for SPT-CL J0142-5032 because of its large measurement error in temperature exceeding 40%. This has a negligible impact on our results. Most of those X-ray measurements were originally presented in McDonald et al. (2013), and they were largely acquired through a Chandra X-ray Visionary Project (PI: Benson). This sample is now supplemented with observations of eight high-redshift z > 1.2 clusters (McDonald et al. 2017). We refer the reader to these references for the details of the X-ray analysis.

The X-ray data products entering this analysis are (i) lookup tables of the total gas mass, ${M}_{\mathrm{gas}}$, within an outer radius ranging from 80 to 2000 kpc (calculated using a fiducial cosmology), allowing interpolation of ${M}_{\mathrm{gas}}$ within any realistic value of r₅₀₀; and (ii) spectroscopic temperatures, T_X, in the 0.15r₅₀₀–1.0r₅₀₀ aperture. All X-ray measurements were remade for this work using the Chandra calibration CALDB v4.7.7. Note that this calibration does not change the results from dH16.

We use WL measurement for 32 clusters in our sample. Of these, 19 were observed with Magellan/Megacam at redshifts 0.29 ≤ z ≤ 0.69 (D19), and 13 at redshifts 0.576 ≤ z ≤ 1.132 with the Advanced Camera for Surveys on board HST (S18). Details on the data reduction and analysis methods can be found in these works.

The data products from these works used in our analysis are the reduced tangential shear profiles in angular coordinates, corrected for contamination by cluster galaxies, and the estimated redshift distributions of the selected source galaxies. These are the observable quantities, which are independent from cosmology, whereas mass estimates or shear profiles in physical coordinates depend on cosmology through the redshift–distance relation and the cosmology dependence of the NFW profile. Our approach ensures a clean separation between the actual measurements and their modeling.

In addition to our cluster data set, we will also consider external cosmological probes. We use measurements of primary CMB anisotropies from Planck and focus on the TT+lowTEB data combination from the 2015 analysis (Planck Collaboration et al. 2016a). We use angular diameter distances as probed by baryon acoustic oscillations (BAOs) by the 6dF Galaxy Survey (Beutler et al. 2011), the SDSS Data Release 7 Main Galaxy Sample (Ross et al. 2015), and the BOSS Data Release 12 (Alam et al. 2017). We also use measurements of luminosity distances from Type Ia supernovae (SNe Ia) from the Pantheon sample (Scolnic et al. 2018).

We consider three cluster mass proxies: the unbiased SZ significance ζ, the X-ray ${Y}_{{\rm{X}}}$, and the WL mass M_WL. We parameterize the mean observable–mass relations as

$$\begin{eqnarray}\begin{array}{rcl}\langle \mathrm{ln}\zeta \rangle & = & \mathrm{ln}{A}_{\mathrm{SZ}}\ +{B}_{\mathrm{SZ}}\ \mathrm{ln}\left(\displaystyle \frac{{M}_{500c}\,{h}_{70}}{4.3\times {10}^{14}{M}_{\odot }}\right)\\ & & +\ {C}_{\mathrm{SZ}}\ \mathrm{ln}\left(\displaystyle \frac{E(z)}{E(0.6)}\right)\end{array}\end{eqnarray} \tag{ 1 }$$

$$\begin{eqnarray}\begin{array}{rcl}\mathrm{ln}\left(\displaystyle \frac{{M}_{500c}\,{h}_{70}}{8.37\times {10}^{13}{M}_{\odot }}\right) & = & \mathrm{ln}{A}_{{Y}_{{\rm{X}}}}\ +{B}_{{Y}_{{\rm{X}}}}\langle \mathrm{ln}{Y}_{{\rm{X}}}\rangle \\ & & +\ {B}_{{Y}_{{\rm{X}}}}\mathrm{ln}\left(\displaystyle \frac{{h}_{70}^{5/2}}{3\times {10}^{14}{M}_{\odot }\mathrm{keV}}\right)\\ & & +\ {C}_{{Y}_{{\rm{X}}}}\mathrm{ln}E(z)\end{array}\end{eqnarray} \tag{ 2 }$$

$$\begin{eqnarray}&&\langle \mathrm{ln}{M}_{\mathrm{WL}}\rangle =\mathrm{ln}{b}_{\mathrm{WL}}+\mathrm{ln}{M}_{500c}.\end{eqnarray} \tag{ 3 }$$

The ζ–mass and ${Y}_{{\rm{X}}}$–mass relations are equivalent to the ones adopted in dH16, except for replacing h/0.72 by h₇₀ in ${Y}_{{\rm{X}}}$–mass.

The intrinsic scatter in ln ζ, ln ${Y}_{{\rm{X}}}$, and ln M_WL at fixed mass and redshift is described by normal distributions with widths ${\sigma }_{\mathrm{ln}\zeta }$, ${\sigma }_{\mathrm{ln}{Y}_{{\rm{X}}}}$, and ${\sigma }_{\mathrm{WL}}$. These widths are assumed to be independent of mass and redshift. Note that the parameters ${\sigma }_{\mathrm{ln}\zeta }$ and ${\sigma }_{\mathrm{ln}{Y}_{{\rm{X}}}}$ have been called D_SZ and D_X in some previous SPT publications. We allow for correlated scatter between the SZ, X-ray, and WL mass proxies as described by the covariance matrix

$$\begin{eqnarray}&&\begin{array}{l}{{\boldsymbol{\Sigma }}}_{\mathrm{multi} \mbox{-} \mathrm{obs}}\\ =\ \left(\begin{array}{ccc}{\sigma }_{\mathrm{ln}\zeta }^{2} & {\rho }_{\mathrm{SZ}-\mathrm{WL}}{\sigma }_{\mathrm{ln}\zeta }{\sigma }_{\mathrm{WL}}\ & {\rho }_{\mathrm{SZ}-{\rm{X}}}{\sigma }_{\mathrm{ln}\zeta }{\sigma }_{\mathrm{ln}{Y}_{{\rm{X}}}}\\ {\rho }_{\mathrm{SZ}-\mathrm{WL}}{\sigma }_{\mathrm{ln}\zeta }{\sigma }_{\mathrm{WL}}\ & {\sigma }_{\mathrm{WL}}^{2} & {\rho }_{\mathrm{WL}-{\rm{X}}}{\sigma }_{\mathrm{WL}}{\sigma }_{\mathrm{ln}{Y}_{{\rm{X}}}}\\ {\rho }_{\mathrm{SZ}-{\rm{X}}}{\sigma }_{\mathrm{ln}\zeta }{\sigma }_{\mathrm{ln}{Y}_{{\rm{X}}}}\ & {\rho }_{\mathrm{WL}-{\rm{X}}}{\sigma }_{\mathrm{WL}}{\sigma }_{\mathrm{ln}{Y}_{{\rm{X}}}}\ & {\sigma }_{\mathrm{ln}{Y}_{{\rm{X}}}}^{2}\end{array}\right)\end{array}\end{eqnarray} \tag{ 4 }$$

with correlation coefficients ρ_SZ−X, ρ_SZ−WL, and ρ_WL−X. With this, the full description of the multi-observable–mass relation is

$$\begin{eqnarray}&&\begin{array}{l}P\left(\left[\begin{array}{c}\mathrm{ln}\zeta \\ \mathrm{ln}{M}_{\mathrm{WL}}\\ \mathrm{ln}{Y}_{{\rm{X}}}\end{array}\right]| M,z,{\boldsymbol{p}}\right)\\ =\ { \mathcal N }\left(\left[\begin{array}{c}\langle \mathrm{ln}\zeta \rangle (M,z,{\boldsymbol{p}})\\ \langle \mathrm{ln}{M}_{\mathrm{WL}}\rangle (M,z,{\boldsymbol{p}})\\ \langle \mathrm{ln}{Y}_{{\rm{X}}}\rangle (M,z,{\boldsymbol{p}})\end{array}\right],\ {{\boldsymbol{\Sigma }}}_{\mathrm{multi} \mbox{-} \mathrm{obs}}\right).\end{array}\end{eqnarray} \tag{ 5 }$$

All parameters of the observable–mass relations are listed in Table 2.

While our default X-ray observable is ${Y}_{{\rm{X}}}$, we also consider the X-ray gas mass ${M}_{\mathrm{gas}}$. Note that both observables share the same ${M}_{\mathrm{gas}}$ data, and so we do not use them simultaneously. We define a relation for the gas mass fraction ${f}_{\mathrm{gas}}$ ≡ ${M}_{\mathrm{gas}}$/M_500c,

$$\begin{eqnarray}\begin{array}{rcl}\langle \mathrm{ln}{f}_{\mathrm{gas}}\rangle & = & \mathrm{ln}\left(\displaystyle \frac{{A}_{{M}_{{\rm{g}}}}}{{h}_{70}^{3/2}}\right)+({B}_{{M}_{{\rm{g}}}}-1)\mathrm{ln}\left(\displaystyle \frac{{M}_{500c}\,{h}_{70}}{5\times {10}^{14}{M}_{\odot }}\right)\\ & & +\ {C}_{{M}_{{\rm{g}}}}\mathrm{ln}\left(\displaystyle \frac{E(z)}{E(0.6)}\right),\end{array}\end{eqnarray} \tag{ 6 }$$

with which the ${M}_{\mathrm{gas}}$–mass relation becomes

$$\begin{eqnarray}\begin{array}{rcl}\Space{0ex}{0.25ex}{0ex}\langle \mathrm{ln}\left(\displaystyle \frac{{M}_{\mathrm{gas}}}{5\times {10}^{14}{M}_{\odot }}\right)\Space{0ex}{0.25ex}{0ex}\rangle & = & \mathrm{ln}\left(\displaystyle \frac{{A}_{{M}_{{\rm{g}}}}}{{h}_{70}^{5/2}}\right)+{B}_{{M}_{{\rm{g}}}}\mathrm{ln}\left(\displaystyle \frac{{M}_{500c}\,{h}_{70}}{5\times {10}^{14}{M}_{\odot }}\right)\\ & & +\ {C}_{{M}_{{\rm{g}}}}\mathrm{ln}\left(\displaystyle \frac{E(z)}{E(0.6)}\right).\end{array}\end{eqnarray} \tag{ 7 }$$

3.1.1. The SZ ξ–Mass Relation

The observable we use to describe the cluster SZ signal is ξ, the detection signal-to-noise ratio (S/N) maximized over all filter scales. To account for the impact of noise bias, the unbiased SZ significance ζ is introduced, which is the S/N at the true, underlying cluster position and filter scale (Vanderlinde et al. 2010). Following previous SPT work, ξ across many noise realizations is related to ζ as

$$\begin{eqnarray}&&P(\xi | \zeta )={ \mathcal N }(\sqrt{{\zeta }^{2}+3},1).\end{eqnarray} \tag{ 8 }$$

In practice, we only map objects with ζ > 2 to ξ using this relation, but the exact location of this cut has no impact on our results (see also dH16). The validity of this approach and of Equation (8) has been extensively tested and confirmed by analyzing simulated SPT observations of mock SZ maps (Vanderlinde et al. 2010).

The SPT-SZ survey consists of 19 fields that were observed to different depths. The varying noise levels only affect the normalization of the ζ–mass relation and leave ${B}_{\mathrm{SZ}}$, ${C}_{\mathrm{SZ}}$, and ${\sigma }_{\mathrm{ln}\zeta }$ effectively unchanged (dH16). In the analysis presented here, ${A}_{\mathrm{SZ}}$ is rescaled by a correction factor for each of the 19 fields, which then allows us to work with a single SZ observable–mass relation, given by Equation (1). The scaling factors γ_field can be found in Table 1 in dH16.

In a departure from previous SPT analyses, we do not apply informative (Gaussian) priors on the SZ scaling relation parameters. The self-calibration through fitting the cluster sample against the HMF (see, e.g., Majumdar & Mohr 2004), the constraint on the normalization of the observable–mass relations through our WL data, and the constraint on the SZ scatter through the X-ray data are strong enough to constrain all four SZ scaling relation parameters (in νΛCDM; see Table 3). When not including the X-ray data in our fit, however, we apply a Gaussian prior ${\sigma }_{\mathrm{ln}\zeta }$ = 0.13 ± 0.13 as in dH16 (this constraint was extracted from mock observations of hydrodynamic simulations from Le Brun et al. 2014).

We discuss possible limitations in our description of the ξ–mass relation that would lead to systematic biases in the recovered cosmological constraints. Because of our empirical weak-lensing mass calibration and the parameterization of the SZ scaling relation by power laws and lognormal scatter with free parameters, any bias in the SZ–mass relation that can be described by a power law and/or lognormal scatter would only lead to parameter shifts in the SZ scaling relation but would not affect cosmological parameter constraints. Therefore, important systematics would be from potential contaminants that would lead to an additional, non-lognormal scatter, a mass or redshift dependence in the scatter, or a redshift dependence of the mass slope.

A potential worry might be the dilution of the SZ signal by AGN activity and the presence of dusty star-forming galaxies in the cluster. Various studies have found that emission by dusty star-forming galaxies is negligible compared to the SZ signal (see, e.g., Lin et al. 2009; Sehgal et al. 2010; and the summary in Section 6.4 in Benson et al. 2013). Gupta et al. (2017b) measured the cluster radio luminosity function using an X-ray-selected cluster sample at z ≲ 0.7 and concluded that radio sources obeying this luminosity function would not have a strong impact on the SZ signal. Only a few percent at most of the SPT-SZ clusters would host sufficiently bright radio sources for their SZ signal to drop below the selection threshold, and this is within the Poisson uncertainty of our sample. At higher redshifts, it has been previously measured that the radio fraction in optically selected clusters somewhat decreases at z > 0.65 (Gralla et al. 2011). This result is consistent with simulations of the microwave sky from Sehgal et al. (2010), which predicted that the amount of radio contamination in SZ surveys was either flat or falling at z > 0.8. Using tests against mocks, we find, for example, that to cause a shift in w by more than Δ(w) = −0.3, the level of SZ contamination would have to be strong enough to remove more than ∼30% of all cluster detections at redshifts z ≳ 1, which by far exceeds the measurement by Gupta et al. (2017b). In conclusion, none of the discussed sources of potential SZ cluster contamination have an impact that is strong enough to introduce large biases in our cosmological constraints.

Another approach to testing the robustness of the SZ observable–mass relation is to compare it with other cluster mass proxies and try to find deviations from the simple scaling relation model. Note that if such a deviation were found, it would be hard to discern which observable is behaving in an unexpected way, but importantly, one would learn that the multi-observable model needs an extension. At low and intermediate redshifts z ≲ 0.8, comparisons with cluster samples selected through optical and X-ray methods have shown that the cluster populations can be described by power-law observable–mass scaling relations with lognormal intrinsic scatter (Vikhlinin et al. 2009a; Mantz et al. 2010a, 2016; Saro et al. 2015, 2017). At higher redshifts, the subset of the SPT-selected sample with available X-ray observations from Chandra and XMM-Newton exhibits scaling relations in X-ray T_X, ${Y}_{{\rm{X}}}$, ${M}_{\mathrm{gas}}$, and L_X, as well as in stellar mass galaxies, that are consistent with power-law relations in mass and redshift with lognormal intrinsic scatter (Chiu et al. 2016, 2018; Hennig et al. 2017; Bulbul et al. 2019). When a redshift-dependent mass slope parameter has been included in the analyses of these data sets, the parameter constraints have been statistically consistent with 0 in all cases (see Table 4 in Bulbul et al. 2019).

In conclusion, our description of the ξ–mass relation has been confirmed by various independent techniques, especially for redshifts z ≲ 1. Note that these tests are harder to perform at higher redshifts, where non-SZ-selected samples are small and more challenging to characterize. Our expectation is that as the cluster sample grows larger and the mass calibration information improves we will be able to characterize the currently negligible departures from our scaling relation model. At that point, we will need to extend our observable–mass relation to allow additional freedom.

3.1.2. The Weak-lensing Observable–Mass Relation

The WL modeling framework used in this work is introduced in D19, and we refer the reader to their Section 5.2 for details.

The WL observable is the reduced tangential shear profile g_t(θ), which can be analytically modeled from the halo mass M_200c, assuming an NFW halo profile and using the redshift distribution of source galaxies (Wright & Brainerd 2000). Miscentering, halo triaxiality, large-scale structure along the line of sight, and uncertainties in the concentration–mass relation introduce bias and/or scatter. As introduced in Equation (3), we assume a relation ln M_WL = ln(${b}_{\mathrm{WL}}$M_true) and use numerical simulations to calibrate the normalization ${b}_{\mathrm{WL}}$ and the scatter about the mean relation. Our WL data set consists of two subsamples (Megacam and HST) with different measurement and analysis schemes. We expect some systematics to be shared among the entire sample, while others will affect each subsample independently.

We model the WL bias as

$$\begin{eqnarray}\begin{array}{rcl}{b}_{\mathrm{WL},i} & = & {b}_{\mathrm{WL}\mathrm{mass},i}\\ & & +\ {\delta }_{\mathrm{WL},\mathrm{bias}}\,{\rm{\Delta }}{b}_{\mathrm{WL}\mathrm{mass}\mathrm{model},i}\\ & & +\ {\delta }_{i}\,{\rm{\Delta }}{b}_{\mathrm{measurement}\mathrm{systematics},i},\\ i & \in & \{\mathrm{Megacam},{\text{}}{HST}\},\end{array}\end{eqnarray} \tag{ 9 }$$

where b_{WL mass} is the mean bias due to WL mass modeling, Δb_{WL mass model} is the uncertainty on b_{WL mass}, and Δb_{measurement systematics} is the systematic measurement uncertainty due to multiplicative shear bias and uncertainties in the determination of the source redshift distribution; δ_WL,bias, δ_Megacam, and δ_HST are free parameters in our likelihood. With this parameterization, we apply Gaussian priors ${ \mathcal N }(0,1)$ on the three fit parameters. The numerical values of the different components of the WL bias are given in Table 1.

Table 1.  WL Modeling Parameters (S18; D19)

Effect	Parameter	Impact on Mass
		Megacam	HST
Intrinsic scatter	σintrinsic	0.214	0.26–0.42
Δ(Intrinsic scatter)	Δσintrinsic	0.04	0.021–0.055
Uncorrelated LSS scatter	σLSS	9 × 1013 M⊙	8 × 1013 M⊙
Δ(Uncorrelated LSS scatter)	ΔσLSS	1013 M⊙	1013 M⊙
WL mass bias	bWL mass	0.938	0.81–0.92
Mass modeling uncertainty	ΔbWL mass model	4.4%	5.8%–6.1%
Systematic measurement uncertainty	Δbmeasurement systematics	3.5%	7.2%
Total systematic uncertainty	N/A	5.6%	9.2%–9.4%

Note. The WL mass bias and the (lognormal) intrinsic scatter are calibrated against N-body simulations. Among other effects, they also account for the uncertainty and the scatter in the c(M) relation. This is done separately for each cluster in the HST sample leading to a range of values; here we report the smallest and largest individual values. The mass modeling uncertainty accounts for uncertainties in the calibration against N-body simulations and in the centering distribution. The systematic measurement uncertainties account for a multiplicative shear bias and the uncertainty in estimating the redshift distribution of source galaxies. Uncorrelated large-scale structure along the line of sight leads to an additional, Gaussian scatter.

Download table as:  ASCIITypeset image

The width of the (lognormal) scatter that is intrinsic to fitting WL shear profiles against NFW profiles is

$$\begin{eqnarray}\begin{array}{rcl}{\sigma }_{\mathrm{WL},i} & = & {\sigma }_{\mathrm{intrinsic},i}+{\delta }_{\mathrm{WL},\mathrm{scatter}}\,{\rm{\Delta }}{\sigma }_{\mathrm{intrinsic},i},\\ i & \in & \{\mathrm{Megacam},{\text{}}\mathrm{HST}\},\end{array}\end{eqnarray} \tag{ 10 }$$

where σ_intrinsic and Δσ_intrinsic are the mean intrinsic scatter and the error on the mean, respectively (given in Table 1); δ_WL,scatter is a free parameter in our likelihood on which we apply a Gaussian prior ${ \mathcal N }(0,1)$.

Finally, the width of the (normal) scatter due to uncorrelated large-scale structure is

$$\begin{eqnarray}\begin{array}{rcl}{\sigma }_{\mathrm{WL},\mathrm{LSS},i} & = & {\sigma }_{\mathrm{LSS},i}+{\delta }_{\mathrm{WL},\mathrm{LSS},i}\,{\rm{\Delta }}{\sigma }_{\mathrm{LSS},i},\\ i & \in & \{\mathrm{Megacam},{\text{}}\mathrm{HST}\},\end{array}\end{eqnarray} \tag{ 11 }$$

with the mean scatter σ_LSS and the error on the mean Δσ_LSS given in Table 1, and where we apply a Gaussian prior ${ \mathcal N }(0,1)$ on the fit parameters ${\delta }_{\mathrm{WL},{\mathrm{LSS}}_{\mathrm{Megacam}}}$ and ${\delta }_{\mathrm{WL},{\mathrm{LSS}}_{{\text{}}\mathrm{HST}}}$.

For reference, the total systematic error in the WL calibration is 5.6% for the Megacam sample (D19) and 9.2%–9.4% for the HST sample (S18). Given the small sample size of 19 and 13 clusters, our WL mass calibration is still dominated by statistical errors.

The analysis pipeline used in this work evolved from the code originally used in a previous SPT analysis (Bocquet et al. 2015). Since then, we have updated it to the full 2500 deg² survey, included the handling of WL data and the ability to account for correlated scatter among all observables, and modified the X-ray analysis (see Section 3.2.2). The pipeline is written as a module for CosmoSIS (Zuntz et al. 2015) and was also used for other WL scaling relation studies of SPT-SZ clusters (D19; Stern et al. 2019).

We start from a multi-observable Poisson log-likelihood

$$\begin{eqnarray}\begin{array}{rcl}\mathrm{ln}{ \mathcal L }({\boldsymbol{p}}) & = & \displaystyle \sum _{i}\mathrm{ln}\displaystyle \frac{{dN}(\xi ,{Y}_{{\rm{X}}},{g}_{{\rm{t}}},z| {\boldsymbol{p}})}{d\xi {{dY}}_{{\rm{X}}}\ {{dg}}_{{\rm{t}}}{dz}}\left|{}_{{\xi }_{i},{Y}_{{{\rm{X}}}_{i}},{g}_{{{\rm{t}}}_{i}},{z}_{i}}\right.\\ & & -\ \iint \iint d\xi \,{{dY}}_{{\rm{X}}}{{dg}}_{{\rm{t}}}\,{dz}\\ & & \times \ \left[\displaystyle \frac{{dN}(\xi ,{Y}_{{\rm{X}}},{g}_{{\rm{t}}},z| {\boldsymbol{p}})}{d\xi {{dY}}_{{\rm{X}}}\ {{dg}}_{{\rm{t}}}{dz}}{{\rm{\Theta }}}_{{\rm{s}}}\right]\\ & & +\ \mathrm{const},\end{array}\end{eqnarray} \tag{ 12 }$$

where the sum runs over all clusters i in the sample, and Θ_s is the survey selection function; in our case Θ_s = Θ(ξ > 5, z > 0.25).

As discussed in Bocquet et al. (2015) and explicitly shown in their Appendix, we rewrite the first term in Equation (12) as $P({Y}_{{\rm{X}}},{g}_{{\rm{t}}}| {\xi }_{i},{z}_{i},{\boldsymbol{p}})\left|{}_{{Y}_{{{\rm{X}}}_{i}},{g}_{{{\rm{t}}}_{i}}}\right.\times \tfrac{{dN}(\xi ,z| {\boldsymbol{p}})}{d\xi {dz}}\left|{}_{{\xi }_{i},{z}_{i}}\right.$. The second term in Equation (12) represents the total number of clusters in the survey, which are selected in ξ and z (and without any selection based on the follow-up observables). Therefore, this term reduces to $\int d\xi {dz}{{\rm{\Theta }}}_{{\rm{s}}}{dN}(\xi ,z| {\boldsymbol{p}})/d\xi {dz}$. With these modifications, and after explicitly setting the survey selection, the likelihood function becomes

$$\begin{eqnarray}\begin{array}{rcl}\mathrm{ln}{ \mathcal L }({\boldsymbol{p}}) & = & \displaystyle \sum _{i}\mathrm{ln}\displaystyle \frac{{dN}(\xi ,z| {\boldsymbol{p}})}{d\xi {dz}}\left|{}_{{\xi }_{i},{z}_{i}}\right.\\ & & -\ {\displaystyle \int }_{{z}_{\mathrm{cut}}}^{\infty }{dz}{\displaystyle \int }_{{\xi }_{\mathrm{cut}}}^{\infty }d\xi \displaystyle \frac{{dN}(\xi ,z| {\boldsymbol{p}})}{d\xi {dz}}\\ & & +\ \displaystyle \sum _{j}\mathrm{ln}P({Y}_{{\rm{X}}},{g}_{{\rm{t}}}| {\xi }_{j},{z}_{j},{\boldsymbol{p}})\left|{}_{{Y}_{{{\rm{X}}}_{j}},{g}_{{{\rm{t}}}_{j}}}\right.\end{array}\end{eqnarray} \tag{ 13 }$$

up to a constant. The first sum runs over all clusters i in the sample, and the second sum runs over all clusters j with ${Y}_{{\rm{X}}}$ and/or WL g_t measurements.

The first two terms in Equation (13) can be interpreted as the likelihood of the abundance (or number counts) of SZ clusters, while the third term represents the information from follow-up mass calibration. These two components are also visualized in the analysis flowchart in Figure 3: the number counts on the lower left side use the distribution of clusters in (ξ, z) space, and the mass calibration on the lower right also uses all available WL and X-ray follow-up data.

We note that the subsamples of clusters that were targeted for follow-up WL and/or X-ray data were selected at random within some cuts in ξ and redshift. Importantly, the selection was not made on WL and/or X-ray measurements. Therefore, the likelihood function presented above is complete; importantly, it does not suffer from biases from WL and/or X-ray selections.

3.2.1. Implementation of the Likelihood Function

We compute the individual terms in Equation (13) as follows:

$$\begin{eqnarray}\begin{array}{rcl}\displaystyle \frac{{dN}(\xi ,z| \,{\boldsymbol{p}})}{d\xi {dz}} & = & \iint {dM}\,d\zeta \left[\right.P(\xi | \zeta )P(\zeta | M,z,{\boldsymbol{p}})\\ & & \times \ \displaystyle \frac{{dN}(M,z| {\boldsymbol{p}})}{{dMdz}}{\rm{\Omega }}(z,{\boldsymbol{p}})\left.\right],\end{array}\end{eqnarray} \tag{ 14 }$$

where ${\rm{\Omega }}(z,{\boldsymbol{p}})$ is the survey volume and ${dN}(M,z| {\boldsymbol{p}})/{dMdz}$ is the HMF. We evaluate Equation (14) in the space (ξ, z) by convolving the HMF with the intrinsic scatter in $P(\zeta | M,z,{\boldsymbol{p}})$ and the measurement uncertainty in $P(\xi | \zeta )$.

The first term in Equation (13) is computed by evaluating Equation (14) at each cluster's measured (ξ_i, z_i), marginalizing over photometric redshift errors where present. The second term is a simple two-dimensional integral over Equation (14).

Our cluster sample contains 22 SZ detections for which no optical counterparts were found; these were assigned lower redshift limits z_lim in Bleem et al. (2015). We used simulations to determine the expected false-detection rate dN_false(ξ)/dξ given survey specifics (see Section 2.2 and Table 1 in dH16). For each unconfirmed cluster candidate, we evaluate a modified version of the first term in Equation (13),

$$\begin{eqnarray}\begin{array}{rcl}\displaystyle \frac{{{dN}}_{\mathrm{unconf}.\mathrm{cand}.}(\xi ,z| {\boldsymbol{p}})}{d\xi {dz}} & = & \displaystyle \frac{{{dN}}_{\mathrm{cluster}}(\xi ,z| {\boldsymbol{p}})}{d\xi {dz}}\\ & & +\ \displaystyle \frac{{{dN}}_{\mathrm{false}}(\xi )}{d\xi },\end{array}\end{eqnarray} \tag{ 15 }$$

and marginalize over the candidate's allowed redshift range ${z}_{\mathrm{lim}}\lt z\lt \infty$. Note that the total expected number of false detections $\int d\xi {{dN}}_{\mathrm{false}}(\xi )/d\xi$ is independent of ${\boldsymbol{p}}$ and is therefore neglected in Equation (13). The expected number of false detections in the SPT-SZ survey is 18 ± 4, which is consistent with our 22 unconfirmed candidates (dH16). In practice, we obtain essentially unchanged results if we simply discard the 22 optically unconfirmed SZ detections from the catalog. There are nine clusters that are detected in the overlap region between adjacent SPT fields. We follow dH16 and double-count these clusters in our analysis. Accounting for only one object of each pair of these clusters instead does not change our results in any significant way.

The mass calibration term in Equation (13) is computed as

$$\begin{eqnarray}&&\begin{array}{l}P({Y}_{{\rm{X}}}^{\mathrm{obs}},{g}_{{\rm{t}}}^{\mathrm{obs}}| \xi ,z,{\boldsymbol{p}})\\ =\ \int \int \int \int {dM}\,d\zeta \,{{dY}}_{{\rm{X}}}{{dM}}_{\mathrm{WL}}\\ \ \ \times \ \left[P({Y}_{{\rm{X}}}^{\mathrm{obs}}| {Y}_{{\rm{X}}})P({g}_{{\rm{t}}}^{\mathrm{obs}}| {M}_{\mathrm{WL}})P(\xi | \zeta )\right.\\ \ \ \times \ \left.P(\zeta ,{Y}_{{\rm{X}}},{M}_{\mathrm{WL}}| M,z,{\boldsymbol{p}})P(M| z,{\boldsymbol{p}})\right],\end{array}\end{eqnarray} \tag{ 16 }$$

with the HMF $P(M| z,{\boldsymbol{p}})$ and the multi-observable scaling relation $P(\zeta ,{Y}_{{\rm{X}}},{M}_{\mathrm{WL}}| M,z,{\boldsymbol{p}})$ that includes the effects of correlated scatter. Computing this multidimensional integral in the (ζ, ${Y}_{{\rm{X}}}$, ${M}_{\mathrm{WL}}$) space is expensive. We minimize the computational cost of this step by (i) only considering parts of the (ζ, ${Y}_{{\rm{X}}}$, ${M}_{\mathrm{WL}}$) space that have non-negligible probability densities (we estimate this subspace from the measurements and ${\boldsymbol{p}}$), (ii) using fast Fourier transform convolutions, and (iii) only performing this computation for clusters that actually have both follow-up measurements ${Y}_{{\rm{X}}}$ and ${M}_{\mathrm{WL}}$; otherwise, we restrict the computation to the much cheaper two-dimensional (${Y}_{{\rm{X}}}$, ζ) or (${M}_{\mathrm{WL}}$, ζ) spaces. The mass calibration term does not need to be computed at all for clusters that have no X-ray or WL follow-up data.

3.2.2. Update of the X-Ray Analysis Scheme

The X-ray observable is a measurement of the radial ${Y}_{{\rm{X}}}$ profile. The scaling relation, on the other hand, predicts a value of the observable integrated out to r₅₀₀ for a given M₅₀₀. In a self-consistent analysis, the likelihood should be extracted by comparing the data and the model prediction at the same radius.

In previous SPT analyses, a ${Y}_{{\rm{X}}}$ value was extracted from the profile by iteratively solving for the radius r_iter at which the measured ${Y}_{{\rm{X}}}$ and the X-ray scaling relation prediction from Equation (2) match (the scaling relation is evaluated at ${M}_{500}\equiv 4\pi /3{r}_{\mathrm{iter}}^{3}500{\rho }_{c}$). This iteration was repeated for each set of parameters ${\boldsymbol{p}}$, but within a fixed reference cosmology. However, this method introduces a bias, because r_iter is not equal to the radius r₅₀₀ at which the scaling relation $P(\zeta ,{Y}_{{\rm{X}}},{M}_{\mathrm{WL}}| {M}_{500},z,{\boldsymbol{p}})$ in Equation (16) is evaluated.

We choose a different approach and evaluate both the (integrated) measured profile and the model prediction at a fixed fiducial radius r_fid. We define r_fid for each cluster by computing r_500,fid from its SZ significance ξ using a fiducial set of SZ scaling relation parameters and setting r_fid = r_500,fid. Then, for each set of parameters ${\boldsymbol{p}}$ in the analysis, we convert the model prediction ${Y}_{{\rm{X}}}$(r₅₀₀) from radius r₅₀₀ to r_fid. We use the fact that the radial profiles are well approximated by power laws in radius

$$\begin{eqnarray}&&\displaystyle \frac{{Y}_{{\rm{X}}}(r)}{{Y}_{{\rm{X}}}({r}_{500})}={\left(\displaystyle \frac{r}{{r}_{500}}\right)}^{d\mathrm{ln}{Y}_{{\rm{X}}}/d\mathrm{ln}r},\end{eqnarray} \tag{ 17 }$$

where r₅₀₀ is derived from M_500c. In our analysis, we assume isothermality (see Section 2.2), and so d ln ${Y}_{{\rm{X}}}$/d ln r equals the radial slope in gas mass d ln M_g/d ln r. From our sample we measure

$$\begin{eqnarray}&&d\mathrm{ln}{M}_{{\rm{g}}}/d\mathrm{ln}r=1.12\pm 0.23.\end{eqnarray} \tag{ 18 }$$

We are now able to make a model prediction at r_fid, starting from the scaling relation prediction ${Y}_{{\rm{X}}}$(r₅₀₀):

$$\begin{eqnarray}&&{Y}_{{\rm{X}}}({r}_{\mathrm{fid}})={Y}_{{\rm{X}}}({r}_{500}){\left(\displaystyle \frac{{r}_{\mathrm{fid}}}{{r}_{500}}\right)}^{d\mathrm{ln}{M}_{{\rm{g}}}/d\mathrm{ln}r}.\end{eqnarray} \tag{ 19 }$$

In the analysis, we marginalize over the uncertainty in d ln M_g/d ln r, which shows negligible correlation with any other parameter. Note that this prescription for the model prediction ${Y}_{{\rm{X}}}$(r_fid) contains an additional dependence on r₅₀₀ and thus on M₅₀₀.

We note that a similar approach was adopted by other groups (e.g., Mantz et al. 2010a, 2015). We have shown through tests against mock catalogs that the new analysis scheme is unbiased and that the previous method biased ${B}_{{Y}_{{\rm{X}}}}$ low at a level that is comparable to the uncertainty on that parameter, while the effect on other parameters was very small.

We assume the HMF fit by Tinker et al. (2008). This approach assumes universality of the HMF across the cosmological parameter space considered in this work and uses a fitting function that was calibrated against N-body simulations. In principle, the HMF is also affected by baryonic effects. However, hydrodynamic simulations suggest that these have negligible impact for clusters with masses as high as those considered here (Velliscig et al. 2014); this was explicitly tested for a simulated and idealized SPT-SZ cluster survey (Bocquet et al. 2016). Finally, note that the Tinker et al. (2008) fit applies to mean spherical overdensities in the range 200 ≤ Δ_mean ≤ 3200, and we thus convert to Δ_500crit using Δ_mean(z) = 500/${{\rm{\Omega }}}_{{\rm{m}}}$(z). As the HMF fit is only calibrated up to Δ_mean = 3200, we require ${{\rm{\Omega }}}_{{\rm{m}}}$(z) ≥ 500/3200 = 0.15625 for all redshifts z ≥ 0.25 relevant for our cluster sample.

We have run extensive tests to ensure that our analysis pipeline is unbiased at a level that is much smaller than our total error budget. The primary approach is testing against mock catalogs. Of course, such tests are only useful if producing mocks is easier and more reliable than the actual analysis. In our case, the analysis is challenging mainly because of the computation of multidimensional integrals. To create one of our mocks, on the other hand, one has to compute the HMF, apply the observable–mass relations, draw random deviates, and compute WL shear profiles. Using the same code to compute the HMF for the mocks and the analysis would undercut the usefulness of the testing, and so we also created mocks using HMFs computed with an independent code. For the same reason, the mock shear profiles were created using an independent code. We typically create mock catalogs that contain an order of magnitude more clusters and calibration data than our real sample. We created and analyzed sets of mocks using different random seeds and different sets of input parameters (notably, some with w $\ne$ −1). No test indicated any biases in our analysis pipeline at the level relevant for our data set.

We characterize the agreement between constraints obtained from pairs of probes (e.g., clusters and primary CMB anisotropies) by quantifying whether the difference between the two posterior distributions is consistent with zero difference. We draw representative samples $[{{\boldsymbol{x}}}_{1}]$ and $[{{\boldsymbol{x}}}_{2}]$ from the posteriors of the two probes ${P}_{1}({\boldsymbol{x}})$ and ${P}_{2}({\boldsymbol{x}})$, compute the difference between all pairs of points ${\boldsymbol{\delta }}\equiv {{\boldsymbol{x}}}_{1}-{{\boldsymbol{x}}}_{2}$, and then construct the probability distribution D from the ensemble $[{\boldsymbol{\delta }}]$. The probability value (or p-value) that the two distributions represent the same underlying quantity is

$$\begin{eqnarray}&&p={\int }_{D\lt D(0)}d{\boldsymbol{\delta }}D({\boldsymbol{\delta }}),\end{eqnarray} \tag{ 20 }$$

where $D(0)$ is the probability of zero difference. The p-value can be converted into a significance assuming Gaussian statistics. This measure can be applied to one-dimensional and multidimensional parameter spaces. The code is publicly available.⁵⁴

In our cosmological fits, we assume spatial flatness and allow the sum of neutrino masses to vary. The comparison of our results with constraints from primary CMB anisotropies is of prime interest—notably, the comparison of constraints on ${\sigma }_{8}$. For primary CMB anisotropies, ${\sigma }_{8}$ is strongly degenerate with $\sum {m}_{\nu }$, and so the latter should be a free parameter of the model to avoid artificially tight constraints. We refer to the flat ΛCDM model with a varying sum of neutrino masses as νΛCDM and to its extension with a free dark energy equation-of-state parameter as νwCDM.

In the νΛCDM cosmology, we vary the cosmological parameters ${{\rm{\Omega }}}_{{\rm{m}}}$, ${{\rm{\Omega }}}_{\nu }{h}^{2}$, ${{\rm{\Omega }}}_{{\rm{b}}}{h}^{2}$, A_s, h, n_s; ${\sigma }_{8}$ is a derived parameter. Our cluster data primarily constrain ${{\rm{\Omega }}}_{{\rm{m}}}$ and ${\sigma }_{8}$, and we marginalize over flat priors on the other parameters. The parameter ranges for ${{\rm{\Omega }}}_{{\rm{b}}}{h}^{2}$ and n_s are chosen to roughly match the 5σ credibility interval of the Planck constraints; h is allowed to vary in the range 0.55–0.9. We assume two massless neutrinos and one massive neutrino and allow ${{\rm{\Omega }}}_{\nu }{h}^{2}$ to vary in the range 0–0.01; this corresponds to a range in $\sum {m}_{\nu }$ of 0–0.93 eV. We note that the minimum allowed sum of neutrino masses from oscillation experiments is $\sum {m}_{\nu }\gt 59.5\pm 0.5$ meV (Tanabashi et al. 2018). In a departure from previous SPT analyses, we do not apply a BBN prior on ${{\rm{\Omega }}}_{{\rm{b}}}{h}^{2}$ or constraints from direct measurements of H₀. We remind the reader that the implementation of the theory HMF leads to an effective, hard prior ${{\rm{\Omega }}}_{{\rm{m}}}$(z) ≳ 0.16 for all redshifts z > 0.25 relevant to our survey (see Section 3.3); however, this prior does not affect our results. All parameters and their priors are summarized in Table 2.

Table 2.  Summary of Cosmological and Astrophysical Parameters Used in Our Fiducial Analysis

Parameter	Prior
Cosmological
${{\rm{\Omega }}}_{{\rm{m}}}$	${ \mathcal U }(0.05,0.6),{{\rm{\Omega }}}_{{\rm{m}}}(z\gt 0.25)\gt 0.156$
${{\rm{\Omega }}}_{{\rm{b}}}{h}^{2}$	${ \mathcal U }(0.020,0.024)$
${{\rm{\Omega }}}_{\nu }{h}^{2}$	${ \mathcal U }(0,0.01)$
Ωk	Fixed (0)
As	${ \mathcal U }({10}^{-10},{10}^{-8})$
h	${ \mathcal U }(0.55,0.9)$
ns	${ \mathcal U }(0.94,1.00)$
w	Fixed (−1) or ${ \mathcal U }(-2.5,-0.33)$
Optical Depth to Reionization
τ	Fixed or ${ \mathcal U }(0.02,0.14)$
SZ Scaling Relation
${A}_{\mathrm{SZ}}$	${ \mathcal U }(1,10)$
${B}_{\mathrm{SZ}}$	${ \mathcal U }(1,2.5)$
${C}_{\mathrm{SZ}}$	${ \mathcal U }(-1,2)$
${\sigma }_{\mathrm{ln}\zeta }$	${ \mathcal U }(0.01,0.5)$ ($\times { \mathcal N }(0.13,{0.13}^{2}))$
X-ray ${Y}_{{\rm{X}}}$ Scaling Relation
${A}_{{Y}_{{\rm{X}}}}$	${ \mathcal U }(3,10)$
${B}_{{Y}_{{\rm{X}}}}$	${ \mathcal U }(0.3,0.9)$
${C}_{{Y}_{{\rm{X}}}}$	${ \mathcal U }(-1,0.5)$
${\sigma }_{\mathrm{ln}{Y}_{{\rm{X}}}}$	${ \mathcal U }(0.01,0.5)$
d ln Mg/d ln r	${ \mathcal U }(0.4,1.8)\times { \mathcal N }(1.12,{0.23}^{2})$
WL Modeling
δWL,bias	${ \mathcal U }(-3,3)\times { \mathcal N }(0,1)$
δMegacam	${ \mathcal U }(-3,3)\times { \mathcal N }(0,1)$
δHST	${ \mathcal U }(-3,3)\times { \mathcal N }(0,1)$
δWL,scatter	${ \mathcal U }(-3,3)\times { \mathcal N }(0,1)$
${\delta }_{\mathrm{WL},{\mathrm{LSS}}_{\mathrm{Megacam}}}$	${ \mathcal U }(-3,3)\times { \mathcal N }(0,1)$
${\delta }_{\mathrm{WL},{\mathrm{LSS}}_{{\text{}}\mathrm{HST}}}$	${ \mathcal U }(-3,3)\times { \mathcal N }(0,1)$
Correlated Scatter
ρSZ−WL	${ \mathcal U }(-1,1)$
ρSZ−X	${ \mathcal U }(-1,1)$
ρX−WL	${ \mathcal U }(-1,1)$
	$\det ({{\boldsymbol{\Sigma }}}_{\mathrm{multi} \mbox{-} \mathrm{obs}})\gt 0$

Note. The Gaussian prior on ${\sigma }_{\mathrm{ln}\zeta }$ is only applied when no X-ray data are included in the fit. The parameter ranges for ${{\rm{\Omega }}}_{{\rm{b}}}{h}^{2}$ and n_s are chosen to roughly match the 5σ interval of the Planck ΛCDM results. w is fixed to −1 for ΛCDM and is allowed to vary for $w\mathrm{CDM}$. The optical depth to reionization τ is only relevant when Planck data are included in the analysis. The WL modeling systematics are presented in Table 1.

Download table as:  ASCIITypeset image

The likelihood sampling is done within CosmoSIS using the Metropolis (Metropolis et al. 1953) and MultiNest (Feroz et al. 2009) samplers. We confirmed that they produce consistent results.

From the cluster abundance measurement of our SPTcl (SPT-SZ +WL+${Y}_{{\rm{X}}}$) data set we obtain our baseline results

$$\begin{eqnarray}&&{{\rm{\Omega }}}_{{\rm{m}}}=0.276\,\pm \,0.047\end{eqnarray} \tag{ 21 }$$

$$\begin{eqnarray}&&{\sigma }_{8}=0.781\,\pm \,0.037\end{eqnarray} \tag{ 22 }$$

$$\begin{eqnarray}&&{\sigma }_{8}{\left({{\rm{\Omega }}}_{{\rm{m}}}/0.3\right)}^{0.2}=0.766\,\pm \,0.025.\end{eqnarray} \tag{ 23 }$$

The remaining cosmological parameters (including $\sum {m}_{\nu }$; see Figure 9) are not or only weakly constrained by the cluster data. Constraints on scaling relation parameters can be found in Table 3. We note that applying priors on ${{\rm{\Omega }}}_{{\rm{b}}}{h}^{2}$ and H₀ from BBN and direct measurements of H₀ and/or fixing the sum of neutrino masses to 0.06 eV, approximately the lower limit predicted from terrestrial oscillation experiments, does not affect our constraints on ${{\rm{\Omega }}}_{{\rm{m}}}$ and ${\sigma }_{8}$ in any significant way (see Figure 15 in Appendix A for the impact of fixing the sum of the neutrino masses).

4.1.1. Goodness of Fit

In Figure 4, we compare the measured distribution of clusters as a function of their redshift and SPT detection significance with the model prediction evaluated for the recovered parameter constraints. This figure does not suggest any problematic feature in the data.

Zoom In Zoom Out Reset image size

For a more quantitative discussion, we bin our confirmed clusters into a grid of 30 × 30 in redshift and detection significance, and we confront this measurement with the expected number of objects in each two-dimensional bin. The expected (and measured) numbers in each bin are too small to apply Gaussian χ² statistics, and we estimate the goodness of fit using a prescription for the Poisson statistic (Kaastra 2017).⁵⁵ This approach is similar to our likelihood analysis, which applies Poisson statistics within infinitesimally small bins, instead of the larger bins we assume here. Adopting the maximum-posterior νΛCDM parameters, we compute the expected number of clusters in each of the 30 × 30 bins and follow Kaastra (2017) to evaluate the test statistic C. We obtain an expected mean C_e and variance C_v,

$$\begin{eqnarray}&&{C}_{{\rm{e}}}=439.8;\,{C}_{{\rm{v}}}={26.8}^{2}.\end{eqnarray} \tag{ 24 }$$

For samples that contain at least a few hundred objects—like ours—the statistic C is well approximated by a Gaussian with mean C_e and variance C_v (Kaastra 2017). The data statistic for our sample is

$$\begin{eqnarray}&&{C}_{{\rm{d}}}=449.3,\end{eqnarray} \tag{ 25 }$$

in full agreement with the range expected for C_e, indicating that the model provides an adequate fit to the data.

4.1.2. Comparison with Previous SPT Results

As discussed in the Introduction, this work uses the same SPT-SZ cluster sample (Bleem et al. 2015, now with updated photometric redshifts; see Section 2.1) that was analyzed in dH16, and the key update is the inclusion of WL data. In dH16, the amplitude of the observable–mass relation was set by a prior on the X-ray normalization ${A}_{{Y}_{{\rm{X}}}}$, which in turn was informed by external WL data sets (CCCP and WtG; Applegate et al. 2014; von der Linden et al. 2014; Hoekstra et al. 2015). Gaussian priors were applied to the remaining SZ and X-ray scaling relation parameters, which we dropped for this analysis. In Figure 5, we compare our constraints on ${{\rm{\Omega }}}_{{\rm{m}}}$–${\sigma }_{8}$ with the ones presented in dH16. We recover very similar results; in ${{\rm{\Omega }}}_{{\rm{m}}}$–${\sigma }_{8}$ space, the agreement is p = 0.86 (0.2σ). Since the key difference between dH16 and this work is the inclusion of WL data, this agreement indirectly confirms that our internal WL mass calibration agrees with the external priors adopted previously. This is expected because the X-ray prior adopted in previous work agrees well with the measurement enabled by our own WL data set (D19).

Zoom In Zoom Out Reset image size

4.1.3. Comparison with External Probes

In Figure 6, we show a comparison of our results with constraints from Planck (TT+lowTEB) and from combined analyses of cosmic shear, galaxy–galaxy lensing, and galaxy clustering from the Kilo Degree Survey and the Galaxies And Mass Assembly survey (KiDS+GAMA; van Uitert et al. 2018) and the Dark Energy Survey (DES) Year 1 results (Abbott et al. 2018). We also compare our results with another cluster study that used internal WL mass calibration but a sample based on X-ray selection (Weighing the Giants, or WtG; Mantz et al. 2015). Overall, the constraining power of all probes is roughly similar in this plane. There is good agreement among all probes, as the 68% contours all overlap. In particular, the cluster-based constraints yield very similar ${{\rm{\Omega }}}_{{\rm{m}}}$, but WtG favors a somewhat higher ${\sigma }_{8}$. Interestingly, the degeneracy axis of WtG is slightly tilted with respect to SPTcl, which we attribute to the different redshift and mass ranges spanned by the two samples.

Zoom In Zoom Out Reset image size

We pay particular attention to a comparison with Planck (TT+lowTEB). Our constraint on ${\sigma }_{8}{\left({{\rm{\Omega }}}_{{\rm{m}}}/0.3\right)}^{0.2}$ = 0.766 ± 0.025 is lower than the one from Planck (${\sigma }_{8}{\left({{\rm{\Omega }}}_{{\rm{m}}}/0.3\right)}^{0.2}$ = ${0.814}_{-0.020}^{+0.041}$); the agreement between the two measurements is p = 0.28 (1.1σ). In the two-dimensional ${{\rm{\Omega }}}_{{\rm{m}}}$–${\sigma }_{8}$ space, the agreement is p = 0.13 (1.5σ).

We note that the latest analysis of the cluster sample selected by the Planck satellite is qualitatively in agreement with our constraint, as shown in Figure 32 of Planck Collaboration et al. (2018b). Notably, the 95% contour of their result, calibrated using CMB lensing, encompasses the Planck primary CMB result in the ${{\rm{\Omega }}}_{{\rm{m}}}$–${\sigma }_{8}$ plane.

4.1.4. Impact of X-Ray Follow-up Data

We compare our baseline results from SPTcl (SPT-SZ+WL+${Y}_{{\rm{X}}}$) with the ones obtained from the SPT-SZ+WL data combination, in which no X-ray follow-up data are included. In this case, we apply an informative Gaussian prior to the SZ scatter ${\sigma }_{\mathrm{ln}\zeta }$. As in all of this work, no informative priors are applied on the remaining three SZ scaling relation parameters and on the X-ray scaling relation parameters. A figure showing constraints on all relevant parameters can be found in Appendix B (Figure 16, compare blue and red contours), and Table 3 summarizes parameter constraints. Both data combinations, with and without X-ray data, provide very similar constraints on cosmological and scaling relation parameters. Without informative priors on the X-ray amplitude, mass slope, or redshift evolution, the inclusion of X-ray data does not enable tighter constraints. The use of X-ray data does, however, enable constraints on the SZ and X-ray scatters ${\sigma }_{\mathrm{ln}\zeta }$ and ${\sigma }_{\mathrm{ln}{Y}_{{\rm{X}}}}$, with flat priors applied to both.

Note that our data set is sensitive to the SZ-to-X-ray relation. As defined in Section 3.1, our model consists of two observable–mass relations that each relate one observable to mass. This implies that the amplitudes, mass slopes, and redshift evolutions of the two scaling relations are degenerate, as shown in Figure 7. The degeneracy between ${\sigma }_{\mathrm{ln}\zeta }$ and ${\sigma }_{\mathrm{ln}{Y}_{{\rm{X}}}}$ is particularly interesting: while the marginalized posterior of either or both parameters has substantial mass near 0 scatter (see Figure 16), the bottom right panel of Figure 7 shows that 0 total scatter is clearly ruled out.

Zoom In Zoom Out Reset image size

Our data set is not able to constrain any of the coefficients describing the correlated scatter among the observables. The visual impression of a constraint in Figure 16 stems from the requirement that the matrix describing the multi-observable scatter must be a valid non-degenerate covariance matrix that prevents combinations of extreme correlation coefficients.

4.1.5. Constraints on X-Ray Scaling Relation Parameters

Without any informative priors on the X-ray scaling relation parameters, we can use the SPTcl data set to constrain the ${Y}_{{\rm{X}}}$–mass relation. The recovered amplitude

$$\begin{eqnarray}&&{A}_{{Y}_{{\rm{X}}}}=6.35\,\pm \,0.69\end{eqnarray} \tag{ 26 }$$

is very close to the WL-informed prior (Applegate et al. 2014; von der Linden et al. 2014; Hoekstra et al. 2015; Mantz et al. 2015) that was used in our previous cosmology analysis (${A}_{{Y}_{{\rm{X}}}}$ = 6.38 ± 0.61; dH16). We constrain the redshift evolution of the ${Y}_{{\rm{X}}}$–mass relation to

$$\begin{eqnarray}&&{C}_{{Y}_{{\rm{X}}}}=-{0.31}_{-0.21}^{+0.14}.\end{eqnarray} \tag{ 27 }$$

The self-similar expectation ${C}_{{Y}_{{\rm{X}}}}$ = −0.4 is well within 1σ. Our measurement of the ${Y}_{{\rm{X}}}$ scatter

$$\begin{eqnarray}&&{\sigma }_{\mathrm{ln}{Y}_{{\rm{X}}}}=0.18\,\pm \,0.09\end{eqnarray} \tag{ 28 }$$

is higher than but consistent at the 1σ level with the prior 0.12 ± 0.08 adopted in previous SPT analyses. It closely matches the measurement 0.182 ± 0.015 from Mantz et al. (2016), although with larger uncertainty.

The recovered ${Y}_{{\rm{X}}}$ mass slope

$$\begin{eqnarray}&&{B}_{{Y}_{{\rm{X}}}}=0.514\,\pm \,0.037\end{eqnarray} \tag{ 29 }$$

is lower than the self-similar evolution ${B}_{{Y}_{{\rm{X}}}}$ = 0.6 and the measurements ${B}_{{Y}_{{\rm{X}}}}$ = 0.57 ± 0.03 from Vikhlinin et al. (2009a) and ${B}_{{Y}_{{\rm{X}}}}$ = 1/(1.61 ± 0.04) = 0.621 ± 0.015 from Mantz et al. (2016).⁵⁶ From our data, the consistency of ${B}_{{Y}_{{\rm{X}}}}$ with the self-similar value is p = 0.021, corresponding to 2.3σ. Our data constrain ${B}_{{Y}_{{\rm{X}}}}$ through its degeneracy with the SZ mass slope ${B}_{\mathrm{SZ}}$ (Figure 7), which in turn is constrained through the process of fitting the cluster abundance against the HMF. This subject was already discussed in dH16, where a prior on ${B}_{{Y}_{{\rm{X}}}}$ was adopted from the measurement by Vikhlinin et al. (2009a).

As a cross-check, and because other groups have used the X-ray gas mass as their low-scatter mass proxy, we repeat the analysis replacing the ${Y}_{{\rm{X}}}$ data with ${M}_{\mathrm{gas}}$ measurements. We apply no informative priors on the four parameters of the ${M}_{\mathrm{gas}}$ scaling relation of Equation (7). We then analyze this SPT-SZ+WL+${M}_{\mathrm{gas}}$ data set. The constraints on the SZ scaling relation parameters and cosmology are very similar to the results from the fiducial SPT-SZ+WL+${Y}_{{\rm{X}}}$ analysis, and again we observe an X-ray mass slope that disagrees with the self-similar evolution. We measure

$$\begin{eqnarray}&&{A}_{{M}_{{\rm{g}}}}=0.116\,\pm \,0.011\ \end{eqnarray} \tag{ 30 }$$

$$\begin{eqnarray}&&{B}_{{M}_{{\rm{g}}}}=1.22\,\pm \,0.07\end{eqnarray} \tag{ 31 }$$

$$\begin{eqnarray}&&{C}_{{M}_{{\rm{g}}}}=-0.05\,\pm \,0.17\end{eqnarray} \tag{ 32 }$$

$$\begin{eqnarray}&&{\sigma }_{\mathrm{ln}{M}_{{\rm{g}}}}=0.11\,\pm \,0.04.\end{eqnarray} \tag{ 33 }$$

This corresponds to a 2.5σ preference for a slope that is steeper than the self-similar expectation ${B}_{{M}_{{\rm{g}}}}=1$ or the measurement ${B}_{{M}_{{\rm{g}}}}=1.004\pm 0.014$ from Mantz et al. (2016). The measurement ${B}_{{M}_{{\rm{g}}}}=1.15\pm 0.02$ ⁵⁷ from Vikhlinin et al. (2009a) is in between the two results and is 1σ low compared to ours. Because these slopes differ, we compare the measurements of the gas fraction ${A}_{{M}_{{\rm{g}}}}$ at the pivot mass in our relation 5 × 10¹⁴ M_⊙/h₇₀, where we obtain Equation (30). The mean gas fraction at this mass is 0.128 from Mantz et al. (2016) and 0.114 from Vikhlinin et al. (2009a). Both values are contained within the 1σ range of our measurement. Finally, as for ${Y}_{{\rm{X}}}$, our measurement of the redshift evolution encompasses the self-similar evolution (${C}_{{M}_{{\rm{g}}}}$ = 0) within 1σ.

For an extensive discussion of the mass and redshift trends in the ${M}_{\mathrm{gas}}$–mass and ${Y}_{{\rm{X}}}$–mass relations for SPT-selected clusters and how they compare to previously published results, we refer the reader to two recent studies where SZ-based mass information was adopted using the posterior distributions of the SZ ζ–mass relation parameters presented in dH16 (Chiu et al. 2018; Bulbul et al. 2019). Bulbul et al. (2019) used X-ray data from XMM-Newton, while we use data from Chandra; their recovered constraints on the X-ray mass slopes and redshift evolutions are consistent with our findings at the 1σ level, which confirms a consistent X-ray analysis. Here we note that most measurements of X-ray scaling relations have been performed using samples at low redshifts z ≲ 0.5, and so it is of particular interest to examine the mass slopes for the low-redshift half of our sample.

We therefore split our cluster sample (and all follow-up data) into two subsamples above and below z = 0.6, the median redshift of our sample. Constraints on the most relevant parameters are shown in Figure 16 in Appendix B, and Figure 8 shows the constraints on ${B}_{{Y}_{{\rm{X}}}}$. Interestingly, the low-redshift subsample prefers a higher value

$$\begin{eqnarray}&&{B}_{{Y}_{{\rm{X}}}}(0.25\lt z\lt 0.6)={0.583}_{-0.069}^{+0.054}\end{eqnarray} \tag{ 34 }$$

that is closer to the self-similar evolution ${B}_{{Y}_{{\rm{X}}}}$ = 0.6. As expected, the value obtained from the high-z subsample,

$$\begin{eqnarray}&&{B}_{{Y}_{{\rm{X}}}}(z\gt 0.6)={0.503}_{-0.047}^{+0.037},\end{eqnarray} \tag{ 35 }$$

is lower than the one obtained from the full sample. However, note that the low-redshift and high-redshift constraints on ${B}_{{Y}_{{\rm{X}}}}$ only differ with p = 0.44 (0.8σ).

Zoom In Zoom Out Reset image size

We perform the same splits in redshift using the SPT-SZ+WL+${M}_{\mathrm{gas}}$ data set. Here as well, our measurement using the low-redshift subsample,

$$\begin{eqnarray}&&{B}_{{M}_{{\rm{g}}}}(0.25\lt z\lt 0.6)=1.12\,\pm \,0.09,\end{eqnarray} \tag{ 36 }$$

is closer to the self-similar evolution, while the high-redshift half yields a steeper slope

$$\begin{eqnarray}&&{B}_{{M}_{{\rm{g}}}}(z\gt 0.6)=1.36\,\pm \,0.11.\end{eqnarray} \tag{ 37 }$$

To capture a possible redshift dependence of the slope of the X-ray scaling relations, we analyzed models with an extended scaling relation model of the form

$$\begin{eqnarray}\begin{array}{rcl}\mathrm{ln}{{ \mathcal O }}_{{\rm{X}} \mbox{-} \mathrm{ray}} & = & \mathrm{ln}A+B\mathrm{ln}\left(\displaystyle \frac{{M}_{500c}\,{h}_{70}}{5\times {10}^{14}{M}_{\odot }}\right)\\ & & +\ C\mathrm{ln}\left(\displaystyle \frac{E(z)}{E(0.6)}\right)\\ & & +\ E\mathrm{ln}\left(\displaystyle \frac{E(z)}{E(0.6)}\right)\mathrm{ln}\left(\displaystyle \frac{{M}_{500c}\,{h}_{70}}{5\times {10}^{14}{M}_{\odot }}\right)\end{array}\end{eqnarray} \tag{ 38 }$$

that allows for additional freedom and the mass and redshift dependences. However, we do not observe any significant departure in E from 0, in agreement with Bulbul et al. (2019).

Having quantified the consistency between our cluster data set and Planck in Section 4.1.3, we proceed and combine the two probes. The SPTcl+Planck data set yields

$$\begin{eqnarray}&&{{\rm{\Omega }}}_{{\rm{m}}}=0.353\,\pm \,0.027\end{eqnarray} \tag{ 39 }$$

$$\begin{eqnarray}&&{\sigma }_{8}=0.761\,\pm \,0.033\end{eqnarray} \tag{ 40 }$$

$$\begin{eqnarray}&&{\sigma }_{8}{\left({{\rm{\Omega }}}_{{\rm{m}}}/0.3\right)}^{0.2}=0.786\,\pm \,0.025\end{eqnarray} \tag{ 41 }$$

$$\begin{eqnarray}&&\sum {m}_{\nu }=0.39\,\pm \,0.19\,\mathrm{eV}\end{eqnarray} \tag{ 42 }$$

$$\begin{eqnarray}&&\sum {m}_{\nu }\ \lt 0.74\,\mathrm{eV}\,(95 \% \,\mathrm{upper}\ \mathrm{limit}).\end{eqnarray} \tag{ 43 }$$

Compared to constraints from Planck alone, the combination with SPTcl shrinks the errors on ${{\rm{\Omega }}}_{{\rm{m}}}$, ${\sigma }_{8}$, and ${\sigma }_{8}{\left({{\rm{\Omega }}}_{{\rm{m}}}/0.3\right)}^{0.2}$ by 3%, 12%, and 20%. By breaking parameter degeneracies (notably between ${\sigma }_{8}$ and $\sum {m}_{\nu }$; see Figure 9), the addition of cluster data to the primary CMB measurements by Planck affects the inferred sum of neutrino masses. If interpreted as a Gaussian probability distribution (i.e., ignoring the hard cut $\sum {m}_{\nu }\gt 0$), our joint measurement corresponds to a 2.0σ preference for a nonzero sum of neutrino masses.

Zoom In Zoom Out Reset image size

The Planck Collaboration recently presented an updated analysis of primary CMB anisotropies (Planck Collaboration et al. 2018a). Most notably, the optical depth decreased to τ = 0.054 ± 0.007. As the updated Planck likelihood code is not available yet, we estimate the impact of the updated Planck analysis on our results and especially our constraint on $\sum {m}_{\nu }$ by analyzing the Planck 2015 TT data (without lowTEB) with a prior on $\tau \sim { \mathcal N }(0.054,{0.007}^{2})$. We analyze the joint SPTcl+Planck TT+τ prior data set and obtain

$$\begin{eqnarray}&&\sum {m}_{\nu }=0.35\,\pm \,0.21\,\mathrm{eV}.\end{eqnarray} \tag{ 44 }$$

The recovered constraint is lower than our fiducial constraint using the (SPTcl+Planck TT+lowTEB) data set, and the 95% credible interval runs against the hard prior $\sum {m}_{\nu }=0$. The preference for a nonzero sum of neutrino masses reduces to 1.7σ. We caution the reader that this result is only preliminary owing to the way it depends on the prior on τ that we adopted. The full analysis will require analyzing our cluster sample jointly with the latest Planck analysis.

We explain the shift in $\sum {m}_{\nu }$ toward lower values as follows. In ΛCDM, the relationship between A_s and ${\sigma }_{8}$ is essentially fixed. However, in νΛCDM, the additional degree of freedom $\sum {m}_{\nu }$ allows for different values of ${\sigma }_{8}$ at a fixed A_s. In any joint analysis of Planck+low-redshift growth-of-structure probe such as SPTcl, $\sum {m}_{\nu }$ is constrained to accommodate the Planck measurement of A_s with the low-redshift measurement of ${\sigma }_{8}$. As has been pointed out many times, the Planck15 measurement of A_s implies a higher ${\sigma }_{8}$ in ΛCDM than obtained from local measurements, which leads to an apparent detection of $\sum {m}_{\nu }$ in νΛCDM. Meanwhile, CMB temperature fluctuations are sensitive to the combination A_se^−2τ—i.e., A_s and τ are positively correlated in TT parameter constraints—so imposing a τ prior with a lower central value results in a lower inferred value of A_s. In ΛCDM, this shifts the Planck-inferred ${\sigma }_{8}$ to lower values. Finally, when analyzing Planck TT+τ+SPTcl in νΛCDM, ${\sigma }_{8}$ is dominated by the local constraint from SPTcl, and the lower A_s implies that $\sum {m}_{\nu }$ need not be as high as in our fiducial analysis.

We further test the impact of using only the low-redshift half of our cluster sample. The SPTcl(0.25 < z < 0.6)+Planck data set yields

$$\begin{eqnarray}&&\sum {m}_{\nu }={0.29}_{-0.29}^{+0.09}\,\mathrm{eV}.\end{eqnarray} \tag{ 45 }$$

The probability distribution in $\sum {m}_{\nu }$ runs against the hard prior $\sum {m}_{\nu }\gt 0$, which shifts the mean recovered value away from the mode; the 68% credible interval starts at $\sum {m}_{\nu }=0$. In conclusion, all preference for a nonzero sum of neutrino masses vanishes when only considering the low-redshift half of our cluster sample. Figure 10 shows the constraints on $\sum {m}_{\nu }$ as obtained in our fiducial analysis, the analysis with the τ prior, and the analysis where we only use the low-redshift cluster data.

Zoom In Zoom Out Reset image size

The sum of neutrino masses is degenerate with the amplitude of the SZ scaling relation ${A}_{\mathrm{SZ}}$ with a correlation coefficient ${\rho }_{{A}_{\mathrm{SZ}}-\sum {m}_{\nu }}=0.83$; see Figure 11. Therefore, an improved (WL) mass calibration will improve the constraints on $\sum {m}_{\nu }$. Also note that the effect of massive neutrinos on the HMF depends (weakly) on mass and redshift (Ichiki & Takada 2012). Therefore, an improved mass calibration covering the entire cluster sample will in principle allow for measurements of the sum of neutrino masses from clusters alone.

Zoom In Zoom Out Reset image size

We consider an extension to modeling dark energy as a cosmological constant by allowing for an equation-of-state parameter w that is different from w = −1. This modification impacts the expansion history of the universe E(z) and the growth of structure; both affect the cluster abundance. Therefore, as noted in, e.g., Haiman et al. (2001), measuring the abundance over a range of redshifts allows for a measurement of w. Using our SPTcl data set, we obtain

$$\begin{eqnarray}&&w=-1.55\,\pm \,0.41\ \end{eqnarray} \tag{ 46 }$$

$$\begin{eqnarray}&&{{\rm{\Omega }}}_{{\rm{m}}}=0.299\,\pm \,0.049\end{eqnarray} \tag{ 47 }$$

$$\begin{eqnarray}&&{\sigma }_{8}=0.766\,\pm \,0.036,\end{eqnarray} \tag{ 48 }$$

as shown by the blue contours in Figure 12. Constraints on scaling relation parameters can be found in Table 3. The consistency of our recovered constraint on w with a cosmological constant w = −1 has a p-value 0.076 (1.8σ). Note that the SPTcl contours in the ${{\rm{\Omega }}}_{{\rm{m}}}$–${\sigma }_{8}$–w space close.

Zoom In Zoom Out Reset image size

Our constraint on w is in general agreement with the result obtained from the SPT-SZ+${Y}_{{\rm{X}}}$+X-ray priors data combination w = −1.28 ± 0.31 as presented in dH16. In that earlier analysis, informative (Gaussian) priors were applied on the scaling relation parameters ${A}_{\mathrm{SZ}}$, ${B}_{\mathrm{SZ}}$, ${C}_{\mathrm{SZ}}$, ${\sigma }_{\mathrm{ln}\zeta }$, ${A}_{{Y}_{{\rm{X}}}}$, ${B}_{{Y}_{{\rm{X}}}}$, ${C}_{{Y}_{{\rm{X}}}}$, ${\sigma }_{\mathrm{ln}{Y}_{{\rm{X}}}}$, whereas we marginalize over flat priors and use our internal WL mass calibration. However, even when analyzing the same data combination used in dH16 (without WL data) and applying the same priors, our analysis pipeline gives a more negative value of w = $-{1.53}_{-0.25}^{+0.36}$. As described in Section 3.4, we have extensively tested our analysis pipeline, including tests against mock catalogs with input values of w $\ne$ −1. The analysis pipeline used in dH16 was not subjected to that test. Using our internal WL mass calibration shifts the constraints on w toward even more negative values. Finally, the cluster photometric redshifts were updated since the dH16 analysis (see Section 2.1), with the net impact being a shift in w toward less negative values of similar magnitude to the shift due to our WL mass calibration. In the end, some of these shifts in w partially cancel out, and the final constraint we present here is 0.7σ low in comparison to that in dH16.

We proceed and analyze the joint SPTcl+Planck data set. The cluster data break some of the Planck parameter degeneracies shown in Figure 12, and we measure

$$\begin{eqnarray}&&w=-1.12\,\pm \,0.21\ \end{eqnarray} \tag{ 49 }$$

$$\begin{eqnarray}&&{{\rm{\Omega }}}_{{\rm{m}}}=0.347\,\pm \,0.039\end{eqnarray} \tag{ 50 }$$

$$\begin{eqnarray}&&{\sigma }_{8}=0.761\,\pm \,0.027\end{eqnarray} \tag{ 51 }$$

$$\begin{eqnarray}&&\sum {m}_{\nu }=0.50\,\pm \,0.24\,\mathrm{eV}.\end{eqnarray} \tag{ 52 }$$

Interestingly, while the individual constraints on w are both centered on w ≈ −1.5, the joint analysis provides a constraint that is offset closer toward w = −1. This is due to the different orientations of the w − ${\sigma }_{8}$ degeneracies in Figure 12 that overlap close to w = −1. Compared to the results obtained in νΛCDM, the constraints on ${\sigma }_{8}$ and ${\sigma }_{8}{\left({{\rm{\Omega }}}_{{\rm{m}}}/0.3\right)}^{0.2}$ do not degrade. However, the constraining power on the remaining cosmological parameters weakens (see Table 3).

Figure 12 further shows the constraints obtained from BAOs and SNe Ia. Neither of the two are affected by ${\sigma }_{8}$ and $\sum {m}_{\nu }$.⁵⁸ However, they both exhibit narrow parameter degeneracies that cut through the region of parameter space that is allowed by Planck. Therefore, the joint analyses of Planck+BAO and Planck+SN Ia allow for constraints on νwCDM that are tighter than the ones from Planck+SPTcl (see Figure 13).

Zoom In Zoom Out Reset image size

Finally, we analyze the joint Planck+BAO+SN Ia+SPTcl data set (see constraints in Table 3). In comparison to Planck+BAO+SN Ia, the addition of the SPTcl data set leads to a shift Δ${\sigma }_{8}$ = −0.031. The constraints on ${{\rm{\Omega }}}_{{\rm{m}}}$, h, and w are negligibly affected. However, note that the 95% upper limit on $\sum {m}_{\nu }$ from Planck+BAO+SN Ia increases by 63% when adding SPTcl. A similar effect was seen in the DES 3 × 2 pt analysis (Abbott et al. 2018), where the upper limit on $\sum {m}_{\nu }$ from Planck+BAO+SN Ia increased by a similar amount when adding the DES data. Both effects are due to the lower clustering amplitude measured by SPTcl and DES relative to the prediction by Planck+BAO+SN Ia.

4.3.1. νwCDM: Robustness of Our Results to Data Cuts

In Appendix C (Figure 17), we show the parameter constraints that we recover when cutting our cluster sample in half at redshift 0.6, or when choosing a higher SZ selection threshold ξ > 6.5. There are no significant departures from our fiducial results for any data subset. However, both the low-redshift half of the data and the subsample above ξ > 6.5 yield constraints on w that are closer to the cosmological constant w = −1:

$$\begin{eqnarray}&&w(0.25\lt z\lt 0.6)=-{1.01}_{-0.25}^{+0.41}\end{eqnarray} \tag{ 53 }$$

$$\begin{eqnarray}&&w(\xi \gt 6.5)=-{1.21}_{-0.29}^{+0.42}.\end{eqnarray} \tag{ 54 }$$

Conversely, the high-redshift half of the data gives

$$\begin{eqnarray}&&w(z\gt 0.6)=-1.58\,\pm \,0.46.\end{eqnarray} \tag{ 55 }$$

We note that the constraints on w from the full sample are quite similar to this constraint from the high-redshift half of the data.

Figure 17 further shows a strong degeneracy between w and the redshift evolution parameters of the scaling relations ${C}_{\mathrm{SZ}}$ and ${C}_{{Y}_{{\rm{X}}}}$. To tighten the dark energy constraints in future analyses, it will therefore be important to improve the mass calibration over the entire redshift range of the cluster sample.

We consider another extension to ΛCDM where we do not alter the background expansion but change the growth of structure. Clusters have been used to constrain modified structure growth by, for example, fitting for the growth index γ, which is defined by the relation $d\mathrm{ln}\delta /d\mathrm{ln}a\equiv {{\rm{\Omega }}}_{{\rm{m}}}^{\gamma }(a)$ (e.g., Peebles 1980). A value of γ ≈ 0.55 corresponds approximately to the growth rate in ΛCDM, and clusters allow for constraints on γ at the ∼40% level (e.g., Rapetti et al. 2013; Bocquet et al. 2015; Mantz et al. 2015).

Instead of modeling linear deviations from GR via the growth index, we pursue a different route and constrain the growth of structure by directly measuring the linear amplitude of the density fluctuations, ${\sigma }_{8}$, as a function of redshift. We can then compare the measured ${\sigma }_{8}(z)$ with predictions from νΛCDM, νwCDM, and more exotic models. This approach is nonparametric in that it does not assume a specific description for modified growth of structure, but rather assumes a νΛCDM model (with its parameters allowed to vary) within each redshift bin.

We start from the νΛCDM model and modify the amplitude of the linear matter power spectrum P(k, z) within different redshift bins. We introduce an additional model parameter σ₈(z_i) in each bin and normalize P(k, z) within each redshift bin i to match σ₈(z_i). The HMF is then computed from the modified P(k, z) in the usual way. We define four redshift bins such that all bins contain approximately equal numbers of SPT clusters. We choose bin limits (z = 0.25, 0.45, 0.6, 0.75, 1.7). We include Planck primary CMB data in the fit. By construction, in our model the Planck data only constrain the background cosmology (expansion history E(z)), but they do not contribute to constraining ${\sigma }_{8}$(z). For simplicity, we do not use any X-ray data here and thus use the joint SPT-SZ+WL+Planck data set; we apply the simulation prior on the SZ scatter.

We explore two scenarios: (i) We assume our fiducial SZ scaling relation model across the entire redshift range. This means that the mass calibration will be correlated across the four redshift bins. (ii) We additionally introduce independent normalizations of the SZ scaling relation A_SZ,i in each redshift bin i. This way, the amplitude of the SZ scaling relation is independently determined in each redshift bin (up to the shared WL systematics that are, however, subdominant here given the low number of clusters with WL constraints). We call the first scenario "coupled" and the second "decoupled," in reference to the treatment of the normalization of the SZ observable–mass relation.

Constraints on ${A}_{\mathrm{SZ}}$ and ${\sigma }_{8}$(z) for the coupled and decoupled analyses are shown in Table 4. In the coupled analysis, the four measurements of ${\sigma }_{8}$(z) are quite correlated with correlation coefficients ρ(${\sigma }_{8}$(z)) = 0.55–0.60 because they are limited by the uncertainty in the observable–mass scaling relation that is shared across the entire redshift range. In the decoupled analysis, however, the ${\sigma }_{8}$(z) parameters are much less correlated (ρ(${\sigma }_{8}$(z)) = 0.06–0.12), as mass calibration in each bin is done almost independently, and each ${\sigma }_{8}$(z_i) is mostly degenerate with the corresponding normalization parameter ${A}_{\mathrm{SZ}}$(z_i). As expected, the decoupled analysis leads to weaker constraints.

Table 4.  Constraints on ${\sigma }_{8}$ and the SZ Scaling Relation Normalization ${A}_{\mathrm{SZ}}$ as a Function of Redshift, Measured in Four Redshift Bins

Parameter	Coupled Analysis	Decoupled Analysis	Planck(TT+lowTEB)
${A}_{\mathrm{SZ}}$	${5.40}_{-1.15}^{+0.80}$	⋯	⋯
${A}_{\mathrm{SZ}}$(0.25 < z < 0.45)	⋯	${4.90}_{-1.09}^{+0.88}$	⋯
${A}_{\mathrm{SZ}}$(0.45 < z < 0.6)	⋯	${10.29}_{-4.23}^{+2.56}$	⋯
${A}_{\mathrm{SZ}}$(0.6 < z < 0.75)	⋯	${7.29}_{-4.47}^{+1.66}$	⋯
${A}_{\mathrm{SZ}}$(0.75 < z < 1.7)	⋯	${10.63}_{-3.08}^{+2.19}$	⋯
${\sigma }_{8}$(z = 0.35)	0.592 ± 0.031	0.609 ± 0.028	0.656 ± 0.029
${\sigma }_{8}$(z = 0.525)	0.543 ± 0.029	0.484 ± 0.034	0.597 ± 0.028
${\sigma }_{8}$(z = 0.675)	0.519 ± 0.026	0.505 ± 0.046	0.555 ± 0.026
${\sigma }_{8}$(z = 1.225)	0.415 ± 0.023	0.371 ± 0.020	0.432 ± 0.021
ρ(${\sigma }_{8}$(z))	0.55–0.60	0.06–0.12	⋯

Note. In the coupled analysis we assume a single ${A}_{\mathrm{SZ}}$ parameter. In the decoupled analysis we fit for ${A}_{\mathrm{SZ}}$ separately in each redshift bin; this decorrelates the measurements of ${\sigma }_{8}$(z) as evidenced by ρ(${\sigma }_{8}$(z)).

Download table as:  ASCIITypeset image

In Figure 14, we show measurements of ${\sigma }_{8}$ as a function of redshift. The red band shows the prediction for ${\sigma }_{8}$(z) assuming νΛCDM and Planck cosmological parameters. Blue and orange data points show measurements of ${\sigma }_{8}$(z) using our clusters (with Planck priors on the background cosmology). The cluster measurements are all slightly lower than the predictions using Planck data, which simply reflects the difference in ${\sigma }_{8}$ discussed above for the νΛCDM model (see Figure 9). We emphasize that this offset is roughly constant throughout the entire range in redshift. In particular, the two bins above z > 0.6 that are leading to some shifts in cosmology and scaling relations as described in earlier sections do not seem to provide constraints that are qualitatively different from those obtained from the low-redshift bins.

Zoom In Zoom Out Reset image size

Our measurements of ${\sigma }_{8}$(z) are limited by the determination of ${A}_{\mathrm{SZ}}$, especially in the "decoupled" analysis. The three low-redshift bins will benefit from including cluster WL data from the Dark Energy Survey (M. Paulus et al. 2019, in preparation). The highest-redshift bin can only be improved with deep, high-resolution WL data, e.g., from HST, or with lensing information from the CMB (e.g., Baxter et al. 2018). On the other hand, our cluster sample, together with this technique, allows us to place competitive constraints on the growth of structure over a wide range in redshifts.

For the νΛCDM and νwCDM analyses discussed above, Table 3 also presents constraints on the SZ scaling relation parameters. These, together with Equations (1) and (8), allow one to compute mass estimates $P({M}_{500c}| \,\xi ,z,{\boldsymbol{p}})$ for each cluster in our sample. Moreover, the scaling relation parameter constraints provide another point of comparison with past analyses.

The results in Table 3 exhibit a range of parameters across the six different analyses, but in no case are the parameter differences statistically significant. This indicates that the best estimates of the cluster masses are consistent among the different combinations of data and within the different cosmological models. As an example, the addition of the Planck data set as an external prior leads to preferred values of the amplitude parameter ${A}_{\mathrm{SZ}}$ that are lower, corresponding to ∼8% and ∼4% higher masses at the pivot in νΛCDM and νwCDM, respectively. These mass shifts are smaller than those presented by Bocquet et al. (2015), where the impact of external priors was first discussed. Interestingly, the redshift slope ${C}_{\mathrm{SZ}}$ prefers higher values in the νwCDM model, which corresponds to high-redshift masses that are smaller relative to clusters with the same ξ at low redshift. In the νwCDM model these same values of ${C}_{\mathrm{SZ}}$ ∼ 1 are preferred with or without Planck priors but shift back to a lower value when BAO+SN Ia constraints are added. In comparison, in the results for ΛCDM presented in Table 3 of dH16, the amplitude parameter for the SPTcl+Planck+BAO analysis was ${A}_{\mathrm{SZ}}$ = 3.53 ± 0.27, which is significantly lower than the values presented here. Note, however, that massive neutrinos were not marginalized over in the baseline analysis in dH16.

Given the consistency in the implied masses across all six analyses presented here, we adopt the νΛCDM results for the baseline SPTcl data set in calculating mean masses and mass uncertainties (Table 3, Column (3)). The mass uncertainties include the ξ measurement and intrinsic scatter uncertainties (together these correspond to ∼20% uncertainty for a cluster near our selection threshold), as well as marginalization over the posterior parameter distributions for ${A}_{\mathrm{SZ}}$, ${B}_{\mathrm{SZ}}$, ${C}_{\mathrm{SZ}}$, and ${\sigma }_{\mathrm{ln}\zeta }$ and over the cosmological parameters (this corresponds to an additional ∼15% uncertainty due to the remaining uncertainties in the mass calibration of the SPT-SZ sample). These masses are calculated by sampling the distribution

$$\begin{eqnarray}&&P(M| \xi ,z,{\boldsymbol{p}})=\iint {dMd}\zeta P(\xi | \zeta )P(\zeta | M,z,{\boldsymbol{p}})P(M| z,{\boldsymbol{p}})\end{eqnarray} \tag{ 56 }$$

at each step in the likelihood analysis.

Table 5 contains a list of all cluster candidates in our sample with the associated sky location, SPT detection significance ξ, redshift, and halo mass M_500c. In addition, we present the mass M_200c for each system, assuming the concentration–mass relation from Duffy et al. (2008).

Table 5.  Galaxy Cluster Candidates with ξ > 4.5 in the SPT-SZ Survey

SPT ID	R.A.	Decl.	ξ	θc	YSZ	Redshift	M500c	M200c	${M}_{500{\rm{c}}}^{\mathrm{no}\,\mathrm{syst}.}$
	(J2000)	(J2000)		(arcmin)	(10−6 arcmin−2)		(1014 M⊙/h70)	(1014 M⊙/h70)	(1014 M⊙/h70)
SPT-CL J0000-4356	0.0663	−43.9494	5.92	0.25	77 ± 21	0.73 ± 0.05*	${3.55}_{-0.83}^{+0.65}$	${5.56}_{-1.34}^{+1.04}$	${4.01}_{-0.74}^{+0.67}$
SPT-CL J0000-5748 ${}^{\mathrm{WL},{\rm{X}}}$	0.2499	−57.8064	8.49	0.50	83 ± 13	0.702	${4.33}_{-0.86}^{+0.65}$	${6.79}_{-1.37}^{+1.04}$	${4.88}_{-0.71}^{+0.59}$
SPT-CL J0001-4024	0.3610	−40.4108	5.42	0.75	60 ± 13	0.83 ± 0.03*	${3.19}_{-0.83}^{+0.65}$	${5.03}_{-1.33}^{+1.04}$	${3.62}_{-0.73}^{+0.64}$
SPT-CL J0001-4842	0.2768	−48.7132	5.69	1.25	78 ± 13	0.33 ± 0.02*	${3.89}_{-0.90}^{+0.75}$	${5.92}_{-1.40}^{+1.16}$	${4.48}_{-0.83}^{+0.75}$
SPT-CL J0001-5440	0.4059	−54.6697	5.69	1.00	60 ± 13	0.82 ± 0.08*	${3.47}_{-0.88}^{+0.64}$	${5.46}_{-1.42}^{+1.02}$	${3.88}_{-0.75}^{+0.67}$
SPT-CL J0001-6258	0.4029	−62.9808	4.69	1.50	52 ± 16	0.21 ± 0.02	${3.32}_{-0.88}^{+0.66}$	${4.97}_{-1.34}^{+1.00}$	${3.89}_{-0.87}^{+0.65}$
SPT-CL J0002-5224	0.6433	−52.4092	4.67	1.00	48 ± 13	>0.71	⋯	⋯	⋯
SPT-CL J0002-5557	0.5138	−55.9621	5.20	0.25	55 ± 18	1.15 ± 0.10	${2.84}_{-0.81}^{+0.58}$	${4.55}_{-1.32}^{+0.95}$	${3.20}_{-0.68}^{+0.56}$
SPT-CL J0003-4155	0.7842	−41.9307	4.75	0.25	59 ± 20	>0.76	⋯	⋯	⋯
SPT-CL J0007-4706	1.7514	−47.1159	4.55	0.75	58 ± 13	0.46 ± 0.04*	${3.14}_{-0.85}^{+0.62}$	${4.80}_{-1.33}^{+0.96}$	${3.65}_{-0.83}^{+0.58}$

Note. The subsample with ξ ≥ 5 and z ≥ 0.25 is used in the cosmological analysis. Clusters with follow-up WL and/or X-ray data that we use in this work are marked with "WL" and/or "X." The positions ξ, core radii θ_c, and Y_SZ are the same as presented in Bleem et al. (2015), while the redshifts marked with an asterisk have been updated. Spectroscopic redshifts are quoted without uncertainties. The mean mass estimates and mass uncertainties take the intrinsic and measurement scatter into account. We quote redshift lower limits for unconfirmed SZ detections. The mass estimates M_500c and M_200c are derived from the SPTcl data set in the νΛCDM model (Table 3, Column (3)) and are fully marginalized over cosmology and scaling relation parameter uncertainties. The estimates M_500c^no syst. are computed assuming a fixed cosmology and using the best-fit scaling relation parameters obtained from fitting the SPT-SZ number counts against that fixed cosmology (this approach was also adopted in Bleem et al. 2015). The catalog can also be found at https://pole.uchicago.edu/public/data/sptsz-clusters.

Only a portion of this table is shown here to demonstrate its form and content. A machine-readable version of the full table is available.

Download table as:  DataTypeset image