The XXL Survey: XXXVIII. Scatters and correlations of X-ray proxies in the bright XXL cluster sample

TL;DR

This study measures the intrinsic scatter and correlations among X-ray properties of galaxy clusters, revealing that gas mass is the least scattered proxy and that properties are positively correlated, consistent with self-similar evolution.

Contribution

It introduces a self-consistent method to measure the scatter covariance matrix of X-ray cluster properties without external calibration or mass knowledge.

Findings

Gas mass has the lowest intrinsic scatter (~8%).

Temperature is less scattered (~20%) but has larger measurement uncertainties.

X-ray properties are positively correlated, influenced by dynamical state and merger history.

Abstract

Context. Scaling relations between cluster properties embody the formation and evolution of cosmic structure. Intrinsic scatters and correlations between X-ray properties are determined from merger history, baryonic processes, and dynamical state. Aims. We look for an unbiased measurement of the scatter covariance matrix between the three main X-ray observable quantities attainable in large X-ray surveys -- temperature, luminosity, and gas mass. This also gives us the cluster property with the lowest conditional intrinsic scatter at fixed mass. Methods. Intrinsic scatters and correlations can be measured under the assumption that the observable properties of the intra-cluster medium hosted in clusters are log-normally distributed around power-law scaling relations. The proposed method is self-consistent, based on minimal assumptions, and requires neither the external calibration by weak lensing, dynamical, or hydrostatic masses nor the knowledge of the mass completeness. Results. We analyzed the 100 brightest clusters detected in the XXL Survey and their X-ray properties measured within a fixed radius of 300 kpc. The gas mass is the less scattered proxy (~8%). The temperature (~20%) is intrinsically less scattered than the luminosity (~30%) but it is measured with a larger observational uncertainty. We found some evidence that gas mass, temperature and luminosity are positively correlated. Time-evolutions are in agreement with the self-similar scenario, but the luminosity-temperature and the gas mass-temperature relations are steeper. Conclusions. Positive correlations between X-ray properties can be determined by the dynamical state and the merger history of the halos. The slopes of the scaling relations are affected by radiative processes.

Tables

Table 1. Table 1: Observed scaling relations. Conventions and units are as in Section  2 . Cols. 1-2: variables of the regression procedure. Cols. 3, 4, and 5: intercept, slope, and time evolution of the scaling relation. Cols. 6-7: scatter of X 𝑋 X and its time-evolution. Cols. 8-9: scatter of Y 𝑌 Y and its time-evolution. Cols. 10-11: correlation between the scatters and its time evolution. If the evolution is not considered ( γ σ X | Z = γ σ Y | Z = γ ρ X ​ Y | Z = 0 subscript 𝛾 subscript 𝜎 conditional 𝑋 𝑍 subscript 𝛾 subscript 𝜎 conditional 𝑌 𝑍 subscript 𝛾 subscript 𝜌 conditional 𝑋 𝑌 𝑍 0 \gamma_{\sigma_{X|Z}}=\gamma_{\sigma_{Y|Z}}=\gamma_{\rho_{XY|Z}}=0 ), scatters and correlations are meant as weighted averages over the redshift range. If the evolution is free, scatters and correlations are intended as local values at z ref = 0.01 subscript 𝑧 ref 0.01 z_{\text{ref}}=0.01 . Col. 12: self-similar slope ( β ss subscript 𝛽 ss \beta_{\text{ss}} ) for observable measured within fixed radii. Values in square brackets correspond to parameters kept fixed in the regression. We report the medians of the marginalized posterior distributions and the 68.3% probability range. Quoted values are biweight estimators of the marginalized distributions.

		Scaling						Intrinsic scatter in $X$				Intrinsic scatter in $Y$				Scatter correlation
		Intercept		Slope		Evolution		Local scatter		Evolution		Local scatter		Evolution		Local correlation		Evolution
$x$	$y$	${\alpha }_{Y|Z}$		${\beta }_{Y|Z}$		${\gamma }_{Y|Z}$		${\sigma }_{X|Z,0}$		${\gamma }_{{\sigma }_{X|Z}}$		${\sigma }_{Y|Z,0}$		${\gamma }_{{\sigma }_{Y|Z}}$		${\rho }_{XY|Z,0}$		${\gamma }_{{\rho }_{XY|Z}}$		${\beta }_{\mathrm{ss}}$
$t$	$l$	$-$2.15$\pm$	0.18	2.78$\pm$	0.42	$-$0.48$\pm$	1.25	0.03$\pm$	0.04	[0]		0.29$\pm$	0.08	[0]		0.20$\pm$	0.48	[0]		2
$t$	$l$	$-$2.13$\pm$	0.18	2.70$\pm$	0.42	$-$0.50$\pm$	1.29	0.03$\pm$	0.03	0.39$\pm$	0.34	0.29$\pm$	0.08	0.28$\pm$	0.28	0.28$\pm$	0.48	0.32$\pm$	0.66	2
$t$	$m_{\text{g}}$	$-$2.29$\pm$	0.06	1.53$\pm$	0.20	0.53$\pm$	0.53	0.10$\pm$	0.02	[0]		0.03$\pm$	0.03	[0]		0.64$\pm$	0.48	[0]		1
$t$	$m_{\text{g}}$	$-$2.27$\pm$	0.05	1.47$\pm$	0.17	0.67$\pm$	0.52	0.09$\pm$	0.02	0.23$\pm$	0.25	0.02$\pm$	0.02	0.45$\pm$	0.34	0.46$\pm$	0.52	$-$0.05$\pm$	0.63	1
$m_{\text{g}}$	$l$	2.32$\pm$	0.18	1.91$\pm$	0.09	$-$0.44$\pm$	0.49	0.04$\pm$	0.02	[0]		0.12$\pm$	0.05	[0]		0.89$\pm$	0.20	[0]		2
$m_{\text{g}}$	$l$	2.43$\pm$	0.14	1.97$\pm$	0.07	$-$0.73$\pm$	0.42	0.02$\pm$	0.01	0.42$\pm$	0.33	0.06$\pm$	0.03	0.42$\pm$	0.33	0.31$\pm$	0.58	0.04$\pm$	0.63	2

Table 2. Table 2: Scaling relations and intrinsic scatters from the multi-response regression as described in Section  3 . The gas mass m g subscript 𝑚 g m_{\text{g}} , the temperature t 𝑡 t , and the luminosity l 𝑙 l act as y 1 subscript 𝑦 1 y_{1} , y 2 subscript 𝑦 2 y_{2} , and y 3 subscript 𝑦 3 y_{3} , respectively. Y 1 subscript 𝑌 1 Y_{1} is an unbiased but scattered proxy of Z 𝑍 Z . Scatters and correlations are meant as weighted averages over the redshift range. Values in square brackets are kept fixed in the regression. The last column report the self-similar expectations. Quoted values are biweight estimators of the marginalized distributions.

Description	Relation	Parameter	Results		S-S
intercept	$M_{\text{gas,300kpc}}$-$M_{\text{gas,300kpc}}$	${\alpha }_{m_{\text{g}}^{\left(Y\right)}|m_{\text{g}}^{\left(Z\right)}}$	[0]		–
intercept	$T_{\text{300kpc}}$-$M_{\text{gas,300kpc}}$	${\alpha }_{t^{\left(Y\right)}|m_{\text{g}}^{\left(Z\right)}}$	1.60$\pm$	0.17	–
intercept	$L_{\text{300kpc}}^{\text{XXL}}$-$M_{\text{gas,300kpc}}$	${\alpha }_{l^{\left(Y\right)}|m_{\text{g}}^{\left(Z\right)}}$	2.39$\pm$	0.13	–
slope	$M_{\text{gas,300kpc}}$-$M_{\text{gas,300kpc}}$	${\beta }_{m_{\text{g}}^{\left(Y\right)}|m_{\text{g}}^{\left(Z\right)}}$	[1]		1
slope	$T_{\text{300kpc}}$-$M_{\text{gas,300kpc}}$	${\beta }_{t^{\left(Y\right)}|m_{\text{g}}^{\left(Z\right)}}$	0.71$\pm$	0.09	1
slope	$L_{\text{300kpc}}^{\text{XXL}}$-$M_{\text{gas,300kpc}}$	${\beta }_{l^{\left(Y\right)}|m_{\text{g}}^{\left(Z\right)}}$	1.94$\pm$	0.07	2
time-evolution	$M_{\text{gas,300kpc}}$-$M_{\text{gas,300kpc}}$	${\gamma }_{m_{\text{g}}^{\left(Y\right)}|m_{\text{g}}^{\left(Z\right)}}$	[0]		0
time-evolution	$T_{\text{300kpc}}$-$M_{\text{gas,300kpc}}$	${\gamma }_{t^{\left(Y\right)}|m_{\text{g}}^{\left(Z\right)}}$	0.68$\pm$	0.09	0
time-evolution	$L_{\text{300kpc}}^{\text{XXL}}$-$M_{\text{gas,300kpc}}$	${\gamma }_{l^{\left(Y\right)}|m_{\text{g}}^{\left(Z\right)}}$	0.15$\pm$	0.55	0
intrinsic scatter	$M_{\text{gas,300kpc}}$—$Z$	${\sigma }_{m_{\text{g}}^{\left(Y\right)}|Z}$	0.03$\pm$	0.01	–
intrinsic scatter	$T_{\text{300kpc}}$—$Z$	${\sigma }_{t^{\left(Y\right)}|Z}$	0.09$\pm$	0.02	–
intrinsic scatter	$L_{\text{300kpc}}^{\text{XXL}}$—$Z$	${\sigma }_{l^{\left(Y\right)}|Z}$	0.05$\pm$	0.03	–
scatter correlation	$T_{\text{300kpc}}$-$M_{\text{gas,300kpc}}$— $Z$	${\rho }_{m_{\text{g}}^{\left(Y\right)}{t}^{\left(Y\right)}|Z}$	0.35$\pm$	0.52	–
scatter correlation	$L_{\text{300kpc}}^{\text{XXL}}$-$M_{\text{gas,300kpc}}$— $Z$	${\rho }_{m_{\text{g}}^{\left(Y\right)}{l}^{\left(Y\right)}|Z}$	0.40$\pm$	0.43	–
scatter correlation	$L_{\text{300kpc}}^{\text{XXL}}$-$T_{\text{300kpc}}$— $Z$	${\rho }_{t^{\left(Y\right)}{l}^{\left(Y\right)}|Z}$	0.07$\pm$	0.70	–

