Cosmological Constraints from DES Y1 Cluster Abundances and SPT Multi-wavelength data

The left panels of Figure 4 compare the abundances of the DES Y1 RM clusters (boxes) with the corresponding mean model predictions (markers). The right panels show the residuals between the data and the model expectations for the two scatter models and data combinations considered.

Starting with our baseline data set DES-NC+SPT-OMR, the SPT multi-wavelength data carry the information to constrain the observable–mass relation parameters, while the DES Y1 RM abundance data, thanks to the SPT-OMR calibrated richness–mass relation, constrain the cosmological parameters. Specifically, the richness–mass relation parameters are constrained through the calibration of the $\xi$–mass scaling relation, which in turn is primarily informed by the weak lensing data. The X-ray data mainly affect the constraints on the intrinsic scatter parameters [see also 7]. We explicitly verified that when dropping the X-ray data, we obtain perfectly consistent results for all parameters but for the scatters $D_{SZ}$ and $D_{\lambda}$ whose mean values increase and decrease by $\sim 0.1$ ($\sim 1\sigma$), respectively.

The further inclusion of the SPT-NC data bring additional cosmological information which slightly improves the $\sigma_{8}$ and $\Omega_{{\rm m}}$ constraints — by $30\%$ and $20\%$, respectively — while shifting their confidence contours along the $S_{8}$ degeneracy direction (black dashed and green contours in figure 3). The shift of the $\sigma_{8}$ posterior can be understood by looking at figure 5 which compares the SPT-SZ number count data with predictions from the DES-NC+SPT-OMR and DES-NC+SPT-[OMR,NC] analyses. The larger $\sigma_{8}$ value preferred by the DES-NC+SPT-OMR data tend to over-predict the number of SPT-SZ clusters above $z>0.65$. Consequently, when included, the SPT number count data shift $\sigma_{8}$ towards lower values to recover the correct number of SPT-SZ clusters (see also orange contours in figure 7). Concurrently, to counterbalance the lower $\sigma_{8}$ mean value and thus keep roughly unvaried the predictions for the DES Y1 RM cluster counts, $\Omega_{m}$, $B_{\lambda}$ and $C_{\lambda}$ move toward larger values following the corresponding degeneracy directions with $\sigma_{8}$. We will further comment on the origin of this shift in section IV.1.4. Finally, the SPT abundance data improve the constraints on $B_{\rm SZ}$ and $C_{\rm SZ}$ thanks to the sensitivity of the SPT-NC likelihood to the SZ–mass scaling relation.

Moving to the modeling of $P(\lambda^{\rm ob}|\lambda)$, we find consistent results between the two models adopted for the observational noise on $\lambda$ (BKG with orange and black contours and PRJ with blue and green contours; see section III.1), albeit the PRJ model prefers a slightly lower $\Omega_{{\rm m}}$ value, driven by a shallower ($B_{\lambda}=0.86\pm 0.04$) and redshift independent ($C_{\lambda}=-0.02\pm 0.34$) richness–mass scaling relation compared to the BKG results. This result can be understood as follows: the PRJ model, which accounts also for projection and masking effects, tends to bias high the richness estimates and introduces a larger scatter between $\lambda^{\rm ob}$ and $\lambda$ compared to the BKG model. As a consequence, for a given set of cosmological and scaling relation parameters, the slope of the $\lambda^{\rm ob}$–mass relation increases, as well as the predicted cluster counts for DES Y1 RM. Given the strong degeneracy between $A_{\lambda}-A_{\rm SZ}$ and $D_{\lambda}-D_{\rm SZ}$, and the tight constraints on SZ parameters provided by the SPT-OMR data, $B_{\lambda}$ is the only parameter which can compensate for such effects by moving its posterior to lower values. Similarly, the preference for a non-evolving $\lambda$–mass scaling relation is explained by the redshift dependent bias and scatter intrinsic to the PRJ model, which is a consequence of the worsening of the photo-z accuracy with increasing redshift. These findings are consistent with those obtained in DES20, where it is shown the robustness of the cosmological posteriors to different model assumptions for $P(\lambda^{\rm ob}|\lambda)$.

As for the correlation coefficients between scatters in all the four cases analyzed the posteriors are prior dominated. We note, however, that while the posteriors of the correlation coefficients between SZ and WL and WL and $\lambda$ peak around zero, the $\rho({\rm SZ},\lambda)$ posterior always has its maximum at $\sim-0.2$, suggesting an anti-correlation between the two observables (see figure 12 in appendix A).

