Detection of anti-correlation of hot and cold baryons in galaxy clusters

TL;DR

This paper reports the observational detection of an anti-correlation between hot gas and cold stellar matter in galaxy clusters, supporting simulation predictions about baryon distribution.

Contribution

It provides the first observational evidence of hot and cold baryon anti-correlation in galaxy clusters, using a novel hierarchical Bayesian analysis of multi-wavelength data.

Findings

Detection of anti-correlation consistent with simulations

Use of hierarchical Bayesian model for analysis

Analysis of 41 galaxy clusters across multiple wavelengths

Abstract

The largest clusters of galaxies in the Universe contain vast amounts of dark matter, plus baryonic matter in two principal phases, a majority hot gas component and a minority cold stellar phase comprising stars, compact objects, and low-temperature gas. Hydrodynamic simulations indicate that the highest-mass systems retain the cosmic fraction of baryons, a natural consequence of which is anti-correlation between the masses of hot gas and stars within dark matter halos of fixed total mass. We report observational detection of this anti-correlation based on 4 elements of a $9\times9$ element covariance matrix for nine cluster properties, measured from X-ray, optical, infrared and millimetre wavelength observations of 41 clusters from the Local Cluster Substructure Survey. These clusters were selected using explicit and quantitative selection rules that were then encoded in our hierarchical Bayesian model. Our detection of anti-correlation is consistent with predictions from contemporary hydrodynamic cosmological simulations that were not tuned to reproduce this signal.

Tables

Table 1. Table 1 : Non-selection elements of the property vector of the LoCuSS cluster sample. Each is an integated observed quantity or composite thereof. The full analysis includes the original selection property, L X , RASS subscript 𝐿 𝑋 RASS L_{X,\rm RASS} , described in the companion paper [ 19 ] . The subscript 500 below is the multiple of the universe’s critical density, ρ crit ​ ( z ) subscript 𝜌 crit 𝑧 \rho_{\rm crit}(z) ) employed in the the enclosed density condition used to define physical size of cluster, M 500 ≡ M ( < R 500 ) = ( 4 ​ π 3 ) ​ 500 ​ ρ crit ​ ( z ) ​ R 500 3 subscript 𝑀 500 annotated 𝑀 absent subscript 𝑅 500 4 𝜋 3 500 subscript 𝜌 crit 𝑧 superscript subscript 𝑅 500 3 M_{500}\equiv M(<R_{500})=(\frac{4\pi}{3})500\rho_{\rm crit}(z)R_{500}^{3} . Size estimates are derived from weak-lensing (WL) measurements. Ref. [ 19 ] provides details of the instruments and methods used for each property.

Element	Unit	Description
$L_{X,\mathrm{ce}}$$E^{-1}\left(z\right)$	${10}^{44}\mathrm{erg}{s}^{-1}$	core-excised, bolometric X-ray luminosity
$k_B$ $T_{X,\mathrm{ce}}$	keV	core-excised intracluster medium (ICM) thermal temperature
$M_{\mathrm{gas}}$$E\left(z\right)$	${10}^{14}{M}_{\odot }$	ICM gas mass within WL $R_{500}$
$Y_X$$E\left(z\right)$	${10}^{14}{M}_{\odot }\mathrm{keV}$	ICM (spherical) X-ray thermal energy within WL $R_{500}$
$Y_{\mathrm{SZA}}$$E\left(z\right)$	${10}^{-5}{\mathrm{Mpc}}^2$	ICM (spherical) SZ amplitude
$Y_{\mathrm{Pl}}$$E\left(z\right)$	${10}^{-5}{\mathrm{Mpc}}^2$	ICM (cylindrical) SZ amplitude
$L_{K,\mathrm{BCG}}$$E\left(z\right)$	${10}^{12}{L}_{\odot }$	Brightest Cluster Galaxy (BCG) $K$-band luminosity
$L_{K,\mathrm{tot}}$$E\left(z\right)$	${10}^{12}{L}_{\odot }$	total $K$-band luminosity within WL $R_{500}$
$\lambda$$E\left(z\right)$	none	redMaPPer richness (count of galaxies)
$M_{\mathrm{WL}}$ $E\left(z\right)$	${10}^{14}{M}_{\odot }$	weak-lensing estimate of mass, $M_{500}$

Table 2. Table 2 : Lower Triangle : Posterior median and 68 th percentile range of property pair correlation coefficients at fixed weak-lensing mass. Color encodes the magnitude and sign of the correlation coefficient, with red (blue) showing positive (negative) values. Upper Triangle : Statistical significance (p-value) of the sign of the estimated property correlation, calculated as the cumulative posterior probability of having positive (negative) correlation values if the median is negative (positive). Diagonal : Posterior median and 68 th percentile range of the intrinsic scatter of each property at fixed weak-lensing mass.

	$L_{X,\mathrm{ce}}$	$M_{\mathrm{gas}}$	$T_{X,\mathrm{ce}}$	$Y_X$	$Y_{\mathrm{Pl}}$	$Y_{\mathrm{SZA}}$	$L_{K,\mathrm{BCG}}$	$\lambda$	$L_{K,\mathrm{tot}}$