Table 3. Table 3: Intrinsic distributions of the observed samples modeled as Gaussian functions. Conventions and units are as in Sec.  2 . Quoted values are biweight estimators of the marginalized distributions.

	Mean						Dispersion
Sample	${\mu }_{Z,0}$		${\gamma }_{{\mu }_Z,D}$		${\gamma }_{{\mu }_Z,F_z}$		${\sigma }_{Z,0}$		${\gamma }_{{\sigma }_Z}$
$l$ -$t$	$-$0.22$\pm$	0.11	0.44$\pm$	0.10	$-$0.43$\pm$	0.56	0.15$\pm$	0.02	$-$0.72$\pm$	0.32
$m_{\text{g}}$-$t$	$-$0.25$\pm$	0.10	0.46$\pm$	0.08	$-$0.32$\pm$	0.49	0.10$\pm$	0.02	$-$0.80$\pm$	0.23
$l$ -$m_{\text{g}}$	$-$2.68$\pm$	0.12	0.70$\pm$	0.10	0.05$\pm$	0.54	0.16$\pm$	0.02	$-$0.81$\pm$	0.22

Equations

$$X_Z$$

$$Y_Z$$

$$P(x_j,y_j|X_j,Y_j)\propto {\mathcal{N}}^{\text{(2)}}(\{X_j,Y_j\},{\mathbf{V}}_{\delta ,j})\mathcal{U}(y_{\text{th},j},\infty ),$$

$$P(y_{\text{th},j}|y_{\text{th,obs},j})=\mathcal{N}(y_{\text{th,obs},j},{\delta }_{y_{\text{th},j}}^2),$$

$$P(X_j,Y_j|X_{Z,j},Y_{Z,j})\propto {\mathcal{N}}^{\text{(2)}}(\{X_{Z,j},Y_{Z,j}\},{\mathbf{V}}_{\sigma ,j})\mathcal{U}(Y_{\text{th},j},\infty ),$$

$$P(Y_{\text{th},j}|y_{\text{th},j})=\mathcal{N}(y_{\text{th},j},{\delta }_{y,j}^2).$$

$${\sigma }_{X|Z}(z)$$

$${\sigma }_{Y|Z}(z)$$

$${\rho }_{XY|Z}(z)$$

$$P(Z)=\mathcal{N}\left({\mu }_Z,{\sigma }_Z^2\right).$$

$${\mu }_Z(z)={\mu }_{Z,0}+{\gamma }_{{\mu }_Z}\log F_z+{\gamma }_{{\mu }_Z,D}\log D_z,$$

$${\sigma }_Z(z)={\sigma }_{Z,0}{F}_z^{{\gamma }_{{\sigma }_Z}}.$$

$${\alpha }_{Y|Z},{\mu }_{Z,0}\sim \mathcal{U}(-1/ϵ,1/ϵ),$$

$${\rho }_{XY|Z,0}\sim \mathcal{U}(-1,1).$$

$${\beta }_{Y|Z},{\gamma }_{Y|Z},{\gamma }_{{\mu }_Z},{\gamma }_{{\mu }_Z,D}\sim t_1.$$

$${\gamma }_{{\sigma }_{X|Z}},{\gamma }_{{\sigma }_{Y|Z}}\sim \mathcal{U}(0,1),$$

$${\gamma }_{{\rho }_{XY|Z}}\sim \mathcal{U}(-1,1).$$

$${\gamma }_{{\sigma }_Z}\sim \mathcal{U}(-1,1).$$

$$1/{\sigma }_{X|Z,0}^2,1/{\sigma }_{Y|Z,0}^2,1/{\sigma }_{Z,0}^2\sim \Gamma (r,\lambda ),$$

$${\gamma }_{{\sigma }_{X|Z}}={\gamma }_{{\sigma }_{Y|Z}}={\gamma }_{{\rho }_{XY|Z}}=0.$$

$$Y_{Z,i}={\alpha }_{Y_i|Z}+{\beta }_{Y_i|Z}Z+{\gamma }_{Y_i|Z}\log F_z.$$

$$P(Y_{1,j},Y_{2,j},....|Y_{Z,1j},Y_{Z,2j},....)\propto {\mathcal{N}}^{(n)}\left(\{Y_{Z,1j},Y_{Z,2j},....\},{\mathbf{V}}_{\sigma }\right)$$

$$\times \prod _i\mathcal{U}(Y_{\text{th},ij},\infty ),$$

$$P(y_{1,j},y_{2,j},...|Y_{1,j},Y_{2,j},...)$$

$$\times \prod _i\mathcal{U}(y_{\text{th},ij},\infty ),$$

$$P(y_{\text{th}ij}|y_{\text{th,obs,}ij})=\mathcal{N}\left(y_{\text{th,obs,}ij},{\delta }_{y_{\text{th,ij}}}^2\right),$$

$$P(Y_{\text{th,ij}}|y_{\text{th,ij}})=\mathcal{N}\left(y_{\text{th,ij}},{\delta }_{y,ij}^2\right),$$

$${\mathbf{V}}_{\sigma }^{-1}\sim \mathbf{W}(\mathbf{S},d),$$

$${\mathbf{S}}_{aa}\sim \Gamma (ϵ,ϵ).$$

$$m_{\text{g}}$$

$$t$$

$$l$$

$$y_{\text{th},j}=l_{\text{th},j}.$$

$$M_{\text{g}}\propto M.$$

Full text

*Aims. *We look for an unbiased measurement of the scatter covariance matrix between the three main X-ray observable quantities attainable in large X-ray surveys – temperature, luminosity, and gas mass. This also gives us the cluster property with the lowest conditional intrinsic scatter at fixed mass.

*Methods. *Intrinsic scatters and correlations can be measured under the assumption that the observable properties of the intra-cluster medium hosted in clusters are log-normally distributed around power-law scaling relations. The proposed method is self-consistent, based on minimal assumptions, and requires neither the external calibration by weak lensing, dynamical, or hydrostatic masses nor the knowledge of the mass completeness.

*Results. *We analyzed the 100 brightest clusters detected in the XXL Survey and their X-ray properties measured within a fixed radius of 300 kpc. The gas mass is the less scattered proxy ($\sim 8\%$). The temperature ($\sim 20\%$) is intrinsically less scattered than the luminosity ($\sim 30\%$) but it is measured with a larger observational uncertainty. We found some evidence that gas mass, temperature and luminosity are positively correlated. Time-evolutions are in agreement with the self-similar scenario, but the luminosity-temperature and the gas mass-temperature relations are steeper.

*Conclusions. *Positive correlations between X-ray properties can be determined by the dynamical state and the merger history of the halos. The slopes of the scaling relations are affected by radiative processes.

1 Introduction

The physics of baryons and dark matter can be assessed with scaling relations between cluster properties (Pratt et al., 2009; Arnaud et al., 2010; Giodini et al., 2013). Ongoing and upcoming large surveys are measuring a wealth of cluster properties, e.g. optical richness, X-ray luminosity, Sunyaev-Zel’dovich (SZ) flux (Laureijs et al., 2011; Planck Collaboration et al., 2014a; Bleem et al., 2015; Pierre et al., 2016; Melchior et al., 2017; Maturi et al., 2019).