The four analyses perform similarly well in fitting the DES Y1 abundance data. The model predictions are all consistent within $2\sigma$ with the data but for the highest richness/redshift bin, where all the models over-predict the number counts by $\sim 35\%$ (see right panels of figure 4). Notably, while the SPT-OMR data is only available for clusters above $\lambda^{\rm ob}\gtrsim 40$, the scaling relation extrapolated at low richness provides a good fit to the DES Y1 abundance data. Our composite likelihood model and parameter degeneracies do not allow us to apply a $\chi^{2}$ statistic to assess the goodness of the fit. The same tensions between predictions and DES Y1 RM abundance data was observed in DES20, where the authors verified that dropping the highest-$\lambda/z$ bin from the data does not affect their results, but improve the goodness of the fit. Here we use the posterior predictive distribution to asses the likelihood of observing the highest-$\lambda/z$ data point given our models [see e.g. 46, section 6.3]. The method consists of drawing simulated values from the posterior predictive distribution of replicated data and comparing these mock samples to the observed data. The posterior predictive distribution is defined as:

where $y$ is the observed data vector, $y^{\rm rep}$ the replicated one, and $\theta$ the model parameters. In practice, we generate our replicated data for the highest-$\lambda/z$ by sampling the posterior distribution, $P(\theta|y)$, and drawing for each sampled $\theta$ a value from the multivariate normal distribution defined by equation 10 and covariance matrix ${\bm{C}}$. We draw $500$ samples for each of the four analyses, and fit the distributions with a Gaussian to easily quantify the likelihood of the observed data point. As can been seen in figure 6 for the two models and data combinations considered here the observed data lies within the $3\sigma$ region (dashed and dotted vertical lines), thus we conclude that the highest-$\lambda/z$ data point is not a strong outlier of the predicted distribution and our model suffices to describe it.

Similarly for the SPT-SZ abundance data, the models retrieved from the posteriors of the DES-NC+SPT-[OMR,NC] and DES-NC+SPT-[OMR,NC]+PRJ analyses provide a good fit to the SPT number counts but for the highest $\xi$ bin, where the model predictions lie at the edge of the $\sim 2\sigma$ region (see lower panel of figure 5).

As for the SPT-OMR data we inspect the goodness of the fit of the derived $P(\lambda^{\rm ob}|\xi)$ distributions against the cross-matched sample. Specifically, we verified that all the data points lie within the $3\sigma$ region of the posterior predictive distributions independently from the data combination and model assumed for the observational scatter on $\lambda^{\rm ob}$ (see figure 2).

To determine whether our data sets prefer one of the two models adopted for $P(\lambda^{\rm ob}|\lambda)$ — BKG and PRJ — we rely on the Deviance Information Criterion [hereafter DIC; 47]. Specifically, for a given model $M$ the DIC is computed from the mean $\chi^{2}$ over the posterior volume and the maximum posterior $\chi^{2}$ as:

The model with the lower DIC value either fits better the data — lower $\langle\chi^{2}\rangle$ — or has a lower level of complexity — lower $\left(\langle\chi^{2}\rangle-\chi^{2}_{\rm MaxP}\right)$. For the data combination DES-NC+SPT-OMR we obtain $\Delta{\rm DIC}={\rm DIC}(\rm PRJ)-{\rm DIC}(\rm BKG)=3.5$, while for the full data set $\Delta{\rm DIC}=-3.8$. Adopting the Jeffreys’ scale to interpret the $\Delta{\rm DIC}$ values, the DES-NC+SPT-OMR data combination has a ”positive” ($|\Delta{\rm DIC}|\in[2,5]$) – even though not ”strong” ($|\Delta{\rm DIC}|\in[-5,-10]$) or ”definitive” ($|\Delta{\rm DIC}|>10$) – preference for the BKG model, while the full data combination has a ”positive” preference for the PRJ model. Additional follow-up data extending to lower richness — as the one soon available from the combination of DES Y3 and Y6 data with the full SPT surveys or eROSITA — will help to identify the model which better describes the data.

Figure 7 compares the $\sigma_{8}$-$\Omega_{m}$ posteriors derived in this work for a $\Lambda$CDM +$\sum m_{\nu}$ cosmology including (lower panels) or excluding (upper panels) the PRJ calibration, to other results from the literature. To assess the consistency of two data sets $A$ and $B$ in the $\sigma_{8}$-$\Omega_{m}$ plane we test the hypothesis $\bm{p}_{A}-\bm{p}_{B}=0$ [see method ‘3’ in 49], where $\bm{p}_{A}$ and $\bm{p}_{B}$ are the $\sigma_{8}$-$\Omega_{m}$ posterior distributions as constrained by data sets $A$ and $B$, respectively.