$L_{X,\mathrm{ce}}$	${0.36}_{-0.04}^{+0.05}$	$0.001$	$0.001$		$0.007$	$0.03$	$0.09$	$0.09$	$0.018$
$M_{\mathrm{gas}}$	${0.76}_{-0.13}^{+0.09}$	${0.15}_{-0.03}^{+0.03}$	$0.28$	$0.008$	$0.011$	$0.37$	$0.46$	$0.30$	$0.06$
$T_{X,\mathrm{ce}}$	${0.49}_{-0.16}^{+0.13}$	${0.13}_{-0.22}^{+0.20}$	${0.19}_{-0.03}^{+0.04}$		$0.22$	$0.06$	$0.36$	$0.16$	$0.07$
$Y_X$	${0.85}_{-0.09}^{+0.07}$	${0.59}_{-0.19}^{+0.14}$	${0.67}_{-0.14}^{+0.11}$	${0.32}_{-0.05}^{+0.05}$	$0.003$	$0.04$	$0.10$	$0.33$	$0.04$
$Y_{\mathrm{Pl}}$	${0.54}_{-0.16}^{+0.12}$	${0.56}_{-0.18}^{+0.13}$	${0.15}_{-0.20}^{+0.18}$	${0.60}_{-0.16}^{+0.12}$	${0.28}_{-0.04}^{+0.04}$	$0.47$	$0.21$	$0.14$	$0.43$
$Y_{\mathrm{SZA}}$	${0.41}_{-0.20}^{+0.16}$	${0.08}_{-0.25}^{+0.21}$	${0.36}_{-0.22}^{+0.20}$	${0.39}_{-0.21}^{+0.17}$	$-{0.02}_{-0.22}^{+0.21}$	${0.29}_{-0.07}^{+0.08}$	$0.40$	$0.32$	$0.25$
$L_{K,\mathrm{BCG}}$	${0.19}_{-0.14}^{+0.13}$	$-{0.02}_{-0.19}^{+0.18}$	${0.06}_{-0.17}^{+0.17}$	${0.21}_{-0.16}^{+0.15}$	${0.13}_{-0.16}^{+0.16}$	${0.05}_{-0.21}^{+0.20}$	${0.34}_{-0.04}^{+0.05}$	$0.17$	$0.43$
$\lambda$	$-{0.24}_{-0.17}^{+0.18}$	$-{0.11}_{-0.21}^{+0.21}$	$-{0.24}_{-0.23}^{+0.24}$	$-{0.08}_{-0.20}^{+0.19}$	${0.20}_{-0.19}^{+0.16}$	$-{0.12}_{-0.25}^{+0.25}$	$-{0.17}_{-0.18}^{+0.18}$	${0.25}_{-0.04}^{+0.05}$	$0.014$
$L_{K,\mathrm{tot}}$	$-{0.59}_{-0.22}^{+0.27}$	$-{0.56}_{-0.28}^{+0.36}$	$-{0.48}_{-0.29}^{+0.32}$	$-{0.53}_{-0.26}^{+0.30}$	${0.05}_{-0.32}^{+0.29}$	$-{0.26}_{-0.39}^{+0.37}$	${0.06}_{-0.32}^{+0.32}$	${0.77}_{-0.27}^{+0.16}$	${0.09}_{-0.03}^{+0.05}$

Table 3. Table 3 : The projection covariance models employed in the Supplementary Material to assess the effect of the projections on the inferred anti-correlation signal. The diagonal elements are the scatter and the off-diagonal elements are the correlation coefficients.

	Mass	Optical	X-ray
Mass	0
Optical	0	0.15
X-ray	0	r (input)	0.15

Table 4. (a) Input – Σ input subscript Σ input \Sigma_{\rm input}

	Mass	Optical	X-ray
Mass	0
Optical	0	0.15
X-ray	0	r (input)	0.15

Table 5. (b) Projection Model 1 – Σ projection , 1 subscript Σ projection 1 \Sigma_{\rm projection,1}

	Mass	Optical	X-ray
Mass	0.10
Optical	0.50	0.10
X-ray	0.00	0.00	0.05

Table 6. (c) Projection Model 2 – Σ projection , 2 subscript Σ projection 2 \Sigma_{\rm projection,2}

	Mass	Optical	X-ray
Mass	0.20
Optical	0.50	0.20
X-ray	0.00	0.00	0.10

Table 7. (d) Projection Model 3 – Σ projection , 3 subscript Σ projection 3 \Sigma_{\rm projection,3}

	Mass	Optical	X-ray
Mass	0.20
Optical	0.50	0.20
X-ray	0.25	0.25	0.10

Table 8. (e) Projection Model 4 – Σ projection , 4 subscript Σ projection 4 \Sigma_{\rm projection,4}

	Mass	Optical	X-ray
Mass	0.20
Optical	0.75	0.40
X-ray	0.00	0.00	0.10

Equations

$$⟨s_a|\mu ,z⟩={\pi }_a+{\alpha }_a\mu ,$$

$$C_{a,b}=\frac{N}{N-1}\sum _{i=1}^n\text{leavevmode}\delta s_{a,i}\text{leavevmode}\delta s_{b,i},$$

$$r_{a,b}=\frac{C_{a,b}}{\sqrt{C_{a,a}{C}_{b,b}}}.$$

$$r_{\ln \lambda ,\mu |M_{\mathrm{WL}}}=\frac{{\sigma }_{\mu |M_{\mathrm{WL}}}/{\sigma }_{\mu |\lambda }-r}{{\left(1-r^2+({\sigma }_{\mu |M_{\mathrm{WL}}}/{\sigma }_{\mu |\lambda }-r{)}^2\right)}^{1/2}},$$

$${\sigma }_{\mu |\lambda }\ge {\sigma }_{\mu |M_{\mathrm{WL}}}/r=0.2/r.$$

Full text

Introduction

Dark matter, whose nature remains elusive, and ordinary matter described by the Standard Model of particle physics, are the strongly clustered materials of our Universe, with the latter component referred to as baryonic matter, or more simply baryons, by observational cosmologists[1]. The question of how well these two components trace one another, across spatial scales and cosmic time, is central to our understanding of the astrophysics that drives galaxy formation and offers clues to the thermal nature of dark matter and other new physics.

Assuming weak-field (Newtonian) gravitational accretion and collisional shocks under the approximation of spherical symmetry, self-similar solutions [2] emerge in which both collisionless dark matter and collisional baryonic fluids develop similar radial profiles when expressed in terms of a characteristic physical radius, the turn-around radius [3] at which the perturbed Hubble flow is stationary. A key implication of this simple model is that dark matter and baryons exhibit no radial separation. Collapsed structures, referred to as halos, should retain the cosmic mix of these different fluids at all radii.