Gravity is the driving force in structure formation and evolution, and makes clusters self-similar with observable properties following power-law relations in halo mass (Kaiser, 1986; Giodini et al., 2013; Ettori, 2013). Deviations from the self-similar scheme are due to non-gravitational processes, such as feedback and non-thermal mechanisms, which can contribute significantly to the global energy budget (Maughan et al., 2012).

The scaling relations are scattered by underlying processes that can affect different cluster properties to different degrees (Stanek et al., 2010; Truong et al., 2018). Numerical simulations (Stanek et al., 2010; Fabjan et al., 2011; Angulo et al., 2012; Saro et al., 2013) and observational studies (Maughan, 2007; Vikhlinin et al., 2009a) confirm that the properties are log-normally distributed. Broadly speaking, the scatter is related to the regularity of the clusters (Sereno & Ettori, 2015b, a) and to deviations from equilibrium (Fabjan et al., 2011; Saro et al., 2013). Well behaved proxies with small scatters can be used to provide accurate measurement of the mass of galaxy clusters, which is crucial in important branches of cosmology and astrophysics (Ettori et al., 2009; Vikhlinin et al., 2009b; Mantz et al., 2010; Planck Collaboration et al., 2014b).

If we measure scalings, scatters, and correlations of scaling relations, we can study the forces driving cluster formation and evolution. In this paper, we investigate X-ray properties measurable in large surveys. We propose a novel statistical method where scatters are measured exploiting the expected linearity of the relations even without knowing the mass of the clusters. An observable property of a galaxy cluster can be a good proxy if it is easy to measure and well behaved. Intrinsic scatter enables us to determine such a variable. We can view the proxy with the lowest intrinsic conditional scatter with respect to some basic cluster characteristic as ‘optimal’.

We exploit the XXL Survey, the largest completed XMM-Newton project (Pierre et al., 2016, XXL Paper I). Two sky regions for a total of 50 square degrees has been surveyed and several hundreds of galaxy clusters out to redshift $\sim 2$ have been detected (Adami et al., 2018, XXL Paper XX). We consider the three main X-ray global quantities measured by the mission, i.e. temperature, luminosity, and gas mass.

The paper is as follows. The regression method and its extension to multi-response variables are described in Secs. 2 and 3, respectively. The data sample is introduced in Sec. 4. Theoretical expectations are briefly discussed in Sec. 5. Section 6 presents the results. Section 7 is devoted to comparison with previous analyses. Final considerations are in Sec. 8. Appendix A presents a simple recipe to deal with asymmetric errors and logarithmic variables. Additional figures are presented in App. B.

1.1 Notations

The frame-work cosmological model in use in the XXL papers is the flat $\Lambda$CDM universe with density parameter $\Omega_{\text{M}}=0.28$, and Hubble constant $H_{0}=70\leavevmode\nobreak\ \text{km\leavevmode\nobreak\ s}^{-1}\text{Mpc}^{-1}$, as found from the study of the final nine years cosmic microwave background observations of the Wilkinson Microwave Anisotropy Probe satellite (WMAP9), combined with baryon acoustic oscillation measurements and constraints on $H_{0}$ from Cepheids and type Ia supernovae (Hinshaw et al., 2013).

As usual, $H(z)$ is the redshift dependent Hubble parameter and $E_{z}\equiv H(z)/H_{0}$. When $H_{0}$ is not specified, $h$ is the Hubble constant in units of $100\leavevmode\nobreak\ \mathrm{km\leavevmode\nobreak\ s}^{-1}\mathrm{Mpc}^{-1}$.

$O_{\Delta}$ denotes a cluster property measured within the radius $r_{\Delta}$ which encloses a mean over-density of $\Delta$ times the critical density at the cluster redshift, $\rho_{\mathrm{cr}}=3H(z)^{2}/(8\pi G)$.

‘$\log$’ is the logarithm to base 10 and ‘$\ln$’ is the natural logarithm. Results for natural logarithm are quoted as percents, i.e. 100 times the dispersion in natural logarithm.

By intrinsic scatter we mean the standard deviation of the conditional probability, e.g., the probability of the temperature given the mass. If the conditional probability is related to the measurement process, i.e. the probability of the measured output given true input, we name it statistical measurement uncertainty. We do not name as scatter the standard deviation of the marginalized distribution, e.g. the distribution of temperatures. Throughout the paper, we denote the intrinsic scatter as $\sigma$ and the measurement uncertainty as $\delta$.

2 Regression scheme

In this section, we describe the Bayesian fitting procedure used to retrieve the scaling relations and the intrinsic scatters. We assume that the cluster properties (temperature, luminosity, gas mass) are power laws of some basic property, e.g. the mass. Hereafter, we focus on the logarithms of these quantities, which are thus linearly related with each other.

In a nutshell, we first appoint a basic intrinsic cluster feature (denoted by $Z$ as the reasoning would apply to the mass as well as to other choices). For any measurable property, e.g. the temperature, we distinguish three variables, see Fig. 1: i) $X_{Z}$, the latent quantity that is exactly linked to $Z$ through a functional relation $X(Z)$ (Maughan, 2014). ii) $X$, the true quantity that would be measured in an theoretical observation with infinite accuracy and precision (Feigelson & Babu, 2012). $X$ is intrinsically scattered with respect to $X_{Z}$. iii) $x$, the manifest result of the measurement process which shows some observational noise with respect to $X$.

Similarly to $X$, we can consider an additional cluster property $Y$ and the related $y$ and $Y_{Z}$. If the latent variables ($X_{Z}$, $Y_{Z}$) are linearly related to $Z$, they are linearly related with each other. We can then identify the best proxy with the X-ray observable characterized by the lowest conditional intrinsic scatter with respect to $Z$. In this way, the optimal proxy can be established even if we do not know the cluster mass.

Of course, we still need the mass if we want to calibrate the scaling relations and unambiguously identify the best proxy as the best mass proxy. In principle, $Z$ can be any fundamental property of the cluster. Since we do not measure the mass itself, our analysis is based on the fact that $Z$ is the quantity that can best characterize the X-ray properties of the clusters and minimize the scatters. Even though on a theoretical basis the most suitable candidate for this role is the mass, it is not guaranteed that $Z$ is the mass. It could be another quantity, e.g. the optical richness or the SZ signal, or one of the three X-ray properties measured in the present paper. However, this last hypothesis can be discarded if the estimated intrinsic scatters are not null.

2.1 Scaling and distributions

In order to quantify the intrinsic scatter of the X-ray properties, we follow the regression scheme detailed in the CoMaLit series (Comparing Masses in Literature, Sereno & Ettori, 2015b; Sereno et al., 2015; Sereno & Ettori, 2015a, 2017) and in Sereno (2016). This scheme accounts for time-evolution, correlated intrinsic scatters, and selection effects (Malmquist/Eddington biases).

In this section we consider a pair of observables. In the next section, we generalize the procedure to the case of multiple response variables.

As result of the measurement of the $j$-th cluster, the observable $x_{j}$, $y_{j}$, and the related uncertainty covariance matrix $\mathbf{V}_{\delta,j}$ are known. On the other hand, $\{X_{Z,j},Y_{Z,j}\}$, $\{X_{j},Y_{j}\}$, and the covariance matrix of the intrinsic scatters, $\mathbf{V}_{\sigma,j}$, are unknowns to be determined under the assumption of linearity.

In case of a linear relation, $X_{Z}$ and $Y_{Z}$ are related to the same covariate variable, $Z$, as

[TABLE]

where $\alpha$ denotes the normalization, the slope $\beta$ accounts for the dependence with $Z$, and the slope $\gamma$ accounts for the redshift-evolution. $F_{z}$ is the re-normalized Hubble parameter, $F_{z}=E_{z}/E_{z}(z_{\text{ref}})$. Here and in the following, we assume a power-law dependence on $F_{z}$ for the redshift-evolution of the observables.

The relations between $Z$, $X_{Z}$, and $Y_{Z}$ are deterministic and they are not affected by scatter. We assume that the uncertainty on the spectroscopic redshift $z$ is negligible.

If we do not know the value of $Z$, the two slopes and the two normalizations in Eqs. (1) and (2) are redundant. We then assume that $X$ is an unbiased proxy of $Z$, i.e. we fix $\alpha_{X|Z}=0$, $\beta_{X|Z}=1$, and $\gamma_{X|Z}=0$. Fixing the parameters of the $X$-$Z$ rather than the $Y$-$Z$ relation is just a matter of rescaling which does not affect the analysis of the intrinsic scatters. In absence of a direct measurement of $Z$, any bias between $X$ and $Z$ (i.e. $\alpha_{X|Z}\neq 0$) is degenerate with the estimated overall normalization of the scaling between $Y$ and $Z$. The data analysis can only constrain the relative bias between $X$ and $Y$ (Sereno & Ettori, 2015b, CoMaLit-I).

The measured and the true values of the quantities are related as

[TABLE]