Starting with the simpler scatter model (BKG, upper panels), our baseline data combination (DES-NC+SPT-OMR) is consistent within $2\sigma$ with all the probes considered here. The largest tension ($1.7\sigma$) is found with the results from DES20 (DES-[NC,$M_{\rm WL}$] in figure 7) which combine DES Y1 RM abundances and mass estimates from the stacked weak lensing signal around DES Y1 RM clusters [50]. The tension with DES-[NC,$M_{\rm WL}$] results is not surprising and reflects the different richness–mass scaling relation preferred by the DES Y1 weak lensing calibration (see also section IV.1.4). The consistency of our posteriors with the DES Y1 combined analysis of galaxy clustering and weak lensing [DES 3x2pt 48], Planck CMB data [40], and other cluster abundance studies, seems to confirm the conclusions of DES20: the tension between DES-[NC,$M_{\rm WL}$] and other probes is most likely due to flawed interpretation of the stacked weak lensing signal of redMaPPer clusters in terms of mean cluster mass.

The similar constraining power provided by our data set and SPT-SZ 2500, which combine SPT-OMR data and SPT-SZ cluster counts above $z>0.25$, indicates that the two analyses are limited by the uncertainty in the mass calibration, i.e. the data set they have in common. The lower $\sigma_{8}$ value preferred by the SPT SZ-2500 analysis [7]⁹⁹9Note that at odds with the B19 analysis, here we show results for the SPT-SZ 2500 analysis obtained assuming 3 degenerate massive neutrino species and adopting the massive neutrino prescription for the halo mass function presented in [45], consistently with our analysis. The different massive neutrino scheme and the inclusion of this prescription lowers the $\sigma_{8}$ posterior by $0.024$ (corresponding to $\sim 0.5\sigma$) compared to original results of B19. can be again understood by looking at figure 5: the cosmology preferred by the DES-NC+SPT-OMR data combination over-predict the SPT-NC by a factor of $\sim 2$, and the same trend holds for the SPT abundance data below $z=0.65$ (not shown in the figure). As a consequence, when substituting the DES-NC data with the SPT-SZ cluster number counts, the $\sigma_{8}$ posterior shifts toward lower values to accommodate the model predictions to the new abundance data.

The inclusion of SPT-NC data (DES-NC+SPT-[OMR,NC]) worsens the consistency with the other low-redshift probes considered here by shifting the $\Omega_{{\rm m}}$/$\sigma_{8}$ posteriors towards higher/lower values. In particular, the agreement is degraded with the DES 3x2pt and WtG results, with which the tension in the $\sigma_{8}$-$\Omega_{m}$ plane raises to $1.8\sigma$ and $1.9\sigma$, respectively. Notably, the full data combination is at $1.3\sigma$ tension also with results from SPT-SZ 2500 with which it shares part of the abundance data (SPT-SZ counts above $z=0.65$) and the follow-up data. The fact that the DES-NC+SPT-[OMR,NC] posteriors do not lie in the intersection of the DES-NC+SPT-OMR and SPT-SZ 2500 contours suggests the presence of some — yet not statistically significant — tension between the DES-NC, SPT-OMR and SPT-NC data, possibly driven by an imperfect modeling of the scaling relations¹⁰¹⁰10To exclude the possibility that the tension is driven by SPT-SZ abundance data at low redshift we re-analyze the SPT-SZ 2500 catalog excluding the cluster counts data below $z=0.65$ — i.e. analysing the data combination SPT-[OMR,NC] — finding posteriors fully consistent with SPT-SZ 2500 results [see also figure 16 in 7]..

On the other hand, by turning the $\sigma_{8}$-$\Omega_{m}$ degeneracy direction, the inclusion of the PRJ model (lower panel) improves the agreement of the DES-NC+SPT-OMR posteriors with the SPT-SZ 2500 results (from $1\sigma$ to $0.5\sigma$ tension), at the expense of larger, yet not significant ($1.3\sigma$), tension with CMB data (red contours). Also the tension with the DES20 results decreases ($0.7\sigma$) as a consequence of the improved consistency between the richness–mass scaling relations (see section IV.1.4). Similarly, when considering the full data combination, the PRJ model shifts the cosmological posteriors in the intersection of the DES-NC+SPT-OMR and SPT-SZ 2500 contours, solving the above mentioned tension between the three data set. We will go back to this point in the next section.