The most massive dark matter halos host groups and clusters of galaxies. Early X-ray measurements of the hot gas content of clusters upended the standard cold dark matter (CDM) model orthodoxy of a matter dominated universe[4] before observations of Type Ia supernovae ushered in the current $\Lambda$CDM model of a universe dominated by a smooth dark energy component[5]. The argument against a matter dominated universe relied on a fair sampling hypothesis, namely that the mean baryon fraction within clusters (ratio of baryonic mass to total mass) accurately reflects the cosmic mean baryon fraction. A natural consequence of this hypothesis is that the hot and cold fractions of baryons in clusters should be anti-correlated; at fixed total mass, clusters with more cold baryons should have less hot baryons, and vice versa. While this model does not define how baryons are partitioned into these phases, the constancy of the sum implies that a particular system with more hot gas than average must contain a lower stellar mass than average, and vice-versa.

However, this simple model ignores important non-spherical and non-gravitational effects such as hierarchical mergers driven by large-scale filaments and the redistribution of energy, momentum and mass (generically termed feedback) by supernovae and active galactic nuclei (AGN). In low mass halos that host only one bright galaxy like the Milky Way, feedback is energetic enough to vent hot gas phase baryons out of these relatively shallow gravitational potentials[6]. Even smaller halos of dwarf or satellite galaxies suffer severe baryon losses from collective supernova explosions. At the other extreme, the massive halos of rich galaxy clusters have much deeper gravitational potentials that shield them from feedback-driven baryon venting outside of their core regions[7]. Thus, clusters are likely to be closed baryon boxes, unbiased reservoirs of the cosmic baryon fraction.

Studies of mean trends in gas and stellar mass fraction[8, 9] support the expectation that massive clusters are more closed than smaller halos, based on the trend of increasing baryon fraction with halo mass. However, measurements of absolute baryon fractions are currently subject to uncertain biases of ${\cal O}(10\%)$ in estimates of total mass, and this systematic uncertainty limits the reliability of comparison with the cosmic mean baryon fraction. We take a complementary approach based on variance about mean behavior, particularly the covariance of hot gas mass and stellar mass conditioned on total mass. This approach is encouraged by recent findings of strongly negative correlation coefficients ($r\lesssim-0.5$) from a pair of complex, multi-fluid cosmological simulations in which this statistic has been measured[10, 11].

Observational studies have explored baryonic properties conditioned on estimated halo mass, particularly X-ray and thermal Sunyaev-Zel’dovich (SZ) Effect[12] signatures from the hot gas phase and optical/infrared properties of galaxies[13]. While correlations among internal hot gas properties have been measured in a few empirical studies [14, 15, 16, 17], the degree of correlation between hot gas and galactic components has not yet been investigated. The minimum requirement for such an analysis is to obtain high quality observations of both stellar and gas properties for a cluster sample with well-defined selection rules and robust estimates of total cluster mass. Currently, these requirements are only fulfilled by the Local Cluster Substructure Survey (LoCuSS), a multi-wavelength survey of the 41 X-ray brightest galaxy clusters at redshifts of $0.15<z<0.3$.

The LoCuSS sample is selected by applying a redshift-dependent X-ray luminosity ($L_{X}$) cut to clusters identified in the ROSAT All-sky Survey (RASS) catalog at high galactic latitudes. The multi-wavelength observations used in this study, obtained over the period of a decade (2005-2014) by co-authors, includes optical imaging data from the Subaru 8.2-m telescope, infrared data from the 3.8-m United Kingdom Infrared Telescope on Mauna Kea (UKIRT), X-ray observations from the Chandra and XMM-Newton satellites, and millimeter observations of the thermal SZ effect from the Planck satellite and the Sunyaev-Zel’dovich Array (SZA).

Here we report the first observational detection of anti-correlation between the hot and cold baryon contents of galaxy clusters. The key measurements are a subset of posterior estimates of 36 pairwise correlations among nine cluster properties derived from the LoCuSS observations, most measured within a radial scale defined by the weak-lensing estimate of each system’s mass[18]. Details of the galaxy cluster sample and posterior measures of the slope, variance, and $L_{X,\rm RASS}$–property covariance for nine properties are presented in a companion work[19]. Our detection of anti-correlation supports independent evidence that massive galaxy clusters retain close to the cosmic mix of baryons and dark matter, a finding that can underpin improved cluster cosmology from cross-wavelength sample analysis.

Results

Table 1 lists the observable properties employed in this analysis. We model the data using a likelihood based on log-normal property covariance about mean scaling relations behaving as power-laws in halo mass[20]. We employ default redshift scaling behaviors[21], but this assumption is unimportant due to the narrow redshift range of the sample. We assume that, on average, weak gravitational lensing measurements provide unbiased estimates of true cluster masses with 0.2 fractional scatter. To model the selection effect, we employ the threshold selection condition for $L_{X,\rm RASS}$ emission used to define the LoCuSS cluster sample. X-ray emission offers clearer identification of massive halos, being less prone than cluster properties measured at other wavelengths to confusion from additional halos projected along the line of sight. While imprecise models of sample selection can bias scaling parameter estimates, we show in the Supplementary Material that inferred correlations among property pairs are insensitive to biases in posterior slope and variance estimates (See Supplementary Figures 4 and 5).

Our analysis takes a hierarchical Bayesian approach that accounts for the effects of the sample selection, measurement error covariance induced by the use of a common sky aperture and other effects, as well as the halo space density as a function of mass for a $\Lambda$CDM cosmology. Uninformative priors are used in the regression; all quoted constraints are derived solely from the sample data. Our model simultaneously constrains the population scaling parameters associated with the multi-wavelength ensemble, the slopes, normalizations, and the mass-conditioned property covariance (see equation (2) in the Methods section). We report pairwise correlation coefficients, i.e. the covariance divided by the intrinsic scatter of each observable, as in equation (3).

An analysis of variance must be cognizant of astrophysical sources of scatter extrinsic to the host halos of the cluster sample. In particular, other halos along the line-of-sight will add correlated noise to some of a cluster’s observed properties[22]. In the Supplementary Figure 6, we show that such sources of systematic error, including projection, tend to dilute the magnitude of an intrinsically anti-correlated property pair. We argue that it is conservative, then, to consider the measured magnitude of an anti-correlation between stellar mass and hot gas mass as effectively a lower limit to the underlying halo population value.