where ${\cal N}^{\text{(2)}}$ is the bivariate Gaussian distribution and ${\cal U}$ is the uniform distribution. $\mathbf{V}_{\delta,j}$ is the uncertainty covariance matrix whose diagonal elements are denoted as $\delta_{x,j}^{2}$ and $\delta_{y,j}^{2}$, and whose off-diagonal elements are denoted as $\rho_{xy,j}\delta_{x,j}\delta_{y,j}$. The proportionality symbol in Eq. (3) indicates that the function on the right hand side is not normalized.

The probability distribution in Eq. (3) is truncated for $y_{j}<y_{\text{th},j}$, which accounts for selection effects when only clusters above an observational threshold (in the response variable) are included in the sample, i.e. the Malmquist bias (Sereno et al., 2015, CoMaLit-II).

To shorten the notation in Eq. (3) and similarly in the following, we drop the explicit dependence on the fixed parameters, i.e. $\mathbf{V}_{\delta,j}$ and $y_{\text{th},j}$, on the left hand side. On the right hand side, we do not express the functional dependence on the random variables $x_{j}$ and $y_{j}$.

The observational threshold $y_{\text{th}}$ may not be exactly known. This is the case when the quantity which the selection procedure is based on differs from the quantity used in the regression. As an example, XXL clusters are selected according to their flux within 1′ whereas we consider the luminosity within 0.3 Mpc in the regression, see Sec. 4. We have then to consider the additional relation

[TABLE]

where $\delta_{y_{\text{th},j}}$ is the uncertainty associated to the measured threshold $y_{\text{th,obs},j}$. This save us from adding new quantities, e.g. the observed flux, to the regression scheme.

We assume that the intrinsic scatters are Gaussian. It is

[TABLE]

where $\mathbf{V}_{\sigma,j}$ is the scatter covariance matrix whose diagonal elements are denoted as $\sigma_{X|Z,j}^{2}$ and $\sigma_{Y|Z,j}^{2}$, and whose off-diagonal elements are denoted as $\rho_{XY|Z,j}\sigma_{X|Z,j}\sigma_{Y|Z,j}$.

$Y_{\text{th},j}$ is the threshold in the response variable,

[TABLE]

Even though the selection procedure is based only on the value of the measured $y$ rather than the value of $Y$, any threshold in $y$ affects all the probability distributions. In fact, we do not sample a generic distribution of clusters but we select them and we have to model the distribution of the sampled objects. Hence, the distribution of $Y$ given $Z$ for a generic cluster from the full population differs from the distribution specific to a selected sample, which follows Eq. (5) and it is truncated.

We assume that the scatters are time-dependent, hence the subscript $j$ in the covariance matrix, but do not depend on $Z$ (Rozo et al., 2010). The time evolution of the intrinsic scatters and of the correlation can be modeled as

[TABLE]

The intrinsic distribution of $Z$ can be approximated with a mixture of Gaussian functions (Kelly 2007, CoMaLit-II, Sereno & Ettori 2015a, CoMaLit-IV). In the simplest but still effective case of one component (Sereno, 2016),

[TABLE]

Most of the parent populations of astronomical quantities, e.g. the halo mass function or the luminosity function, are locally exponential (in log-space), i.e. $P_{\text{parent}}(Z)\sim\exp(-aZ)$. However, here we have to model just the distribution of the clusters included in the sample rather than the full population. Once the parent population is filtered by the selection process, a Gaussian distribution provides a reliable approximation (CoMaLit-IV).

The redshift evolution of the (mean of the) $Z$-distribution can be modeled as (CoMaLit-IV),

[TABLE]

where $\mu_{Z,0}$ is the local mean and $D_{z}$ is the luminosity distance. We renormalize the distances such that $D_{z}(z_{\text{ref}})=1$.

The dispersion of the $Z$-distribution evolves as

[TABLE]

The dependence on $F_{z}$ is enough to account for the redshift evolution of the scaling relations, see Eqs. (1, 2) and of the scatters, see Eqs. (7–9) and Eq. (12). This is justified by theoretical predictions based on the self-similar model, by results of numerical simulations, and by observational fits (CoMaLit-IV). On the other hand, we introduce an explicit dependence on the cosmological distance for the evolution of the covariate distribution, see Eq. (11). In fact, the completeness of a sample selected according to either flux or signal-to-noise depends on the distance (CoMaLit-IV). The redshift dependence in Eq. (11) is general enough to address even more complicated cases.

2.2 Priors

Priors have to be conveniently non-informative (\al@se+et15_comalit_I,ser+al15_comalit_II; \al@se+et15_comalit_I,ser+al15_comalit_II). The priors on the intercept $\alpha_{Y|Z}$ and on the mean $\mu_{Z,0}$ are flat,

[TABLE]

where $\epsilon$ is a small number. In our calculation, $\epsilon=10^{-4}$. The prior on the correlation $\rho_{XY|Z,0}$ is flat too,

[TABLE]

We model the prior probability of the slopes and of the time-evolution as a Student’s $t_{1}$ distribution with one degree of freedom,

[TABLE]

This is equivalent to uniformly distributed direction angles $\arctan\beta$ and $\arctan\gamma$.

As non informative priors for the evolution of scatters and correlations, we consider uniform distributions. Since intrinsic scatters are expected to slightly increases with redshifts (CoMaLit-IV), we assume

[TABLE]

and

[TABLE]

For the evolution of the dispersion of the $Z$-distribution we adopt

[TABLE]

For the precision, i.e. the inverse of the variance, we adopt a nearly scale-invariant Gamma distribution,

[TABLE]

where the rate $r$ and the shape parameter $\lambda$ are fixed to $r=\lambda=\epsilon$ so that the prior spans a considerable range and it is nearly constant in logarithmic bins.

Alternatively to non-informative priors, we consider strong assumptions. The evolution of scatters and correlations is poorly constrained in samples of $\sim 100$ objects (Sereno, 2016). In our reference regression scheme, we only constrain the redshift weighted scatters and we fix

[TABLE]

3 Multi-response regression scheme

The scheme detailed in Section 2 can be generalized to the simultaneous regression of $n(\geq 2)$ observables. In this scheme, $y_{ij}$ is the $j$-th measurement of the $i$-th observable, $Y_{ij}$ is the true value, and $Y_{Z,ij}$ is the latent unscattered quantity which fits the scaling relation. To simplify the notation, in the following we dismiss the subscripts for the time dependence of the scatters. The scaling relation of the $i$-th property is expressed as

[TABLE]

Due to degeneracy, we anchor the scaling parameters of the first response variable, i.e. $\alpha_{Y_{1}|Z}=0$, $\beta_{Y_{1}|Z}=1$, and $\gamma_{Y_{1}|Z}=0$. The reference $Z$ variable is modeled as in Eq. (10). In absence of Malmquist biases, $Z$ is the expected values of $Y_{1}$given $Z$, $\langle Y_{1}|Z\rangle=Z$.

The intrinsic scatters shape the distribution of the true quantities around the model predictions. For the $j$-th cluster

[TABLE]

where ${\cal N}^{(n)}$ is the multivariate Gaussian distribution. $\mathbf{V}_{\sigma}$ is the $n\times n$ scatter covariance matrix whose diagonal elements are the intrinsic variances, $\sigma_{Y_{i}|Z}^{2}$, and whose off-diagonal elements can be expressed in terms of the correlations as $\rho_{Y_{a}Y_{b}|Z}\sigma_{Y_{a}|Z}\sigma_{Y_{b}|Z}$. The threshold for the $j$-th measurement of the $i$-th response variable is $Y_{\text{th},ij}$.

The measured and the true values are related as

[TABLE]

where $\mathbf{V}_{\delta,j}$ is the $n\times n$ uncertainty covariance matrix of the $j$-th cluster. The thresholds for the measured and the true response values are related as

[TABLE]

and

[TABLE]

where $y_{\text{th,obs,}ij}$ is the observational threshold, $\delta_{y_{\text{th,ij}}}$ is the related uncertainty, and $\delta_{y,ij}$ is the uncertainty associated to $y_{ij}$.

Similarly to Sec. 2, Malmquist bias is treated with the inclusion of the (smooth) truncations in Eqs. (22, 23).

3.1 Priors

We express the prior on the (inverse of the) intrinsic scatter matrix in terms of the Wishart distribution,

[TABLE]

where $d$ is the number of degrees of freedom and $\mathbf{S}$ in the $n\times n$ scale matrix. We take $d=n+1$, so that the marginalized prior distribution of the correlation factors is uniform between $-1$ and $1$. In analogy to the variances in Eq. (19), we model $\mathbf{S}$ as a scalar matrix with diagonal elements

[TABLE]

The Wishart prior is widely regarded as non-informative, even though it favors high variance in case of high correlation. Other priors are defined as in Sec. 2.2.

4 The sample