Being constrained by the SPT multi-wavelength data both the SZ and $Y_{X}$ scaling relations derived from the DES-NC+SPT-OMR analysis are perfectly consistent with those obtained in B19. The inclusion of SPT-NC data in our analysis shifts the slope of the SZ relation, $B_{\rm SZ}$, by $1.5\sigma$ towards steeper values to compensate for the larger $\Omega_{{\rm m}}$ value preferred by the full data combination. As mentioned before, the shift of the cosmological posteriors along the $S_{8}$ direction suggests the presence of some inconsistencies between the scaling relations preferred by the different data sets: DES-NC, SPT-OMR and SPT-NC. To pinpoint the source of tension we re-analyze the abundance and multi-wavelengths data independently using as cosmological priors the product of the posterior distributions obtained from the DES-NC+SPT-OMR and SPT-[OMR,NC] analyses (roughly the intersection between the black and pink contours in the upper right panel of figure 7). This test will allow us to understand why that region of the $\sigma_{8}$–$\Omega_{{\rm m}}$ plane is disfavored by the full data combination.

As can been seen in figure 8 the tension between DES-NC+SPT-OMR, SPT-NC and DES-NC+SPT-[OMR,NC] arises from the different amplitude of the richness and SZ scaling relation preferred by the abundance (blue contours) and SPT-OMR data (orange contours). The PRJ model, lowering the $A_{\lambda}$ value preferred by the abundance data (black dot-dashed contours), but leaving almost unaffected the SPT-OMR posteriors (green dot-dashed contours), largely alleviates the tension between data sets. Once we let the cosmological parameters free to vary, the tight correlation between the SZ and richness scaling relation parameters introduced by the SPT-OMR data, along with the different posteriors for the amplitudes preferred by the latter, moves the $\Omega_{{\rm m}}$ posterior of the full data combination towards larger values. The larger shift with respect to the DES-NC+SPT-OMR data combination observed for the BKG analysis can be understood in terms of the larger tension between multi-wavelength and abundance data displayed in figure  8. Despite the better agreement of the $A_{\lambda}$-$A_{\rm SZ}$ posteriors derived assuming the PRJ calibration, the DIC suggests a mild preference for this model only for the full data combination (see section IV.1.1).

Moving to the mass–richness relation, figure 9 compares the scaling relations derived in this work (hatched bands) with other results from the literature. The scaling relation from DES20 originally derived for $M_{200,m}$ has been converted to $\langle M_{500,c}|\lambda^{\rm ob},z\rangle$ imposing the condition $n(M_{500,c})\text{d}M_{500,c}=n(M_{200,m})\text{d}M_{200,m}$ to the Tinker halo mass function. The mean mass-richness relation and its uncertainty are computed from the $\lambda$–mass parameter posteriors through Bayes’ theorem as follows:

Fitting the $\langle M|\lambda^{\rm ob},z\rangle$ relation derived from DES-NC+SPT-OMR to a power law model similar to the one assumed in [50] we get¹¹¹¹11The corresponding mean richness–mass relations, $\langle\lambda^{\rm ob}|M_{500,c},z\rangle$, for both scatter models are reported for completeness in appendix B.:

The DES-NC+SPT-OMR and DES-NC+SPT-[OMR,NC] analyses provide mass–richness relations consistent with each other within one standard deviation (gray and hatched orange bands). These results are also consistent with a similar analysis performed by [27] who calibrate the $\lambda$–mass relation combining cluster counts from both SPT-SZ and SPTpol Extended Cluster Survey, richnesses obtained by matching the SZ sample with the redMaPPer DES Year 3 catalog, and assuming the fiducial cosmology $\sigma_{8}=0.8$ and $\Omega_{m}=0.3$ (magenta band). Also here, the slightly steeper M-$\lambda^{\rm ob}$ relation preferred by our data is due to the different cosmologies preferred by the DES and SPT abundance data. Indeed, when we include the SPT-NC data in our analysis, the $\langle M|\lambda^{\rm ob},z\rangle$ relation totally overlaps with the results from [27] (hatched orange and magenta bands). Similarly, [51] derived a richness–mass relation consistent with ours ($A_{\lambda}=83.3\pm 11.2$ and $B_{\lambda}=1.03\pm 0.10$) analysing the same redMaPPer-SPT matched sample and adopting as priors the results of B19. A consistent slope of the mass–richness relation is also found in the work of [11, $B_{\lambda}=0.99^{+0.06}_{-0.07}\pm 0.04$], who calibrate the richness-mass relation of a X-ray selected, optically confirmed cluster sample through galaxy dynamics. However, a direct interpretation of their results in the context of this analysis is not possible due to the different assumptions on the X-ray scaling relation and the scatter of the richness–mass relation, made in that work.