Table 2 presents the full property covariance matrix at fixed weak-lensing mass derived from the LoCuSS sample. While we report the entire matrix, our focus is mainly on the last two rows and columns that contain optical properties. The lower triangle elements summarize the correlation coefficients of property pairs while the diagonal elements provide standard deviations of each property. Median values from the Markov chain Monte Carlo (MCMC) chains, $\sim 10^{5}$ in length, are quoted, and uncertainties give 68% confidence limits. As explained below in Methods, we impose a minimum value of $0.05$ on the intrinsic scatter in the log of K-band luminosity, $\ln$ $L_{K,\rm tot}$, at fixed halo mass when determining these statistics. The upper triangle gives the odds that each element has a sign opposite to that of its median value.

We first note the physically sensible result that the two independently measured properties reflecting a halo’s stellar mass – total $K$–band luminosity, $L_{K,\rm tot}$, and optical richness $\lambda$ – have a strong positive correlation, $r=0.77^{+0.16}_{-0.27}$. The probability of this value being negative is very small, $1.4\%$. Of these two measures, the quantity $\lambda$ appears noisier, with median intrinsic scatter of $25\%$ compared to only $9\%$ for $L_{K,\rm tot}$, but this may also reflect the different measurement errors quoted for each property. The fractional statistical uncertainties in $\lambda$ are a factor $\sim 3$ smaller than those in $L_{K,\rm tot}$. As we show in the Supplementary Figure 5, bias and/or extra noise relative to the underlying halo population statistics will tend to dilute measured (anti-)correlations, and these effects can explain why some galaxy-hot gas property pairs yield weaker evidence of anti-correlation.

In the bottom two rows of Table 2, the elements linking galaxy and hot gas properties are mainly negative. Figure 1 highlights the correlation coefficients between the galaxy measures and two key measures of hot gas: the core-excised X-ray luminosity and the derived gas mass. Boxes show inner quartiles ($25-75$-percentile) and whiskers encompass the inner $95$ percent of the marginalized posterior distributions. All pairs tend to be negative, as anticipated by the correlations between hot gas mass and stellar mass seen in hydrodynamical simulations [10, 11], shown as background bands in Figure 1. The consistency in sign of hot-cold phase covariance elements between observed proxy measures and their simulation-derived counterparts is an encouraging sign of fidelity in the sophisticated astrophysical treatments employed to model the coupled evolution of multiple baryon components at sub-cluster scales in modern cosmological simulations[23]. A consistent feature of such simulations is that the mean baryon fraction measured within the characteristic $R_{500}$ length scale used in this work approaches the cosmic value as system mass increases. This aspect, along with a reduction in the population variance in baryon component mass fractions, supports the fair sampling hypothesis and allows the most massive clusters to serve as cosmic distance rulers[24].

Due to the modest sample size, the uncertainty on any individual correlation coefficient remains large. Examining the upper triangle of Table 2, we find that the pairing of $L_{K,\rm tot}$ and $L_{X,\rm ce}$ is the strongest indicator of anti-correlation, with only a $1.8\%$ chance of being zero or positive. As noted in Fig. 3 of the Methods section, the odds of a positive $L_{K,\rm tot}$–$L_{X,\rm ce}$ correlation drop below one percent if the intrinsic scatter in $L_{K,\rm tot}$ at fixed halo mass is larger than 0.07. For the $\lambda$ measure of stellar mass, the evidence is somewhat weaker, with a 9% chance that it correlates positively with $L_{X,\rm ce}$ at fixed halo mass.

It has previously been argued[25, 26] that the $K$-band integrated light is a more accurate indicator of a cluster’s total stellar mass than the number of optically-selected galaxies, i.e. optical richness. Our results appear to reinforce this finding, as the anti-correlations for $L_{K,\rm tot}$ and X-ray properties are systematically more negative than those inferred using $\lambda$. But, as noted above, underestimation of the statistical uncertainty in $\lambda$ could also play a role in diluting $\lambda$ correlations[27] (see Supplementary Figures 5 and 6).

The findings of anti-correlation using core-excised X-ray luminosity are reinforced by the derived gas mass, $M_{\rm gas}$. Again, the infrared light provides a tighter constraint, with only a $6\%$ chance of being zero or positive, while the odds rise to $30\%$ when using $\lambda$. In the companion paper, we note that the slope of the $M_{\rm gas}$ scaling with halo mass is $\sim 2.5\sigma$ lower than values derived by previous studies based solely on X-ray observations and also slopes inferred from modern hydrodynamic simulations. A bias in slope could dilute the anti-correlation signal and explain why $L_{X,\rm ce}$ provides more significant evidence of anti-correlation (Supplementary Figure 4).

Figure 2 shows posterior constraints of correlations among different hot gas properties. Our results (shaded bands) of mainly strong positive correlations are consistent with both hydrodynamic simulation expectations [28, 29] and previous empirical measurements [16, 30, 17]. Our constraints are broader in scope, i.e. a larger number of measurements of the ICM including $Y_{X}$, $Y_{\rm SZA}$, and $Y_{\rm Pl}$, and in most cases, more precise than the few existing estimates.

Discussion

Property covariance has been forecast to significantly improve the power of joint, multi-wavelength survey analysis, especially in the case of anti-correlated properties [31]. This work helps set the stage for such analysis by providing initial estimates of stellar and hot gas covariance and refined estimates of a larger number of property correlations. While statistical errors in our correlation estimates are currently large, the coming decades will see an explosion of multi-wavelength cluster data from wide-area surveys such as Large Synoptic Survey Telescope (LSST), Euclid, e-ROSITA, and the Stage-4 ground-based cosmic microwave background experiment (CMB-S4). These upcoming samples will allow a better understanding of the physics and feedback effects that regulate the ICM along with improved cosmological constraints from joint sample analysis.