The XXL survey covers a total area of 50 square degrees with a X-ray sensitivity of $\sim 10^{-14}\mathrm{erg\leavevmode\nobreak\ cm^{-2}s^{-1}}$ in the [0.5-2] keV band for extended sources (Pacaud et al., 2016, XXL Paper II). The survey has uncovered more than 300 galaxy clusters out to redshift $\sim 2$ (XXL Paper XX) over a wide range of nearly two decades in mass (Lieu et al., 2016, XXL Paper IV). The XXL programme has the potential to constrain at the same time cluster scaling relations and cosmological parameters (Pacaud et al., 2018, XXL Paper XXV).

We here consider the sample of the 100 brightest clusters (XXL-100-GC) from the first data release (DR1) (DR1, XXL Paper I, ). The candidate clusters were selected by setting a lower limit of $3\times 10^{-14}\mathrm{erg\leavevmode\nobreak\ cm^{-2}s^{-1}}$ on the source flux in the soft X-ray band within a 1′ aperture (XXL Paper II). We analyse the spectroscopic temperature within 300 kpc, $T_{\mathrm{300kpc}}$, the luminosity in the rest-frame soft band [0.5-2.0] keV within 300 kpc, $L^{\mathrm{XXL}}_{\mathrm{300kpc}}$ (Giles et al., 2016, XXL Paper III) and the gas mass. Here, we consider the gas mass in a sphere of radius equal to 300 kpc, $M_{\mathrm{gas,300kpc}}$, which are computed following the procedure described in Eckert et al. (2016, XXL Paper XIII), which we refer to for details. The full set of measurements is available for 96 clusters. The sample has a mean redshift of $z=0.38$ with standard deviation of $\Delta z=0.24$. The sample covers an extended range in temperature, from small groups at $T_{\mathrm{300kpc}}\la 1\leavevmode\nobreak\ \text{keV}$ to rich clusters at $T_{\mathrm{300kpc}}\la 7\leavevmode\nobreak\ \text{keV}$, even though the most massive and rare clusters from the extreme tail of the cosmological halo mass function are absent due to the limited survey area coverage. For a comprehensive discussion of the sample properties we refer to XXL Paper II.

In the following, we recover the intrinsic scatters and correlations of the X-ray properties without having access to the mass. We use a simplified notation in the log-space. The measured gas mass, temperature, and luminosity in logarithmic units are

[TABLE]

For the latent variables, we cannot strictly follow the convention in Secs. 2 and 3, since X-ray observables in linear space are usually named by capital letters, e.g. $M_{\text{gas,300kpc}}$. Then, if $m_{\text{g}}$ is the measured gas mass ($y$), the hypothetical measurement in absence of noise ($Y$) is $m_{\text{g}}^{(Y)}$ and the unscattered gas mass ($Y_{Z}$) is $m_{\text{g}}^{(Z)}$. The same convention applies to temperature and luminosity.

For our analysis, we follow a data-driven approach and we consider as proxies only X-ray properties measured within fixed physical radii. The use of quantities measured in a given over-density radius as proxy can be ambiguous to some degree. If we know the over-density radius, we know by definition the mass too. Let $O_{\Delta}$ be a generic observable quantity, e.g. temperature or richness, within the over-density radius $r_{\Delta}$. In practice, we use only a part ($r_{\Delta}$) of the full information we already have (the mass $M_{\Delta}$) to get a deteriorated version, i.e. the mass proxy $M_{\Delta}(O_{\Delta})$ calculated through the scattered scaling relation applied to $O_{\Delta}$, of the main information itself (the mass $M_{\Delta}$). In this sense, we lose information. This can be corrected with iterative approaches by determining at the same time $r_{\Delta}$, $O_{\Delta}$, and $M_{\Delta}$. However, we have still to rely on very strong priors (usually, the knowledge of how $r_{\Delta}$ scales with some observable property). This has little effect if the observable is poorly correlated with the radius, e.g. the X-ray luminosity emitted from a very large area, but it is a major problem otherwise, e.g. the gas mass. Valuable mass proxies can be highly correlated with the integration radius.

The use of X-ray properties measured within fixed physical radii also minimize the impact of the assumed cosmological parameters on our results (Sereno & Ettori, 2015b). A different frame-work cosmology would mostly imply a different normalization and a slightly different time-dependence of the relations. Slopes and scatters would be minimally affected. At present, analyses of number counts (XXL Paper XXV) and clustering (Marulli et al., 2018, XXL Paper XVI) of XXL clusters are compatible with our reference cosmological model.

4.1 Covariance uncertainty matrix

The knowledge of the covariance uncertainty matrix is crucial to obtain unbiased estimates of the intrinsic scatters and of their correlations. Measurements of luminosity and temperature are based on the spectroscopic analysis of the core region and are largely independent of the gas mass measurement process, which exploits the photometry and the surface brightness profile in annular regions. However, the photons are the same and the gas mass measurement process uses the estimated temperature to convert the observed surface-brightness profiles into emission-measure profiles. This conversion is largely insensitive to the temperature and metallicity as long as the temperature exceeds $\sim$1.5 keV. In most cases the temperature and the gas measurement are nearly uncorrelated.

To estimate the uncertainty covariance matrix we proceed in the following way. Luminosity and temperature are estimated in a single measurement process. Their correlation is an output of the spectroscopic analysis. We approximate the probability distribution of the observed luminosity and temperature as a bivariate Gaussian.

Measured temperature and gas mass are correlated too. The estimate of the gas mass relies on the conversion of the observed surface brightness into the emission measure. The conversion factor is computed using $T_{\mathrm{300kpc}}$ and simulating a single-temperature absorbed thin-plasma model with the APEC code (XXL Paper XIII).

To estimate the full covariance uncertainty matrix, we extract $10^{5}$ couples of luminosity and temperature for each cluster from the approximated bivariate normal distribution. We then compute the new conversion factor from surface brightness into emission measure for the sampled temperature and we compute the corresponding gas mass by rescaling. We finally extract a new gas mass measurement from a normal distribution centered on the rescaled $M_{\mathrm{gas,300kpc}}$ and with standard deviation equal to the observational uncertainty $\delta M_{\mathrm{gas,300kpc}}$. The final correlation matrix is computed from the sampled values of temperature, luminosity, and gas mass.

4.2 Selection effects

Cosmological studies of number counts and abundance evolution require a very detailed study of the completeness. Observed properties have to be related to the underlying mass function and the selection function can be expressed in terms of the true cluster parameters rather than in terms of their measured counterparts affected by measurement errors.

The selection and validation of the XXL-100-GC sample is described in XXL Paper II. Here we just recall the main features relevant for our analysis. The XXL-100-GC sample was chosen with a flux-limit selection in a fixed angular aperture. The flux limit of the final catalog is $3\times 10^{-14}\mathrm{erg\leavevmode\nobreak\ s^{-1}cm^{-2}}$ in the [0.5-2.0] keV band within a 1′ radius aperture. Assuming temperature equal $3\leavevmode\nobreak\ \mathrm{keV}$, metallicity equal to $0.3$, and redshift $z=0.3$, this is equivalent to a total MOS1+MOS2+PN count rate of 0.0332 cts/s in [0.5-2] keV.

The C1+2 pipeline selection function of the XXL Survey depends on cluster profile, emissivity, and position (XXL Paper II). The dependence on the exposure and background level is encapsulated in the pointing under consideration and the off-axis angle is implicitly given by the sky position (XXL Paper II). Candidate clusters are then confirmed by visual inspection.

When we express the selection function in terms of the true intrinsic parameters, we have to consider that a cluster with a given true count rate can exceed the cut or fail depending on the flux measurement errors, which depend on the local exposure time and therefore the pointing on which the source was detected. Finally, in regions where several pointings overlap, the completeness must follow exactly the order and manner in which the two selection steps are applied to the source candidates.

Whereas the full information of the selection function and of the completeness of the sample is needed in cosmological studies, simpler methods are better suited to studies focused on the scaling relations. The CoMaLit approach can be applied to heterogeneous samples too. If we do not aim at expressing the observed distribution of clusters in terms of the underlying mass function, we have just to properly model the distribution to avoid Eddington/Malmquist biases. The distribution has to be flexible enough to fit the data.

Within this framework, the Malmquist bias and the probability to exceed a given count rate can be simply modeled in terms of the observed flux as a step function rather than as a more complicated position dependent probability of the true flux.

The correction for the Malmquist bias is relevant in the case of the X-ray luminosity as response variable. The flux limit is set within a 1′ aperture. We compute the luminosity threshold extrapolated out to 300 kpc at the cluster redshift, by assuming a $\beta$-profile with core radius $r_{\text{c}}=0.15\ r_{500}$ and $\beta=2/3$ (XXL Paper III). This procedure gives the luminosity threshold $L^{\mathrm{XXL}}_{\mathrm{300kpc,th}}$ for each cluster. According to the notation of Secs. 2 and 4, the threshold for the $j$-th cluster is