A larger than $1\sigma$ tension below $\lambda^{\rm ob}\simeq 60$ is found with the DES20 results which base their mass calibration on the stacked weak lensing analysis of [50] (cyan band in figure 9). As noted in DES20, the weak lensing mass estimates for $\lambda<30$ are responsible for the low values derived for the slope and amplitude of the richness–mass relation compared to the ones preferred by the SPT multi-wavelength data. We stress again here that the SPT-OMR data can actually constrain the richness–mass relation only at $\lambda^{\rm ob}\gtrsim 50$ and the constraints at low richness follow from the power law model assumed for the $\langle\lambda|M\rangle$ relation.

The inclusion of the PRJ calibration, increasing the fraction of low mass clusters boosted to large richnesses, lowers the mean cluster masses compared to the BKG model up to $\sim 25\%$ at $\lambda^{\rm ob}\lesssim 60$ (compare green and yellow with gray and orange bands in figure 9, respectively). Specifically, from the DES-NC+SPT-OMR+PRJ analysis we obtain:

The improved consistency between the scaling relations derived from the analyses adopting the PRJ calibration and DES20 reflects the improved agreement between the corresponding cosmological posteriors due to the lower $\Omega_{{\rm m}}$ value preferred by the former (see figure 7). The fact that the mass-richness relations derived from the two $P(\lambda^{\rm ob}|\lambda^{\rm true})$ models display a larger than $1\sigma$ tension below $\lambda^{\rm ob}\lesssim 50$, but perform equally well in fitting the data (see section IV.1.1), is due to the lack of multi-wavelength data at low richness. Additional follow-up data at $\lambda^{\rm ob}\lesssim 40$ will be fundamental to clearly reject one of the models and thus enable the full exploitation of the cosmological information carried by photometrically selected cluster catalogs.

It is worth noting that at odds with other studies which rely on stacked weak lensing measurements to calibrate the mean scaling relation [e.g. 8, 9], the SPT-OMR data, allowing a cluster by cluster analysis (see equation 8), enable to constrain also the scatter of the richness–mass relation. This is particularly relevant for the analysis of optically–selected cluster samples for which reliable simulation-based priors on the scatter are not available: if constrained only by the abundance data, the scatter parameter becomes degenerate with $\Omega_{m}$ and $\sigma_{8}$, degrading the constraining power of the sample [e.g. see discussion in 8].

To better investigate the implications of the derived scaling relation for low richness objects we compare in figure 10 our predictions for the mean cluster masses in different richness/redshift bins to the mean weak lensing mass estimates from [50] (filled boxes). We also include the weak lensing mass estimates employed in the DES Y1 cluster analysis (hatched boxes) which adopt an updated calibration of the selection bias based on the simulation analysis of [52] [see also appendix D of 9]. Both weak lensing mass estimates and mean mass predictions have been derived assuming $\Omega_{{\rm m}}=0.3$, $h=0.7$ and $\sigma_{8}=0.8$¹²¹²12The larger tension seen in figure 9 between the scaling relations derived in this work and [50] is due to the different cosmology preferred by the two analyses.. The mean mass predictions for the DES-NC+SPT-OMR analysis are in tension with both weak lensing mass estimates. In particular, in the lowest richness bins, $\lambda^{\rm ob}\in[20,30]$, the mean mass predictions are $25\%$ to $40\%$ higher than the weak lensing mass estimates, while they are consistent within 1 $\sigma$ with the lensing masses at $\lambda^{\rm ob}>30$. The inclusion of the PRJ model, lowering the mean mass predictions, largely reduce the tension at low richness with both weak lensing mass estimates, while at $\lambda^{\rm ob}>30$ the model predictions are consistent within 1 $\sigma$ with the weak lensing masses derived adopting the selection effect bias calibration of [52]. These results are consistent with those of DES20: for the DES Y1 cluster cosmology analysis to be consistent with other probes the weak lensing mass estimates of $\lambda^{\rm ob}<30$ systems need to be boosted. Or conversely, the weak lensing mass estimates of $\lambda^{\rm ob}<30$ systems are biased low compared to the mean masses predicted by DES Y1 abundance data alone assuming a cosmology consistent with other probes. As discussed in DES20 this tension might be due to an overestimate of the selection effect correction at low richness, or to another systematic not captured by the current synthetic cluster catalogs. The good agreement of the PRJ mass predictions with the weak lensing masses adopted in DES20 reflects the consistency of our cosmological posteriors with those derived in DES Y1 cluster analysis (see the lower left panel of figure 7). The same conclusions last also for the full data combination analyses (not shown in figure 10), which provide model predictions fully consistent with those obtained from the combination of DES abundances and SPT multi-wavelength data.