For example, the application of cluster gas fraction as a standard ruler[32] could benefit from the inclusion of stellar mass estimates. Just as lightcurve shape and color corrections are used to improve the quality of Type Ia supernovae distances, so could stellar mass measurements be used to derive a lower scatter distance proxy than provided by the gas fraction alone.

On the computational side, the sensitivity of cluster property population statistics to modeling treatments will become more apparent as computational advances in multi-phase plasma astrophysics enable refinements of processes at sub-resolution scales. Synthetic observations of halo populations produced along past lightcones under model-specific conditions, when mapped through survey-specific observational filters, offer a pathway for likelihood testing of increasingly sensitive, multi-wavelength observational surveys. A new era has arrived for studies of clusters of galaxies as a population, one in which astrophysics-dependent population statistics realized from simulations are tested against corresponding multi-wavelength, empirical data, with outcomes driving improvements in next generation models. Our results reinforce the discovery power of applying population statistical analysis to galaxy cluster samples with complete, uniform multi-wavelength observations that probe hot and cold phase baryons and total mass.

Methods

Here we briefly describe the data vector and regression analysis model employed in this study. Further details are provided in the companion paper [19] that discusses mean scaling behaviors and property variance. This paper presents off-diagonal property covariance terms, except for $L_{X,\rm RASS}$ correlations which are presented in the companion paper.

Cosmology and notation. We assume a universe with dimensionless energy densities at the current time in total matter (baryons plus dark matter) $\Omega_{\rm m}=0.3$ and vacuum energy $\Omega_{\Lambda}=0.7$, with Hubble constant $H_{0}=70\,{\rm km\,s^{-1}\,Mpc^{-1}}$. The Hubble expansion rate is normalized via $E(z)\equiv H(z)/H_{0}=\sqrt{\Omega_{\rm m}(1+z)^{3}+\Omega_{\Lambda}}$. For the halo population, we employ a mass scale convention, $M_{500}$, defined as the mass within a sphere, of radius $r_{500}$, within which the mean enclosed density is $500\rho_{\rm crit}(z)$, where $\rho_{\rm crit}(z)=3H(z)^{2}/8\pi G$, is the critical density of the universe[1]. Unless stated otherwise, the weak-lensing determined radius, $r_{500,\rm WL}$, defines the aperture within which integrated observable properties are derived.

A multi-wavelength vector of observables. We study 41 X-ray bright clusters of the LoCuSS sample derived from RASS [33, 34, 35]. The sample is selected by redshift-dependent thresholds of X-ray luminosity, $L_{X,\rm RASS}E^{-1}(z)>4.4\times 10^{44}\,\rm erg\,s^{-1}$ for clusters between $0.15<z<0.24$ and $L_{X,\rm RASS}E^{-1}(z)>7.0\times 10^{44}\,\rm erg\,s^{-1}$ for $0.24<z<0.30$. For each cluster nine additional properties, listed in Table 1, have been measured. Details are provided in the companion paper [19]. The sample is complete in most, but not all, properties, as detailed below. The integrated observables can be grouped into three distinct sets: (i) a weak-lensing mass estimate of total system mass; (ii) quantities associated with the hot intracluster gas, and; (iii) quantities associated with stellar properties. We briefly describe each set as follows.

(i) Weak-Lensing Mass: Using deep, multi-band optical images from Subaru/Suprime-Cam, a mass estimate for each cluster is derived by fitting the shear signal expected from weak gravitational lensing of a projected Navarro–Frenk–White (NFW, ref. [36]) mass density profile to the measured tangential shear pattern [18].

(ii) Hot Gas Content: Properties of the hot gas content of clusters are mostly observable at X-ray and millimeter wavelengths. We use X-ray measurements of the ICM derived in ref. [37], where the selected sample has been observed with either or both of the Chandra and XMM-Newton X-ray observatories.

To avoid contamination from the complex cool core region, measurements of bolometric luminosity, $L_{X}$, and gas temperature, $T_{X}$, are performed in an annulus of $[0.15-1]r_{500,\rm WL}$. The gas mass, $M_{\rm gas}$, is estimated from the observed X-ray emission profile within $r_{500,\rm WL}$. The SZ effect, caused by the inverse Compton scattering of cosmic microwave background (CMB) photons by hot electrons in the ICM[12], is characterized by the parameter $Y$, which is proportional to the integrated electron thermal energy. The SZ effect from CMB intensity maps, $Y_{\rm SZ}$, is measured via interferometry with SZA and independently with spectral filtering of Planck satellite data. A third estimate of the integrated electron thermal energy, $Y_{X}$, is derived from the X-ray observations as the product of gas mass, $M_{\rm gas}$, and temperature, $T_{X}$. This quantity is measured within its own iteratively-defined $r_{500}$, as discussed in ref. [19].

(iii) Stellar Content: We employ two independent measures of the stellar content of clusters, the total near-infrared (NIR) luminosity, $L_{K}$, and a count of red-sequence galaxies, $\lambda$, referred to as optical richness. The NIR luminosity measurements, obtained with the WFCAM instrument on the UKIRT telescope [26], determine the background-subtracted light within the weak-lensing estimated radius, $r_{500}$, of each cluster, as well as $L_{K}$ of the BCG. NIR data is missing for one cluster (Abell2697). The optical richness, $\lambda$, a measure of the number of red-sequence galaxies within the cluster used by the redMaPPer cluster detection algorithm[38], is determined for 33 clusters in the overlap region of the LoCuSS sample and the Sloan Digital Sky Survey (SDSS, ref. [39]).

Regression model. We assume a log-normal probability distribution of cluster properties with mean values that scale as a power-law in halo mass and $E(z)$. Because of the narrow redshift range of the LoCuSS sample, we assume standard, self-similar evolution in redshift. We employ a hierarchical Bayesian inference model that accounts for the sample selection truncation, measurement error covariance and intrinsic property covariance. An additional component of this inference model is a prior function on true halo masses derived from the halo mass function in the reference $\Lambda$CDM cosmology with $\sigma_{8}=0.8$. The performance of this method to recover input scaling relations of synthetic, truncated samples is demonstrated in the companion paper [19].