[TABLE]

The uncertainty $\delta_{y_{\text{th},j}}$ related to the extrapolation procedure is estimated by considering a range of radial profiles with scatter in slope of $\sigma_{\beta}\sim 0.1$, scatter in core radius of $\sigma(r_{\text{c}}/r_{500})\sim 0.1$, and correlation $\rho_{\beta r_{c}}\sim 0.66$ as representative of the sample of the 45 bright nearby galaxy clusters in Mohr et al. (1999). The median $\delta_{y_{\text{th}}}$ is $\sim 3\%$ but the distribution of values shows a long tail at larger values, so that the mean is $\sim 6\%$ and the standard deviation is $\sim 7\%$.

5 Theoretical predictions

The self-similar scenario of cluster formation and evolution was first proposed by Kaiser (1986) and later on extended and integrated (see e.g. Giodini et al., 2013). If gravity is the driving force of structure formation, X-ray quantities follow power-laws.

The relations for quantities within a fixed physical length differ from the canonical ones within the over-density radius. Here, the reference radius $R$ is constant and it does not scale with the mass. On the other hand, the density within $R$ is not constant and it changes with the mass. Since the density is not a constant multiple of the critical one, the time dependence factor $E_{z}^{2}$ connected to the critical density does not enter the relations.

Under the assumptions that clusters are closed boxes and baryons track the total mass,

[TABLE]

If the cluster is near hydrostatic equilibrium, then

[TABLE]

Finally, the X-ray emission in the soft band scales as

[TABLE]

with no appreciable temperature dependence (Ettori, 2015). By definition, luminosity can be written as

[TABLE]

By rearranging Eqs. (32–35), we obtain the self-similar scaling relations for properties measured within fixed radii,

[TABLE]

Reported slopes involving the luminosity are appropriate for the soft band (Ettori, 2015).

Baryonic processes can disrupt the self-similar relations. AGN (Active Galactic Nucleus) feedback or radiative cooling can remove cold, dense gas from the inner regions of low mass clusters, which makes the gas mass vs total mass and luminosity vs total mass relations steeper and the temperature-mass relation shallower (Truong et al., 2018).

Scaling relations can be scattered by a number of processes acting in different directions. Non-thermal sources of gas pressure, temperature inhomogeneity, substructures and clumps, unvirialized bulk motions, and subsonic turbulence play a role (Battaglia et al., 2012; Rasia et al., 2012).

Triaxiality is an additional source of scatter (Limousin et al., 2013; Sereno et al., 2013). Observed signals depend on the orientation of the cluster (Gavazzi, 2005; Oguri et al., 2005; Sereno, 2007; Sereno & Umetsu, 2011; Limousin et al., 2013; Sereno et al., 2013). For systems whose major axis points toward the observer, which are typically over-represented in signal-limited samples, X-ray luminosities and gas mass derived under the standard assumption of spherical symmetry are over-estimated. On the other hand, the majority of randomly oriented clusters are elongated in the plane of the sky and properties can be under-estimated.

A certain degree of correlation between intrinsic scatters is in place and has to be considered in multi-property galaxy cluster statistics to properly model the scaling relations (Evrard et al., 2014; Rozo et al., 2014; Maughan, 2014; Mantz et al., 2015). Correlations can come from internal structure, formation history, orientation, environment, and uncorrelated structure (Angulo et al., 2012).

The X-ray luminosity depends on the assembling history of the clusters (Mantz et al., 2016a). Massive mergers impact the luminosity-mass relation (Torri et al., 2004). The dynamical state of the cluster can cause a positive correlation between luminosity and temperature (Mantz et al., 2016b). Apart from transient shocks, luminosity and temperature can be depressed in merging clusters where energy in bulk motions has not yet virialized. On the other hand, dynamically relaxed, hot clusters show bright cores with higher than average luminosities and approximately average temperatures. In fact, perturbed or relaxed clusters move coherently along the luminosity-temperature relation (Rowley et al., 2004; Hartley et al., 2008).

Intra cluster-medium (ICM) processes impact correlations too. Radiative cooling reduce the amount of gas in the smallest systems, and at the same time their total luminosity (Truong et al., 2018).

Ettori (2015) showed how the normalizations of the scaling relations between the hydrostatic mass and the gas mass, the gas temperature, the X-ray bolometric luminosity, and the integrated Compton parameter depend upon the gas density clumpiness, the gas mass fraction, and the logarithmic slope of the thermal pressure profile. Scatter in the thermal pressure profile can cause positive correlation between the spectroscopic estimate of the temperature and the X-ray luminosity. Clumpiness induces positive correlation between the luminosity and the gas mass estimated under the hypothesis of a smooth distribution. Ettori (2015) also argued that deviations of the observed slopes from the self-similar expectations can be explained with a mass dependence of the gas mass fraction and of the logarithmic slope of the thermal pressure profile.

X-ray quantities in numerical simulations are positively correlated (Stanek et al., 2010; Truong et al., 2018). Stanek et al. (2010) performed a numerical study of the intrinsic covariance of cluster observables using the Millennium Gas Simulations. They adopted two different physical treatments: shock heating driven by gravity only (GO), or cooling and preheating (PH). The results in Stanek et al. (2010) depend on the adopted scheme. They found correlation factors at redshift zero between bolometric luminosity, spectroscopic-like temperature, and gas mass within $r_{500}$ of $\rho_{lt|m}=0.67$ (0.73), $\rho_{lm_{\text{g}}|m}=0.60$ (0.76), $\rho_{m_{\text{g}}t|m}=0.42$ (0.37) for the GO (PH) simulation.

6 Results

As reference analysis, we consider the statistical approach detailed in Section 2 in which X-ray observables are analyzed in pairs. Since the procedure is symmetric and both $X$ and $Y$ are affected by scatter, the results do not change whether we associate one observable to either $X$ or $Y$. The analysis was performed with the R-package LIRA111The package LIRA (LInear Regression in Astronomy) is publicly available from CRAN (Comprehensive R Archive Network) at https://cran.r-project.org/web/packages/lira/index.html.. As reference redshift, we consider $z_{\text{ref}}=0.01$. Results for the scaling parameters and the intrinsic scatters are summarized in Table 1.

As an alternative method, we fitted the three X-ray observables at once, see Section 3. The Bayesian hierarchical multi-response model was sampled with JAGS222JAGS (Just Another Gibbs Sampler) is a program for Bayesian data analysis publicly available at http://mcmc-jags.sourceforge.net/..

A critical aspect in Bayesian data analysis is the computational efficiency of the sampling of the posterior probability distribution. Gibbs sampling or other methods can efficiently constrain the distribution (Kelly, 2007; Mantz, 2016), but problems may arise for non standard distributions. This is the case of the truncated multinormal distributions presented in Sec. 3. To circumvent the problem, we neglected the truncation in Eq. (22) but we still considered the truncation in Eq. (23), which is much easier to sample since the uncertainty covariance matrix is known and fixed. We tested that this approximate scheme can still correct for the main effects of the Malmquist bias.

In our multi-response analysis, the gas mass acts as the $Y_{1}$ variable. As for the analysis in pairs, the choice of the $Y_{1}$ variable does not affect the measurement of the scatters. Results are summarized in Table 2. The agreement between the two alternative regression methods is substantial. Values listed in Tables 1 and 2 are in base 10 logarithm.

6.1 Scaling

Results of the reference regression analysis are reported in Table 1. Measured gas mass, temperature, and luminosity are aligned quite well, see Fig. 2. The marginalized 2D posteriori probabilities are showed in App. B.

The $l$-$t$ and $m_{\text{g}}$-$t$ relations are steeper than the self-similar predictions. This is confirmed by the multi-response analysis, see Table 2. There is no statistical evidence for time evolution, i.e. the $\gamma$ parameters are consistent with zero.

Parameter degeneracies can be analyzed with the 2D probability functions. As well known, slope and intercept are anti-correlated (see top-left panels in Figs. 8, 9, and 10).

Measured slopes are consistent with a prominent role of radiative cooling and AGN feedback in low mass systems, see Sec. 5. Stellar formation consumes the reservoirs of cold gas in the inner regions. Even though AGN feedback balances against over-cooling, it expels gas from the cluster core, which also reduces the luminosity and the total gas supply. As a result of these baryonic processes, which are most effective in low mass systems, the $l$-$t$ and $m_{\text{g}}$-$t$ relations are steeper than the self-similar predictions, in agreement with our results.

6.2 Intrinsic scatters

The gas mass is the less scattered proxy, $\langle\sigma_{m_{\text{g}}^{(Y)}|Z}\rangle=0.04\pm 0.01$ ($8.2\pm 3.0\%$), followed by temperature, $\langle\sigma_{t^{(Y)}|Z}\rangle=0.09\pm 0.02$ ($21.5\pm 4.2\%$) and luminosity, $\langle\sigma_{l^{(Y)}|Z}\rangle=0.19\pm 0.08$ ($43.6\pm 18.5\%$), see Fig. 3. In the reference analysis, the three X-ray quantities were analyzed in pairs, and we have two measurements of each scatter. The scatter values quoted above were obtained as the weighted mean and the standard deviation of the probability density obtained by combining the results from the two regressions summarized in Table 1.