The key element of our model is the conditional joint property likelihood[20], $p(\mathbf{S}\,|\,{M}_{\rm halo},z)$, of a vector of observables, $\mathbf{S}$ (elements in Table 1), given the true mass of the halo, ${M}_{\rm halo}$, at redshift, $z$. For the LoCuSS sample clusters, we assume that the cluster weak-lensing mass, $M_{WL}$, is an unbiased measure of ${M}_{\rm halo}$ with 20% fractional scatter. Our method returns posterior estimates of the intercepts, slopes, and intrinsic variance of each property element as a function of the cluster weak-lensing mass, along with the covariance of pairs of observables. The latter is assumed to be independent of mass and redshift within the narrow ranges probed by the LoCuSS sample. Uninformative priors are used in the analysis.

Using natural logarithms of the properties, $\mathbf{s}=\ln\mathbf{S}$, and mass, $\mu=\ln{M}_{\rm halo}$, the log-mean scaling of observable $a$ at a fixed redshift is linear

[TABLE]

in which $\alpha_{a}$ and $\pi_{a}$ are the slope and normalization of the scaling relation of property $a$.

For a pair of observables, $a$ and $b$, the intrinsic property covariance matrix is

[TABLE]

where $\delta s_{a,i}\equiv s_{a,i}-\alpha_{a}\mu_{i}-\pi_{a}$ is the residual deviation from the mean scaling relation and $N$ is the total number of clusters. Finally the property correlation coefficient is

[TABLE]

This correlation coefficient is the quantity of interest that is studied in this letter. Our method constrains these correlation coefficients and the scaling parameters simultaneously, while including a covariance contribution from the reported measurement errors of the properties.

Statistical significance and scatter in K-band luminosity. In the companion paper, we show that the posterior constraints on the intrinsic scatter in $L_{K,\rm tot}$  are not bounded from below; values near zero are not only allowed by the data but the modal value of the posterior Probability Density Function (PDF) is zero. The correlation coefficients between $L_{K,\rm tot}$ and other properties vary substantially as the scatter in \ln\,$$L_{K,\rm tot} drops to very low values. Very small values of this scatter, $\sigma_{\ln L_{K}|M}$, are not physically reasonable. Cosmological hydrodynamics simulations have found values of $\sigma_{\ln L_{K}|M}>0.10$[11], $0.32$[10], or $0.16$[40]; and a recent observational study estimates a value of $0.22\pm 0.04$[41].

The confidence intervals and statistical significance of the anti-correlation signals reported in the main text employ a lower limit of $\sigma_{\ln L_{K}|M}=0.05$. This choice is a conservative one, two times smaller than the smallest value reported above. We discard any point in the posterior chain with $\sigma_{\ln L_{K}|M}<0.05$. All numbers reported in the main text are based on this truncated posterior distribution. For the sake of symmetry, we also impose the same limit on the richness scatter, $\sigma_{\ln\lambda|M}\geq 0.05$, but this has a much smaller effect as the posterior PDF has very little support in this region.

The statistical significance of the hot-cold baryon phase anti-correlation reported here is sensitive to the choice of minimum value for the stellar mass scatter. Figure 3 illustrates the odds of having a positive correlation for an optical and $L_{X}$ observable changes as a function of the imposed minimum value of $\sigma_{\ln M_{\rm star}|M}$. For a minimum value of 0.1, the odds of a positive correlation between $L_{K}$ and $L_{X}$ at fixed halo mass are 0.006, or roughly three-$\sigma$ evidence. The odds that both optical measures correlate positively with $L_{X,\rm ce}$ is very small, $0.005$ for our fiducial minimum of $0.05$ in $\sigma_{\ln M_{\rm star}|M}$.

Data availability

The observational data vector employed in this work is available via the companion paper[19] [https://academic.oup.com/mnras/article/484/1/60/5274143]. The full posterior chains that support the findings of this study is available in a figshare repository, [https://doi.org/10.6084/m9.figshare.8001218]. The source data underlying Figs 1, 2, and 3 and Table 2 are reproducible via the posterior chains.

Code Availability

All analyses were performed with custom-built Python scripts. The MCMC is performed with PyMC[42] [https://pymc-devs.github.io/pymc/], a Bayesian stochastic modeling package in Python. The analyses scripts are available from the corresponding author upon request.

Acknowledgements

We are grateful to Adam Mantz and Ben Maughan for useful discussions on population statistics. We thank Arif Babul and members of the LoCuSS, DES and XMM-XXL collaborations for their support. A.F. is supported by a McWilliams Postdoctoral Fellowship. A.E.E. and A.F. acknowledge support from NASA Chandra grants GO8-19107B and GO6-17116B. S.L.M. and G.P.S. acknowledge support from the STFC. G.P.S. acknowledges support from the Royal Society. C.O., D.P.M., Z.A., J.E.C., and CARMA operations were supported by NSF grant AST-1140019. C.P.H. acknowledges support from PRIN INAF 2014. H.B. and P.M. acknowledge financial contribution from the agreement ASI-INAF n.2017-14-H.O and from University of Rome “Tor Vergata” Grant, “Mission: Sustainability” EnClOS.

Author contributions statement

A.F. led the data analysis with contributions from A.E.E., S.L.M. G.P.S., and A.F. S.L.M. and G.P.S led the multi-wavelength observational studies. A.F. implemented, tested, and performed the statistical analysis. A.F. and A.E.E wrote the main body of the manuscript. H.B., R.M. and P.M. contributed in deriving all the observational quantity related to the X-ray observations of the hot gas content and the SZ signal with the Planck satellite. C.P.H. was responsible for the K-band photometry and cluster galaxy selection. N.O. was responsible for measuring the weak-lensing masses. D.P.M. and J.E.C. were responsible for obtaining the SZA observations and D.P.M supervised their analysis. C.O. analyzed interferometric data from SZA to estimate the $Y_{SZ}$ for each cluster. All authors discussed the results, contributed to the discussions in the manuscript, and reviewed the manuscript. G.P.S. proposed and designed the LoCuSS (the principle investigator).

Competing interests

The authors declare no competing interests.

Supplementary Information

Systematic Effects

In this section, we study the systematic effects to support the claims made in the main text. We are primarily interested in constraining the correlation coefficient between hot gas mass and stellar mass of the underlaying halos population, which is expected to be anti-correlated. We examine here sources of bias which could modify the estimated correlation. Throughout this section we assume power-law for the mean relation between mass and observables with log-Normal scatter, unless otherwise mentioned.

Effect of biased scaling relation slope

One potential systematic error arises from a potential bias in the posterior mean scaling relation parameters relative to the underlying truth. To illustrate this effect, we generate a synthetic sample of clusters with a known input scaling parameters, then estimate the correlation coefficient imposing a biased set of parameters. To start, we take the LoCuSS weak-lensing masses and assume a set of $M_{\rm gas}$–$M_{\rm WL}$ and $L_{K,\rm tot}$–$M_{\rm WL}$ relations. For each cluster, we draw a random $M_{\rm gas}$ and $L_{K,\rm tot}$ from a multivariate log-normal distribution with a set of input correlation coefficients and $20\%$ intrinsic scatter. Cluster residuals about the mean relation are estimated assuming a set of biased scaling relations consistent with ref.[19]. Finally, the correlation coefficient is estimated with Equation 3. For each input correlation coefficient, 1,000 realizations of LoCuSS-like cluster sample are generated. Figure 4 illustrates the shift in the estimated correlation coefficient as a function of input correlation coefficient.

While the correlation coefficients based on the gas mass, $M_{\rm gas}$, support negative values for both $L_{K,\rm tot}$ and $\lambda$, but the statistical significance is smaller than X-ray luminosity, $L_{X,\rm ce}$, and stellar mass proxies. In the companion paper, we note that the slope of the $M_{\rm gas}$ scaling with halo mass is $0.77\pm 0.1$, lower than values above one derived by previous observational studies[16] and from modern hydrodynamic simulations[11], which yield values slightly above unity. The purple line in Figure 4 shows that a bias of 0.3 in ${M}_{\rm gas}$ slope would produce an underestimate in the anti-correlation magnitude, potentially helping to explain why $L_{X,\rm ce}$ provides more significant evidence for anti-correlation than $M_{\rm gas}$.

Effect of extrinsic or underestimated statistical errors

Our likelihood model assumes that statistical errors of each property are accurate and that the remaining residuals about the mean scaling relation reflect only intrinsic scatter of the underlying halo population that host the clusters in the LoCuSS sample. Properties that are subject to extrinsic contributions not already incorporated into the statistical error budget represent another potential source of systematic error. For example, a recent study shows that the uncertainties quoted for $\lambda$ by the redMaPPer algorithm underestimate the extrinsic scatter driven by projected galaxies lying outside the host halo[27].

We test the effect of such systematics by realizing synthetic samples of pairs of properties, {$M_{\rm star}$, $M_{\rm gas}$}, subject to additional, uncorrelated sources of uncertainty. We assume $M_{\rm star}$–$M_{\rm WL}$ and $M_{\rm gas}$–$M_{\rm WL}$ relations with unit mean slopes and a two-dimensional, log-normal probability distribution with $20\%$ intrinsic scatter in each component and an input correlation coefficient, $r$. Each halo property is then further perturbed with uncorrelated scatter in each component, of variable magnitude $\sigma_{\rm star}$ and $\sigma_{\rm gas}$, to obtain the observed quantities. For each input correlation, we generate 1,000 realizations of 41-object samples, measuring the property correlation between the observed quantities with our likelihood model. Figure 5 illustrates that such extrinsic scatter leads to inferred correlations that underestimate the magnitude of the underlying population. We note that the difference in posterior correlation coefficient estimates for the $\{{M}_{\rm gas},$ L_{K,\rm tot}$$\} pair and $\{{M}_{\rm gas},\lambda\}$ pair seen in Figure 1 could be understood if $\lambda$ is subject to a larger extrinsic scatter effect than $L_{K,\rm tot}$.

The effects of bias in inferred slope and uncalibrated scatter can explain the systematic shifts in the marginalized posterior distributions seen in Figure 1. With current data, we are not able to assess the contribution of each systematic effect. For example, whether the richness is a noisier measurement of true halo stellar mass with respect to $L_{K,\rm tot}$, or whether the richness measurement uncertainties are underestimated, or both, cannot be addressed with the existing observational data. Detailed simulation and future observational studies are needed to calibrate and understand these effects.

Orientation effects

Massive halos deviate from spherical symmetry, with collisionless components tending toward prolate ellipsoidal shapes with minor to major axis ratio of 0.6 [43]. The collisional nature of the hot gas drives it to a more spherical shape formed by equipotential surfaces. It is reasonable to ask whether correlations between weak-lensing mass and stellar properties driven by orientation might induce an anticorrelation with hot gas properties. To address this question, we apply the log-normal population model of [20] under the assumption of moderate to strong covariance between weak-lensing mass and total K-band light at fixed true mass, and zero intrinsic correlation of either of these properties with hot gas properties. We show here that the correlation of hot gas and K-band luminosity conditioned on weak-lensing mass can be driven to negative values, but only under extreme values of model parameters that are very unlikely.

We consider the likelihood of a vector of properties, $p(\mathbf{s}\,|\,M_{\rm true},z)$, with elements, $\mathbf{s}=\ln\{M_{\rm WL},L_{k},X\}$ ($X=$ a hot gas property such as $L_{X},M_{\rm gas}$, etc.) and let $\mu=\ln M_{\rm true}$. The PDF of true selected mass, $p(\mu\,|\,M_{\rm WL})$, is Gaussian with an assumed width $0.2$. Considering the optical properties, $L_{K}$ or $\lambda$, for a sample selected on weak-lensing mass, the correlation between $\lambda$ (for example) and true halo mass, $\mu$, at fixed weak-lensing mass is given by equation (8) of ref.[20],

[TABLE]

where $\sigma_{\mu\,|\,M_{\rm WL}}=0.2$ is the scatter in true mass at fixed weak-lensing mass, $\sigma_{\mu\,|\,\lambda}$ is the similar measure at fixed $\lambda$, and $r$ is the correlation coefficient of $\ln M_{\rm WL}$ and $\ln\lambda$ at fixed $M_{\rm true}$. Similar expressions apply for $\lambda\rightarrow L_{K}$. We assume that the correlation coefficient of weak-lensing-mass and richness is significant, $r\gtrsim 0.5$. The above expression indicates that an anti-correlation between true mass and $\lambda$ or $L_{K}$ could be induced if the numerator is negative. With the default assumption that a hot gas property, $X$, is uncorrelated with $M_{\rm WL}$ and $\lambda$ at fixed true mass, $\mu$, this condition would produce an anti-correlation of $\lambda$ and $X$ at fixed $M_{\rm WL}$, consistent with our findings.

However, the requirement that Equation 4 be negative implies that the mass scatter at fixed K-band total luminosity must be large in order to suppress the first term in the numerator,

[TABLE]

A previous work that studied orientation and projection effects of different mass proxies in N-body simulations finds $r=0.55$ for weak-lensing and red galaxy count[22], meaning $\sigma_{\mu\,|\,\lambda}\simeq 0.4$ .

For this effect to be responsible for our findings of a large anti-correlation coefficient between hot and cold gas properties, two factors must conspire: i) $M_{\rm WL}$ and $L_{K}$ (or $\lambda$) at fixed true mass must be very tightly tied, $r\sim 0.9$, and, ii) both optical proxies must have a factor $\sim 2$ larger scatter in selected true mass compared to $M_{\rm WL}$. Such a strong correlation between weak-lensing mass and optical richness is unlikely, given that random LSS projections across a broad redshift range contribute at least $0.06$ to the scatter in $M_{\rm WL}$ [44]. Such a large intrinsic scatter in optical richness for halos above $5\times 10^{14}\,{\rm M}_{\odot}$, above which the majority of LoCuSS clusters lie, is not supported by state-of-the-art hydrodynamics simulations[11] and existing observational studies[45].