The reference results (plotted as orange lines in Fig. 3) are fully consistent with the results from the multi-response analysis (plotted as red lines in Fig. 3), when we obtained $\langle\sigma_{m_{\text{g}}^{(Y)}|Z}\rangle=0.04\pm 0.02\leavevmode\nobreak\ (9.4\pm 5.2\%)$, $\langle\sigma_{t^{(Y)}|Z}\rangle=0.09\pm 0.02\leavevmode\nobreak\ (20.9\pm 4.4\%)$ and $\langle\sigma_{l^{(Y)}|Z}\rangle=0.07\pm 0.05\leavevmode\nobreak\ (15.6\pm 11.7\%)$. Due to skewness, the mean and standard deviation slightly differ from the biweight estimators quoted in Table 2.

The probability distribution of the luminosity intrinsic scatter extends over a significantly larger range than gas mass and temperature with a tail at the upper end, see Fig. 3. The standard deviation, i.e. the quoted uncertainty, is then larger than for $m_{\text{g}}$ and $t$.

Furthermore, the analysis of luminosity is more directly affected by Malmquist bias, which we treat as a threshold in the observed luminosity, see Sec. 4.1. Thresholds are expressed in terms of probability distributions, see Eqs. (4, 6, 24, 25), which affects the precision within which regression parameters are recovered.

The intrinsic scatter is best constrained when at least one low scatter proxy, e.g. the gas mass, is included in the fitting procedure. As seen from the comparison with the multivariate analysis in the middle and lower panels of Fig. 3, the $l$-$t$ fitting overestimates the larger scatter, i.e. $\sigma_{l^{(Y)}|Z}$, and underestimates the smaller one, i.e. $\sigma_{t^{(Y)}|Z}$. These two scatters are in fact anti-correlated, see Fig. 8. The posterior probability distributions are however compatible.

The temperature is intrinsically less scattered than the luminosity but is measured with a larger uncertainty. The two effects partially counter-balance and make the two proxies nearly as effective.

The time evolution of scatters is not well constrained. Parameter values are affected by large statistical uncertainties and are consistent with zero.

Our estimates of the scatter can be slightly overestimated since we compare deprojected quantities measured within the sphere, e.g. the gas mass, to projected quantities measured within the cylinder, e.g. the temperature and the luminosity. The dispersion associated to the non-universality of the density profiles slightly inflates the measured scatters.

In our analysis, we assume that the variable $X$ is an unbiased proxy, i.e. Eq. (1) reduces to $X_{Z}=Z$. The role of the $X$ was covered by either the temperature or the gas mass. The gas mass also worked as $Y_{1}$, i.e. the analog of $X$ in the multi-response analysis. These assumptions are in line with our expectations for quantities measured within a fixed radius, when $M_{\text{g}}\propto M$ and $T\propto M$, see Sec. 6.1. However, the determination of the conditional scatters and of the correlation factors is independent of the assumption on $X$. In fact, the estimates of scatters and slope are just weakly correlated, see Figs. 8, 9, and 10 and CoMaLit-IV. We further validate the stability of our results by considering tilted relations between $X$ and $Z$, e.g. $X_{Z}=3\times Z$, or interchanging the roles of $Y$ and $X$. Notwithstanding the different values of $\beta_{X|Z}$, the results for scatter covariance matrix are the same.

6.3 Scatter correlations

Our analysis suggests that the intrinsic scatters of gas mass, temperature, and luminosity are positively correlated, i.e. clusters of a given mass which are over-luminous have high temperature and an excess of gas mass. The probability distributions are peaked towards $\rho_{XY|Z}\sim 1$. Strong correlations are slightly preferred but the evidence is marginal due to the large statistical uncertainties. In fact, positive correlation is inferred at just the 1-$\sigma$ confidence level. Lower values of the correlation cannot be excluded, see the last rows of Figs. 8, 9, and 10. This is confirmed by the multi-response analysis, see Fig. 11.

The intrinsic correlation $\rho_{XY|Z}$ is partially degenerate with the slope of the relation $\beta_{Y|Z}$, see Figs. 8, 9, and 10.

Positive correlation between the properties of the ICM is expected as a results of formation history and dynamical state, see Sec. 5. Dynamically relaxed clusters are usually hotter and more luminous than merging systems of comparable total mass, where the ICM has yet to virialize (Rowley et al., 2004; Hartley et al., 2008; Mantz et al., 2016b). Positive correlation between the luminosity and the gas mass estimated under the hypothesis of a smooth distribution can be induced by clumpiness and the degree of regularity of the system (Ettori, 2015). Luminosity and gas mass can be overestimated in triaxial clusters elongated along the line of sight.

Our results support positive correlations even though the large statistical uncertainties ($\delta\rho\sim 0.5$, see Tale 2) prevent to distinguish between extreme scenarios ($\rho\leavevmode\nobreak\ \la 1$) and mild effects of co-evolution ($\rho\leavevmode\nobreak\ \ga 0$).

6.4 Distribution of the selected clusters

As expected for flux-limited selected samples, we find that the distributions of gas mass and temperature evolve with redshift, with more massive clusters preferentially included at high redshifts, see Table 3 and Fig. 4. These distributions are found as a result of the regression procedure, which does not exploit the knowledge of the mass completeness function obtained by numerical simulations, see Sec. 4.2, but recovers the distributions from the data.

The Gaussian function provides a good approximation to the redshift evolving distributions, see Figs.  5, 6, and 7. The temperature distributions derived from the analysis of the $l$-$t$ and of the $m_{\text{g}}$-$t$ are fully consistent.

As far as the tails of the covariate distributions are accounted for, conditional scatters and parameters of the scaling relations are well recovered (Kelly, 2007; Sereno, 2016). In fact, a simple normal distribution provides reliable results since it can account at the same time for the selection effects which penalizes low signal clusters and for the rarity of high mass clusters (Lima & Hu, 2005; Sereno & Ettori, 2015a). A covariate distribution following the halo mass function would fail to reproduce the suppression at the low mass end.

7 Previous results

Intrinsic scatter and correlation can be mass and time-dependent (Truong et al., 2018), and results from different samples should be compared cautiously. We have also to consider differences in measurements and definitions, mostly when comparing relations for either core-excised or core-included quantities. Furthermore, our sample extends to small groups whereas most of the previous analyses considered more massive clusters. Finally, the present level of statistical uncertainties is too large to make firm conclusions on apparent disagreements. However, the positive correlation found with the analysis of XXL-100-GC is consistent with previous results. Maughan (2014) applied the PICACS model to two X-ray samples of clusters with $T\ga 2\leavevmode\nobreak\ \text{keV}$ with measured core-excised temperatures, gas masses, and either hydrostatic masses or luminosities. Quantities were measured within $r_{500}$. The analysis suggested a positive correlation between the intrinsic scatter of $T$ and $M_{\text{g}}$ ($\rho_{m_{\text{g}}t|m}=0.31\pm 0.30$), and between $T$ and the core excluded bolometric luminosity ($\rho_{tl_{\text{ce}}|m}=0.37\pm 0.30$). A strong and significant correlation between the scatter in $M_{\text{g}}$ and $L_{\text{X}}$ was found ($\rho_{l_{\text{ce}}m_{\text{g}}|m}=0.85\pm 0.14$).

Mantz et al. (2015) constrained the cosmological parameters through the analysis of the mass function of a sample of X-ray selected massive clusters ($T\ga 4\leavevmode\nobreak\ \text{keV}$) detected in the ROSAT All-Sky Survey. They used follow-up measurements of soft band X-ray luminosity, temperature, and gas mass within $r_{500}$ and weak gravitational lensing measurements of a sub-sample of massive clusters. In the process, they assumed the gas mass scatter as being uncorrelated, i.e. they fixed the off-diagonal covariance terms involving $m_{\text{g}}$ to zero, and measured the correlation $\rho_{tl|m}=0.11\pm 0.19$. However, just using a larger amount of follow-up measurements and an updated calibration for X-ray observations, Mantz et al. (2016b) found positive strong correlation in the intrinsic scatters of luminosity and temperature at fixed mass ($\rho_{tl|m}=0.53\pm 0.10$).