The following section includes an explicit test that confirms the required conditions.

Line-of-sight projection effect

The cluster properties we use are integrated quantities in the projected space of sky coordinates and photon frequency. Projection of extrinsic material will perturb the intrinsic quantities and, unlike the extrinsic study above, can do so in a correlated manner[22]. Here we model the effect of projection explicitly as a source of extrinsic covariance, $\Sigma_{\rm total}=\Sigma_{\rm halo}+\Sigma_{\rm proj}$. The covariance estimated in this work, Table 2, is the total covariance, and our findings in the main text implicitly assume $\Sigma_{\rm total}\approx\Sigma_{\rm halo}$. Because $\Sigma_{\rm halo}$ and $\Sigma_{\rm proj}$ are degenerate, we cannot disentangle these two components with the observational data in hand.

We perform a set of explicit simulations similar to what we have done above. We generate synthetic realizations of intrinsic quantities of a set of halos with an input property covariance, specified in Table 3 (a), and four models for projection covariance specified in Table 3 (b)-(e). Model (b) is motivated by the fact that the weak-lensing mass and optical observables are expected to be highly correlated [22]. Model (c) is similar to model (b) with larger extrinsic scatter, which illustrates the diluting effect of such a large scatter. It is expected that the optical and weak-lensing observables to be correlated with hot gas observables, but smaller in magnitude [22]. Thus, we consider a third model – model (d) – which includes a small correlation between the X-ray observable and the optical and the mass observables. Finally, model (e) is motivated by our analytical model of the previous section. This improbable model shows that a large extrinsic optical scatter strongly coupled with the mass proxy may indeed induce negative covariance.

Shaded areas are $68\%$ confidence intervals for inferred correlation derived from 1,000 realizations for each input correlation coefficient. The set examines a broad range of scenarios for the projection effect, and our results confirm that negative correlations are diluted. Unless the scatter due to the 2-halo term (projection effect) and correlation between $M_{\rm WL}$ and $L_{K}$ is very large the projection cannot induce a negative correlation.

As a final test, we revise our inference algorithm to include an explicit projection term with redefining the covariance $\Sigma=\Sigma_{\rm int}+\Sigma_{\rm proj}$. Due to the degeneracy between $\Sigma_{\rm int}$ and $\Sigma_{\rm proj}$, we cannot constrain both quantities simultaneously. We, therefore, employ a fixed $\Sigma_{\rm proj}$, model parameters specified in Table 3 (b), and preform the likelihood to infer $\Sigma_{\rm int}$. We find that $\Sigma_{\rm int}$ posteriors are consistent with our main results. Saying that, it is worth noting that the statistical significance of the anti-correlation between gas mass and K-band luminosity is improved – new p-value $=0.03$, and new $r=-0.70^{+0.34}_{-0.20}$.

We do not attempt to correct for this effect as an estimation of such effect is poorly constrained. Estimating the magnitude of this effect requires comprehensive numerical simulations, making synthetic observations, and performing the measurement algorithms which is beyond the scope of this work and is a subject of our future studies. We note that for a larger sample with better statistical power, this effect need to be calibrated and better understood.

Systematics Test Conclusion

The main conclusion drawn from the above tests is that these systematics cannot not induce an anti-correlation signal. Instead, any systematics which is not accounted for or exists in our dataset would reduce the statistical significance of our results and dilute the intrinsic halo anti-correlation signal. Therefore, our results should be interpreted as a conservative lower bound on the magnitude of the anti-correlation between ${M}_{\rm gas}$ and ${M}_{\rm star}$ at fixed halo mass. As the quality of data improving and the sample size growing, we need to handle the above systematics better. A future direction should, particularly, address the projection effect for a broad range of observables.