Mantz et al. (2016a) studied the thermodynamic quantities of 40 massive clusters with $T\ga 5\leavevmode\nobreak\ \text{keV}$ identified as being dynamically relaxed and hot. Being the cluster relaxed, they identified the hydrostatic mass with the true mass and measured the off-diagonal terms of the scatter covariance matrix. They considered gas mass, core excised temperature, and core-excised or core-included soft-band [0.1-2.4] keV luminosity within $r_{500}$. They found $\rho_{tl|m}=-0.06\pm 0.24$ and $\rho_{tm_{\text{g}}|m}=-0.18\pm 0.28$, consistent with zero, and positive correlation between the core-included luminosity and the gas mass, $\rho_{lm_{\text{g}}|m}=0.43\pm 0.22$. The correlation is even stronger considering core-excised luminosity, $\rho_{l_{\text{ce}}m_{\text{g}}|m}=0.88\pm 0.06$.

Based on the analysis of 12 LoCuSS (Local Cluster Substructure Survey) massive clusters with $T\ga 5\leavevmode\nobreak\ \text{keV}$, Okabe et al. (2010) derived a 68.3 per cent confidence lower limit of $\rho_{tm_{\text{g}}|m}=0.185$, suggesting positive correlation between temperature and gas mass.

We caution that slope and normalization of the scaling relations studied here cannot be straightly compared to most literature results for two main reasons. First, we consider quantities within fixed physical radii rather than within $r_{500}$.

Secondly, we are interested in model variables and we consider scatter in both the $X$ and the $Y$ variables. The luminosity-temperature relation determined e.g. in XXL Paper III is the relation between the intrinsic luminosity and temperature and the related scatter measures the dispersions of the luminosities at a given temperature. The $l$-$t$ studied here is the hypothetical relation we would measure if temperature and luminosity were unscattered. This is the relation to study if we want to find the less scattered proxy since scatters are measured with respect to a basic third property. Similar considerations apply to the gas mass-temperature relation and XXL Paper XIII.

Notwithstanding the previous caveats, our results compare well with other recent studies. The $l$-$t$ and the $m_{\text{g}}$-$t$ relations for the XXL-100-GC sample were analyzed by XXL Paper III and XXL Paper XIII, which we refer to for detailed analysis and review of literature results. Taking into account selection effects, XXL Paper III found a bolometric luminosity-temperature relation steeper than the self-similar expectation. XXL Paper XIII found a gas mass-temperature relation steeper than the self-similar expectation, with a slope in agreement with Arnaud et al. (2007), who analyzed ten nearby relaxed clusters. Lovisari et al. (2015) analyzed XMM-Newton observations for a complete sample of local ($z<0.034$), flux-limited galaxy groups selected from the ROSAT All-Sky. They found a steeper than self-similar luminosity-temperature relation and a luminosity-gas mass relation compatible with expectations. Kettula et al. (2015) investigated groups from the XMM-CFHTLS survey together with high-mass systems from the Canadian Cluster Comparison Project and low-mass systems from the Cosmic Evolution Survey to find a steep luminosity-temperature relation.

8 Conclusions

Advancements in statistical methods applied to cluster physics have made it possible the detailed study of scaling relations. The CoMaLit approach exploited in the present paper shares important features with other Bayesian methods, e.g. Kelly (2007), Okabe et al. (2010), Rozo et al. (2010), Evrard et al. (2014), Andreon & Congdon (2014), Maughan (2014), and Mantz et al. (2015). The common theory behind these approaches comprises: the distinction of measured values, intrinsic scattered values, and model values; the modeling of the scaling relations as conditional probabilities; the modeling of the completeness function of the sample in terms of the intrinsic distributions of the underlying quantities. Some important differences can arise from: the treatment of the time-evolution of scaling and scatter; the treatment of uncertainties, scatters, and covariances; the modeling of the distribution of the covariate variable; the treatment of selection biases; the adopted priors.

Our analysis required neither the knowledge of the mass nor the assumption of hydrostatic equilibrium nor spherical symmetry. We only assumed that the matter halos that host clusters of galaxies have X-ray observable properties that are log-normally distributed around power-law scaling relations in halo mass.

We may think of the basic cluster property $Z$ as a ‘mass substitute’, i.e. the best statistical quantity to use when we have no access to total mass. It should be the quantity with the strongest correlation to the studied properties. In principle, this ideal quantity may be fictitious. Only direct measurements of the cluster masses can eliminate any ambiguity or doubt. However, theoretical considerations based on the self-similar model, numerical simulations, and analyses of cluster samples with measured mass support the simplest hypothesis that $Z$ is the mass.

We studied the scaling relations between X-ray properties of the XXL-100-GC sample and the intrinsic conditional scatters without any external mass calibration. In particular, we considered the spectroscopic temperature, the soft band luminosity, and the gas mass within fixed physical radii. The sample spans one order of magnitude in temperature from small groups at $T\sim\leavevmode\nobreak\ 0.6\leavevmode\nobreak\ \text{keV}$ to more massive clusters at $T\sim 7\leavevmode\nobreak\ \text{keV}$. This probes the lower end of the halo mass function at the group scale, whereas most of previous analyses focused on more massive clusters ($T\ga 4\leavevmode\nobreak\ \text{keV}$). The gas mass confirms itself as an excellent proxy. Even when measured within a fixed physical length, cluster properties can be recovered from gas masses with $\sim 8$ per cent accuracy.

Noteworthily, the gas mass is an equal or even better proxy to the weak lensing mass, which has an intrinsic scatter of $\sigma_{\text{WL}}\sim$ 15 per cent (CoMaLit-I; Mantz et al. 2015) and a better proxy than the hydrostatic mass, $\sigma_{\text{HE}}\sim$ 25 per cent (CoMaLit-I).

We considered only a subsample of X-ray properties. Other proposals as mass proxies can be appealing too. The integrated SZ Compton parameter $Y_{\text{SZ}}$ is expected to be tightly correlated to the energy content and the total mass of the clusters (Sereno & Ettori, 2015a), but its measurements can be elusive for small systems. The product of the temperature and $M_{\text{gas,500}}$, $Y_{\text{X}}$, is viewed as a robust mass indicator with low-scatter (Kravtsov et al., 2006). However, this proxy best performs if the temperature is core-excised, which is not practical in small groups. Furthermore, any positive correlation between intrinsic scatters of temperature and gas mass, as found in this paper, can worsen its performance.

Multi-probe analyses can open new windows on the evolution and formations of structures (Sereno et al., 2018). Halo properties can be better understood in terms of multi-dimensional plans than basic one-to-one relations (Fujita et al., 2018). Generalized scaling laws suitably weighting X-ray observables have to be considered to calibrate the proxy with the minimum scatter (Ettori, 2013).

We retrieved positive correlation between measured gas mass, temperature, and luminosity. The study of covariance between intrinsic scatters is important as it impacts the propagation of selection biases based on one observable to biases on other observable quantities (Maughan, 2014). For example, without taking the covariance between luminosity and gas mass or temperature into account, cluster masses estimated from $M_{\text{g}}$ or $T$ in an X-ray flux-limited sample would be biased high, with implications for cosmological studies (Maughan, 2014).

The simple assumption of underlying power-law relations is enough to estimate the intrinsic scatters of the observed properties and rank them. However, the loop cannot be closed without the mass information. This is needed to calibrate the scaling relation and confirm that the optimal cluster proxy is indeed the optimal mass proxy.

Acknowledgements

XXL is an international project based around an XMM Very Large Programme surveying two 25 deg2 extragalactic fields at a depth of $\sim 5\times 10^{-15}$erg cm*-2* s*-1* in the [0.5-2] keV band for point-like sources. The XXL website is http://irfu.cea.fr/xxl/. Multi-band information and spectroscopic follow-up of the X-ray sources are obtained through a number of survey programmes, summarised at http://xxlmultiwave.pbworks.com/.

MS and SE acknowledge financial contribution from the contract ASI-INAF n.2017-14-H.0. The Saclay group acknowledges long-term support from the Centre National d’Etudes Spatiales (CNES). This research has made use of NASA’s Astrophysics Data System (ADS) and of the NASA/IPAC Extragalactic Database (NED), which is operated by the Jet Propulsion Laboratory, California Institute of Technology, under contract with the National Aeronautics and Space Administration.

Appendix A Asymmetric errors

The likelihood function of the measured spectroscopic temperature is approximately Gaussian in log space (XXL Paper III). We used the method of Andreon (2012) to convert the asymmetric errors on $T$ computed by XSPEC to symmetric errors on $\log T$. If the probability distribution $p(T)$ is approximately log-normal, the standard deviation of the distribution $p(\log T)$ is given by

[TABLE]

where $T_{\text{mode}}$ is the mode and the uncertainties $\delta^{+}$ and $\delta^{-}$ are the points where the likelihood is lower than its maximum by a factor $\exp(-1/2)$.

Under the same assumption, we extended the prescription of Andreon (2012) and we also computed the mean as

[TABLE]

Appendix B 2D posteriors

The marginalized 2D posteriori probabilities for the luminosity-temperature, gas mass-temperature, and luminosity-gas mass relations are plotted in Figs. 8, 9, and 10. The results for the multi-response regression are shown in Fig. 11.