CHANDRA CLUSTER COSMOLOGY PROJECT III: COSMOLOGICAL PARAMETER CONSTRAINTS

Recent accelerated expansion of the universe detected in the Hubble diagram for distant type Ia supernovae (SNe) is one of the most significant discoveries of the past 10 years (Perlmutter et al. 1999; Riess et al. 1998). The acceleration can be attributed to the presence of a significant energy density component with negative pressure, hence the phenomenon is commonly referred to as dark energy. For a recent review of the dark energy discovery and related theoretical and observational issues, see Frieman et al. (2008) and references therein. Perhaps the simplest phenomenological model for dark energy is nonzero Einstein's cosmological constant. The SNe data indicated (and other cosmological data sets now generally agree) that a cosmological constant term currently dominates energy density in the universe.

The next big question is whether dark energy really is the cosmological constant. The properties of dark energy are commonly characterized by its equation-of-state parameter, w, defined as p = wρ, where ρ is the dark energy density and p is its pressure. A cosmological constant in the context of general relativity corresponds to a nonevolving w = −1. It is proposed that departures from the cosmological constant model should be sought in the form of observed w being either ≠−1, or evolving with redshift. Combination of SNe, cosmic microwave background (CMB), and baryonic acoustic oscillations (BAOs) data currently constrain |1 + w| < 0.15 at 95% CL (Komatsu et al. 2009). Observational signatures of such deviations of w from −1 are very small, and hence the measurements are prone to systematic errors. For example, variations of w between −1 and −0.9 change fluxes of z = 0.75 SNe in a flat universe with Ω_M = 0.25 by only 0.03 mag. Therefore, it is crucially important that the dark energy constraints at this level of accuracy are obtained from combination of several independent techniques. This not only reduces systematics but also improves statistical accuracy by breaking degeneracies in the cosmological parameter constraints.

One of the methods that has been little used so far is evolution in the number density of massive galaxy clusters. Evolution of the cluster mass function traces (with exponential magnification) growth of linear density perturbations. Growth of structure and distance–redshift relation are similarly sensitive to properties of dark energy, and also are mutually highly complementary methods (e.g., Linder & Jenkins 2003). Mapping between the linear power spectrum and cluster mass function relies on the model for nonlinear gravitational collapse. This model is now calibrated extensively by N-body simulations (see Section 3). The cluster mass function models also use additional assumptions (e.g., that the mass density is dominated by cold dark matter (CDM) in the recent past, and that the fluctuations have Gaussian distribution). However, corrections due to reasonable departures from these assumptions are negligible compared with statistical uncertainties in the current samples (we discuss these issues further in Section 3). It is important also that the theory of nonlinear collapse is insensitive to the background cosmology. For example, the same model accurately describes the relation between the linear power spectrum and cluster mass function in the Ω_M = 1, Ω_Λ = 0, low-density Ω_M = 0.3, Ω_Λ = 0, and "concordant" Ω_M = 0.3, Ω_Λ = 0.7 cosmologies (Jenkins et al. 2001).

Fitting cosmological models to the real cluster mass function measurements uses not only the growth of the structure but also the distance–redshift information because observed properties for objects of the same mass generally depend on the distance. Therefore, constraints on w derived from the cluster mass function internally make a combination of growth of structure and distance based cosmological tests, and thus potentially can be very accurate and competitive with any other technique (e.g., Albrecht et al. 2006).

Previous attempts to use the evolution of the cluster mass function as a cosmological probe were limited by small sample sizes and either poor proxies for the cluster mass (e.g., the total X-ray flux) or inaccurate measurements (e.g., temperatures with large uncertainties). Despite these limitations, reasonable constraints could still be derived on Ω_M (e.g., Borgani et al. 2001; Henry 2004). However, constraints on the dark energy equation-of-state from such studies are weak. For example, Henry (2004) derived the best-fit w = −0.42, only marginally inconsistent with w = −1, using the temperature function of the Einstein Medium Sensitivity Survey clusters; Mantz et al. (2008) determine w = −1.4 ± 0.55 with a larger sample of distant clusters (MACS survey; see Ebeling et al. 2001) but using the X-ray luminosity as a mass proxy.

The situation with the cluster mass function data has been dramatically improved in the past 2 years. A large sample of sufficiently massive clusters extending to z ∼ 0.9 has been derived from ROSAT PSPC pointed data covering 400 deg² (Burenin et al. 2007, Paper I hereafter). Distant clusters from the 400d sample were then observed with Chandra, providing high-quality X-ray data and much more accurate total mass indicators. Chandra coverage has also become available for a complete sample of low-z clusters originally derived from the ROSAT All-Sky Survey. Results from deep Chandra pointings to a number of low-z clusters have significantly improved our knowledge of the outer cluster regions and provided a much more reliable calibration of the M_tot versus proxy relations than what was possible before. On the theoretical side, improved numerical simulations resulted in better understanding of measurement biases in the X-ray data analysis (Nagai et al. 2007; Rasia et al. 2006; Jeltema et al. 2008). Even more importantly, results from these simulations have been used to suggest new, more reliable X-ray proxies for the total mass (Kravtsov et al. 2006). We discuss all this issues in the previous paper (Vikhlinin et al. 2009, Paper II hereafter). The cluster mass functions derived in this paper are reproduced in Figure 1. Overall, these results are an important step forward in providing observational foundation for cosmological work with the cluster mass functions.

Zoom In Zoom Out Reset image size

In this work, we present cosmological constraints from the data discussed in Paper II. The cosmological information contained in the cluster mass function data and relevant to dark energy constraints can be approximately separated into three quasi-independent components.

• 1.  
  Changes in the comoving number density at a fixed mass threshold constrain a combination of the perturbations growth factor and relative distances between low- and high-z samples; this by itself is a dark energy constraint (Section 8).
• 2.  
  The overall normalization of the observed mass function constrains the amplitude of linear density perturbations at z ≈ 0, usually expressed in terms of the σ₈ parameter. Statistical and systematic errors in the σ₈ measurement are now sufficiently small, and the ratio of σ₈ and the amplitude of the cosmic microwave background (CMB) fluctuations power spectrum gives the total growth of perturbations between z ≈ 1000 and z = 0—a second powerful dark energy constraint (Section 8.1).
• 3.  
  The slope of the mass function measures Ω_M × h; this by itself is not a dark energy probe but can be used to break degeneracies present in other methods.

Our dark energy constraints were derived for the following cases. Assuming a constant w and a flat universe, we measure w₀ = −1.14 ±  0.21 using only cluster data (i.e., evolution of the mass function between our two redshift samples) and the Hubble Space Telescope (HST) prior on h (Section 8.2). Combining cluster and Wilkinson Microwave Anisotropy Probe (WMAP) data, we obtain w₀ = −1.08 ± 0.15 but (w₀ is constrained much more tightly for a fixed Ω_M (Section 8.3). Finally, adding cluster data to the joint SNe + WMAP + BAO constraint, we obtain w₀ = −0.991 ± 0.045 (Section 8.3), significantly reducing statistical and especially systematic (Section 8.4) uncertainties compared with the case without clusters. A large fraction of the extra constraining power comes from contrasting σ₈ with normalization of the CMB power spectrum; this procedure is sensitive to nonzero masses of light neutrinos. Allowing for m_ν > 0, we obtain a new conservative upper limit ∑m_ν < 0.33 eV (95% CL) while still improving the w₀ measurement relative to the SN+WMAP+BAO-only case (w₀ = −1.02 ± 0.055, Section 8.5). Adding clusters also improves the equation-of-state constrains for evolving w in a flat universe (Section 9.1) and constant w in a nonflat universe (Section 9.2)

The paper is organized as follows. We start with a short summary of cluster data and systematic uncertainties (Section 2), discuss issues relevant for computing theoretical mass function models (Section 3), and describe our fitting procedure (Section 4). We then discuss constraints that can be obtained from low-redshift mass function only (Ω_Mh in Section 5 and σ₈ in Section 6). We then consider as an example constraints from the cluster evolution in nonflat ΛCDM model (i.e., w fixed at −1); Ω_Λ > 0 is required with ∼5σ confidence (Section 7). Constraints on the dark energy equation of state are considered in Sections 8–9. Systematic errors are discussed in Section 8.4.

This work is based on two cluster samples, originally compiled from ROSAT X-ray surveys (see Paper II for a complete description of the sample selection and data analysis). The low-redshift sample includes the 49 highest flux clusters detected in the All-Sky Survey at Galactic latitudes |b|>20° and z > 0.025. The effective redshift depth of this sample is z < 0.15. The high-redshift sample includes 37 z > 0.35 objects detected in the 400d survey, with an additional flux cut applied; the redshift depth of this sample is z ≈ 0.9. All the low- and high-z clusters were later observed with Chandra, providing good statistical precision spatially resolved spectral data thus yielding several high-quality M_tot estimators for each object. The combined cluster sample is a unique, uniformly observed data set. The volume coverage and effective mass limits at low and high redshifts are similar (see the estimated mass functions in Figure 1).

Because of the sufficiently high quality of the Chandra data, we employ advanced data analysis techniques going well beyond simple flux estimates and β-model fits commonly used in earlier studies. Cosmological cluster simulations have been used to test for the absence of significant observational biases in reconstructing the basic cluster parameters (Nagai et al. 2007). Using these simulations, we also tested which of the X-ray observables are best proxies for the total cluster mass (Kravtsov et al. 2006; Nagai et al. 2007) and concluded that the best three are the average temperature, T_X, measured in the annulus [0.15–1] r₅₀₀ (thus excluding the central region often affected by radiative cooling and sometimes, by active galactic nucleus (AGN) activity in the central galaxy); the intracluster gas mass integrated within r₅₀₀; and the combination of the two, Y_X = T_X × M_gas. These parameters are low-scatter proxies of the total mass (in particular, Y_X, and M_gas is only slightly worse). Simulations and available data show that the scaling of these proxies with M_tot, including the redshift dependence, is very close to predictions of the simple self-similar model. In a sense, even though we use advanced numerical simulations, which include multiple aspects of the cluster physics to test M_tot versus proxy relations, the role of simulations is limited to providing small corrections to predictions of very basic and hence reliable theory. Application of these corrections as well as practical considerations for deriving T_X, M_gas, and Y_X from the real data are discussed in Paper II.

Paper II also presents an observational calibration of the M_tot versus proxy relations using an extremely well observed sample of low-z clusters. This discussion is crucial for understanding the systematic uncertainties in our cluster mass function measurements, and we urge interested readers to consult Paper II. Table 4 in Paper II gives a summary of the main sources of systematic uncertainties in the derived cluster mass functions. They can be separated into three quasi-independent components. First is the uncertainty in calibration of the absolute cluster mass scale by Chandra hydrostatic mass estimates in a sample of dynamically relaxed, well observed low-z clusters (Vikhlinin et al. 2006); the level of this uncertainty (9%) is estimated from comparison of Chandra masses with two recent weak-lensing studies (Hoekstra 2007; Zhang et al. 2008). Second is uncertainties related to possible departures from standard evolution in M_tot–T_X, M_tot–M_gas, and M_tot–Y_X relations. This uncertainty (∼ 5%–6% between z = 0 and z = 0.5) was estimated from general reliability of numerical models of the cluster formation and from the magnitude of corrections that had to be applied to the data (see Section  4 in Paper II for details). The last major source of uncertainty is the evolution in the L_X–M_tot relation, affecting computations of the 400d survey volume coverage; this uncertainty is mostly measurement in nature because we derive the L_X–M_tot relation internally from the same cluster set. Its effect is negligible for the high-M end of the mass function and becomes comparable to Poisson errors for low-M clusters. A representative compilation of the effects of L_X–M_tot uncertainties on the V(M) function is presented in Figure 15 of Paper II.

The general reliability of our analysis is greatly enhanced by using independent, high-quality X-ray indicators of the total cluster mass—T_X, M_gas, Y_X. Since the masses estimated from these proxies depend differently on the distance to the object, the high-z mass functions estimated with different proxies should agree only if the assumed background cosmology is correct. In principle, this can be used as an additional source of information for the distance–redshift relation and folded into the cosmological fit. However, this method is nearly equivalent to the f_gas(z) test, which is more reliably carried out using direct hydrostatic mass estimates in relaxed clusters (Allen et al. 2008), and therefore we ignore this information. Instead, we use the agreement between different proxies observed for the best-fit cosmology as a comforting indication that there are no serious errors in our results.

In the current paradigm of structure formation, galaxy clusters form via gravitational collapse of matter around large peaks in the primordial density field (Kaiser 1984; Bardeen et al. 1986). Their abundance and spatial distribution in a comoving volume will thus depend on the statistical properties of the initial density field, such as Gaussianity⁹ and power spectrum (and hence the cosmological parameters that determine it), and could depend on the details of nonlinear amplification of the density perturbations by gravity. Indeed, semianalytical models based on the linear primordial density field and a simple ansatz describing nonlinear gravitational collapse of density peaks (Press & Schechter 1974; Bond et al. 1991; Lee & Shandarin 1998; Sheth et al. 2001) have proven to be quite successful in describing results of direct cosmological simulations of structure formation (e.g., Lee & Shandarin 1999; Sheth et al. 2001; Jenkins et al. 2001). The accuracy of the existing models, however, is limited and over the last several years the abundance of collapsed objects was calibrated by fitting appropriate fitting function to the results of direct cosmological simulations (Jenkins et al. 2001; Evrard et al. 2002; Warren et al. 2006). The fitting functions are expressed in the so-called universal form¹⁰ as a function of the variance of the density field on the mass scale M. The fact that such universal expressions exist implies that there is a direct link between the linearly evolving density field and cluster abundance.

In our analysis, we use the most recent accurate calibration of the halo mass function by Tinker et al. (2008), which provides fitting formulae for halo abundance as a function of mass, defined in spherical apertures enclosing overdensities similar to the mass we derive from observational proxies for the observed clusters. The Tinker et al. fitting formulae are formally accurate to better than 5% for the cosmologies close to the concordance ΛCDM cosmology and for the mass and redshift range of interest in our study; at this level, the theoretical uncertainties in the mass function do not contribute significantly to the systematic error budget. Although the formula has been calibrated using dissipationless N-body simulations (i.e., without effects of baryons), the expected effect of the internal redistribution of mass during baryon dissipation on the halo mass function are expected to be less than 5% (Rudd et al. 2008) for a realistic fraction of baryons that condenses to form galaxies.

Similarly to Jenkins et al. (2001) and Warren et al. (2006), the Tinker et al. formulae for the halo mass function are presented as a function of variance of the density field on a mass scale M. The variance, in turn, depends on the linear power spectrum of the cosmological model, P(k), which we calculate as a product of the initial power-law spectrum, kⁿ, and the transfer function for the given mixture of CDM and baryons, computed using the analytic approximations of Eisenstein & Hu (1999). This analytical approximation is accurate to better than 2% for a wide range of cosmologies, including cosmologies with nonnegligible neutrino contributions to the total matter density.

Our default analysis assumes that neutrinos have a negligibly small mass. The only component of our analysis that could be affected by this assumption is when we contrast the low-redshift value of σ₈ derived from clusters with the CMB power-spectrum normalization. This comparison uses evolution of purely CDM+baryons power spectra. The presence of light neutrinos affects the power spectrum at cluster scales; in terms of σ₈, the effect is roughly proportional to the total neutrino density, and is ≈20% for ∑m_ν = 0.5 eV (we calculate the effect of neutrinos using the transfer function model of Eisenstein & Hu 1999). Stringent upper limits on the neutrino mass were reported from comparison of the WMAP and Lyα forest data, ∑m_ν < 0.17 eV at 95% CL (Seljak et al. 2006). If neutrino masses are indeed this low, they would have no effect on our analysis. However, possible issues with modeling of the Lyα data have been noted in the literature (see, e.g., discussion in Section  4.2.8 of Dunkley et al. 2009) and so we experiment also with neutrino masses outside the Lyα forest bounds (Section 8.5).

We obtain parameter constraints using the likelihood function computed on a full grid of cosmological parameters affecting cluster observables (and also those for external data sets). The relevant parameters for the cluster data are those that affect the distance–redshift relation, as well as the growth and power spectrum of linear density perturbations: Ω_M, Ω_Λ, w (dark energy equation-of-state parameter), σ₈ (linear amplitude of density perturbations at the 8 h⁻¹ Mpc scale at z = 0), h, tilt of the primordial fluctuations power spectrum, and potentially, the nonzero rest masses of light neutrinos. This is computationally demanding and we describe our approach below.

The computation of the likelihood function for a single combination of parameters is relatively straightforward. Our procedure (described in Paper II) uses the full information contained in the data set, without any binning in mass or redshift, takes into account the scatter in the M_tot versus proxy relations and measurement errors, and so on. We should note, however, that since the measurement of the M_gas and Y_X proxies depends on the assumed distance to the cluster, the mass functions must be rederived for each new combination of the cosmological parameters that affect the distance–redshift relation—Ω_M, w, Ω_Λ, etc. Variations of h lead to trivial rescalings of the mass function and do not require recomputing the mass estimates. Computation of the survey volume uses a model for the evolving L_X–M_tot relation see Section  5 in Paper II, which is measured internally from the data and thus also depends on the assumed d(z) function. Therefore, we refit the L_X–M_tot relation for each new cosmology and recompute V(M). Sensitivity of the derived mass function to the background cosmology is illustrated in Figure 2. The entire procedure, although equivalent to full reanalysis of the Chandra and ROSAT data, can be organized very efficiently if one stores the derived ρ_g(r) and T(r) computed in some reference cosmology. It takes ≈20 s on a single CPU to re-estimate all masses, refit the L_X–M_tot relation, and recompute volumes for each new combination of the cosmological parameters.

Zoom In Zoom Out Reset image size

The next step is to compute, for each combination of Ω_M, Ω_Λ etc., the likelihood function on a grid of those parameters which do not affect the distance–redshift relation. In our case, these are σ₈, h, and when required, the power spectrum tilt or the neutrino mass. The cluster data are extremely sensitive to σ₈ and so we need a fine grid for this parameter. Fortunately, the mass function codes compute the mass functions for different values of σ₈ with other parameters fixed at almost no extra expense. The sensitivity of the cluster data to h and tilt is much weaker; therefore the likelihood can be computed on a coarse grid for these parameters and then interpolated.

With the acceleration strategies outlined above, it took us ∼9600 CPU hours (or 20 days using multiple workstations) to compute the cluster likelihood functions on full parameter grids for several generic models (nonflat ΛCDM, constant dark energy equation of state in a flat universe, constant w with a nonzero neutrino mass, linearly evolving w in a flat universe, constant w in a nonflat universe). Alternatively, simulating the Markov chains (Lewis & Bridle 2002) with sufficient statistics for all these cases would require approximately the same computing time.

After the cluster likelihood function was computed, we also computed χ² for external cosmological data sets—WMAP (5 year results), BAOs, and SNe Ia bolometric distances. Since we basically use analytic Gaussian priors for these data sets (see Section 8.1 below), these computations are fast and can be made on a fine parameter grid. We also use a Gaussian prior for the Hubble constant, h = 0.72 ± 0.08, based on the results from the HST Key Project (Freedman et al. 2001). This prior is important only when the constraints from the shape of the mass function (Section 5) come into play and when external cosmological data sets are not used in the constraints. When fitting the cluster data, we also keep the absolute baryon density fixed at the best-fit WMAP value, Ω_bh² = 0.0227 (Dunkley et al. 2009). This parameter slightly affects the calculation of the linear power spectrum (Eisenstein & Hu 1998). When we add the WMAP information to the total constraints, we marginalize the WMAP likelihood component over this parameter. If not stated otherwise, our cosmological fits also assume a primordial density fluctuation power spectrum with n = 0.95 (Spergel et al. 2007). Our results are completely insensitive to variations of n within the WMAP measurement uncertainties and even to setting n = 1.

Once the combined likelihood as a function of cosmological parameters is available, we use the quantity −2ln L, whose statistical properties are equivalent to the χ² distribution (Cash 1979), to find the best-fit parameters and confidence intervals.

In addition to statistical uncertainties, we also consider different sources of systematics. We do not include systematic errors in the likelihood function but instead refit parameters with the relations affected by systematics varied within the estimated 1σ uncertainties. This approach allows as not only to estimate how the confidence intervals are expanded from combination of all systematic errors, but also to track the most important source of uncertainty for each case. A full analysis of systematic errors is presented in Section 8.4 for the case of constraints on constant w in a flat universe; in other cases the systematic uncertainties contribute approximately the same fraction of the total error budget. We also verified that in the constant w case, our method of estimating the systematic errors produces the results which are very close to the more accurate procedure using the Markov chain analysis.

The shape of the cluster mass function reflects the shape of the linear power spectrum in the relevant range of scales, approximately 10 h⁻¹ Mpc in our case. This shape, for a reasonable range of parameters in the CDM cosmology is controlled (Bardeen et al. 1986) mostly by the quantity Ω_Mh. It is useful to consider constraints on this combination separately because they are nearly independent of the rest of the cosmological parameters we are trying to measure with the cluster data.

Fixing the primordial power-spectrum index to the WMAP value, n = 0.95, the fit to the local mass function¹¹ gives Ω_Mh = 0.184 ±  0.024 (purely statistical 68% CL uncertainties). The best-fit value is degenerate with the assumed primordial power-spectrum index, and the variation approximately follows the relation ΔΩ_Mh = −0.31Δn. The variations of n within the range constrained by the WMAP data, ±0.015, lead to negligibly small changes in our derived Ω_Mh.

An additional source of statistical uncertainty is that related to the derivation of the L–M relation, since we derive this relation from the same set of clusters. Uncertainties in the L–M relation are translated into those of the survey volume and hence the cluster mass function. Most of our cosmological constraints are primarily sensitive to the cluster number density near the median mass of the sample. This median mass, the V(M) uncertainties are small compared with statistics (see Section 6 in Paper II). The Ω_Mh determination, however, is based on the relative number density of clusters near the high and low mass ends of the sample. Since the volume is a fast-decreasing function at low M's, the V(M) variations are important. The most important parameter of the L–M relation in our case is the power-law slope, α (see Equation (20) in Paper II). Variations of α within the error bars (±0.14) of the best-fit value lead to changes in the derived Ω_Mh of ±0.027. Adding this in quadrature to the formal statistical errors quoted above, we obtain a total uncertainty of ±0.035 (see Table 1). We have verified that other sources of systematics in the Ω_Mh determination are much less important than those related to the L–M relation.

In principle, a nonzero mass of light neutrino has some effect on the perturbation power spectrum at low redshifts. We checked, however, that their effect on the shape of the cluster mass function is negligible for any ∑m_ν within the range allowed by the CMB data (Komatsu et al. 2009). Therefore, neutrinos do not affect our results on Ω_Mh.

Our determination of Ω_Mh = 0.184 ± 0.035 compares well with the previous measurements using cluster data and galaxy power spectra. Of the previous cluster results especially noteworthy is the work of Schuecker et al. (2003) whose constraints are based not only on the shape of the mass function but also on the clustering of low-z clusters. Their value is Ω_Mh = 0.239 ± 0.056 (errors dominated by uncertainties in the conversion of cluster X-ray luminosities into mass; this source of uncertainty is avoided in our work by using high-quality X-ray mass proxies). Ω_Mh is measured accurately also by galaxy redshift surveys. The results from the 2dF and Sloan Digital Sky Survey (SDSS) surveys are Ω_Mh = 0.178 ± 0.016 and 0.223 ± 0.023, respectively (Cole et al. 2005; Tegmark et al. 2004, we rescaled to n = 0.95 their best-fit values reported for n = 1). The individual error bars in galaxy survey results are smaller than those from the cluster data; however, a recent work by Percival et al. (2007c) suggests that the previous galaxy redshift results may be affected by scale-dependent biases on large scales. Indeed, there is a tension between the SDSS and 2dF values at ≃90% CL and the difference is comparable to the error bars of our measurement.

The cluster results can be improved in the future by extending the range of the mass function measurements. Not only can this improve statistical errors in the mass function measurements but can also improve the accuracy of the L–M relation, a significant source of uncertainty in our case. We note that it is more advantageous to increase statistics in the high-M range than to extend the mass function into the galaxy group regime. In addition to greater reliability of the X-ray mass estimates in the high-M systems, the surveys become dominated by cosmic variance approximately below the lower mass cut in our sample (the cosmic variance is estimated in Section 7.1 of Paper II using the prescription of Hu & Kravtsov 2003).

Combined with the HST prior on the Hubble constant, our constraint on Ω_Mh becomes a measurement for the matter density parameter, Ω_M = 0.255 ± 0.043 (stat) ±0.037 (sys), where systematic errors are also dominated by the slope of the L–M relation. This agrees within the errors with other independent determinations, such as a combination of BAO and CMB acoustic scales, Ω_M = 0.256 ± 0.027 (Percival et al. 2007b), and a combination of gas fraction measurements in massive clusters with the average baryon density from big bang nucleosynthesis, Ω_M = 0.28 ± 0.06 (Allen et al. 2008). It also agrees with another independent measurement based on our data, Ω_M = 0.30 ± 0.05 from the evolution of the cluster temperature function (see Section 7 below).

The normalization of the cluster mass function is exponentially sensitive to σ₈, the amplitude of linear perturbations at the length scale 8 h⁻¹ Mpc, approximately corresponding to the cluster mass scale (Frenk et al. 1990). Measuring this parameter with the cluster data has been a popular topic of research, especially using statistics of X-ray clusters (Frenk et al. 1990; Henry & Arnaud 1991; Lilje 1992; White et al. 1993, and many others thereafter). The strong sensitivity of the predicted cluster number density to σ₈ makes the determination of this parameter relatively insensitive to the details of the sample selection. Historically, different studies using very different cluster catalogs yielded similar results, if the data were analyzed uniformly. Determination of σ₈ is more sensitive to calibration of the absolute mass scale. For example, Pierpaoli et al. (2003) show that if M_tot for a fixed value of T_X is varied by a factor of 1.5, σ₈ derived from the local cluster temperature function is changed by Δσ₈ ≈ 0.13. Smaller biases are introduced if the effects of scatter in deriving the mass–luminosity relation are neglected resulting in incorrect computations of the survey volume (Stanek et al. 2006). Our present work includes advances in both of these areas and thus it is worth presenting an updated measurement of σ₈.

Determination of σ₈ from the cluster abundance data usually shows a strong degeneracy with the Ω_M parameter, typically, σ₈ ∝ Ω_M^−0.6 (e.g., Huterer & White 2002). The nature of this degeneracy is such that the mass function determines the rms amplitude of fluctuations at the given M_tot scale. The corresponding length scale is a function of Ω_M (M ∼ Ω_Ml³) and thus the derived σ₈ depends also on Ω_M and more weakly on the local slope of the linear power spectrum (see discussion in White et al. 1993). We need, therefore, to constrain σ₈ and Ω_M jointly. We used a grid of parameters of the flat¹² ΛCDM model (Ω_M, h, σ₈), and computed the cluster likelihood using the mass function for the local sample. We then add the Hubble constant prior (Section 4), and marginalized the combined likelihood over h.

The results are shown in Figure 3¹³. For a fixed Ω_M, the value of σ₈ is constrained to within ±0.012 (statistical). The degeneracy between σ₈ and Ω_M can be accurately described as σ₈ = 0.813(Ω_M/0.25)^−0.47. The Ω_M range along this line is constrained by the shape of the local mass function combined with the HST prior on the Hubble constant (Section 5). Including the high-redshift data, we obtain very similar results. For example, for Ω_M = 0.25, the total sample gives σ₈ = 0.803 ± 0.0105, to be compared with σ₈ = 0.813 ± 0.012 from low-z clusters only. This implies that the σ₈ measurement is dominated by the more accurate local cluster data, as expected.

Zoom In Zoom Out Reset image size

Systematic errors of the σ₈ measurement are dominated by the uncertainties in the absolute mass calibration. To test the effect of these uncertainties, we changed the normalization of the mass versus proxy relations by ±9% (our estimate of systematic errors in the mass scale calibration, see Section 2). The effect, shown by the dotted contour in Figure 3, is to shift the estimated values of σ₈ by ±0.02, just outside the statistical 68% CL uncertainties. This range can be considered as a systematic uncertainty in our σ₈ determination for a fixed Ω_M.

Our cluster constraints on σ₈ are more accurate (for a fixed Ω_M) than any other method, even including systematic errors (Figure 4). It is encouraging that our results are in very good agreement with recent results from other methods. The measurements based on lensing sheer surveys, cluster mass function with M_tot estimated from galaxy dynamics, and WMAP (5 year results assuming flat ΛCDM cosmology) are all within their respective 68% CL uncertainties from our best fit. This independently confirms that our calibration of the cluster mass scale is not strongly biased. Furthermore, the present systematic errors in the cluster analysis are smaller than the statistical accuracy provided by WMAP-5 and other methods. This allows us to effectively use the σ₈ information in the dark energy equation-of-state constraints (Section 8.3).

Zoom In Zoom Out Reset image size

We now move to models where the crucial role is played by the high-redshift cluster mass function data. The first case to consider is combined constraints for Ω_M and Ω_Λ in the nonflat ΛCDM cosmology. To better demonstrate what role the different components of the information provided by the cluster mass function play in the combined constraints, we consider two cases: (1) when the full cluster mass function information is used, and (2) when the shape information is artificially removed thus leaving only the evolutionary information.

In the first case, for each combination of parameters, we compute the full likelihood for the low- and high-z mass functions and add the HST prior on the Hubble constant (this is necessary for effective use of the mass function shape information, see Sections 4 and 5). We then marginalize the combined likelihood over nonessential parameters (σ₈ and h in this case) keeping the primordial power-spectrum index fixed at the WMAP best-fit value, n = 0.95. Removal of the shape information (our second case) is achieved by letting n vary and marginalizing over it. This is approximately equivalent to using a free shape parameter for the CDM power spectra, the approach often used in earlier cluster studies (e.g., Borgani et al. 2001). Constraints for both cases were obtained for mass functions estimated using all our three proxies, T_X, M_gas, and Y_X.

The results are presented in Figures 5 and 6. First, we can easily identify the role of using the mass function shape information (illustrated for the M_gas and Y_X proxies). Clearly, it mostly breaks the degeneracies along the Ω_M axis. The best-fit values and statistical uncertainties for Ω_M are very close to those derived from the shape of the local mass function (and nearly identical to those from the total sample, Section 5).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

For a fixed Ω_M, the observed evolution in the cluster mass function provides a constraint on Ω_Λ. Degeneracies in the Ω_M − Ω_Λ plane provided by different mass proxies applied to the same set of clusters differ because of the different distance dependences of the M_tot estimates via T_X, M_gas, and Y_X (see below). Even without the shape information, evolution in the Y_X and M_gas-based mass functions requires Ω_Λ > 0 at the 85% and 99.7% CL, respectively. Including the shape information, we obtain Ω_M = 0.28 ± 0.04, Ω_Λ = 0.78 ± 0.25 (and Ω_Λ > 0 is required at the 99% CL) from the Y_X-based analysis. The evolution of the M_gas-based mass function gives Ω_M = 0.27 ± 0.04, Ω_Λ = 0.83 ± 0.15, and Ω_Λ > 0 at 99.98% CL. The T_X-based mass function does not strongly constrain Ω_Λ but provides an independent measurement of Ω_M with almost no degeneracy with Ω_Λ: Ω_M = 0.34 ±  0.08, in good agreement with the mass function shape results (and also previous measurements based on the evolution of the cluster temperature function, see Henry 2004). In a flat ΛCDM model (the one with Ω_M + Ω_Λ = 1), the constraint is slightly tighter, Ω_M = 0.30 ± 0.05.

Systematic uncertainties of the Ω_Λ measurements are dominated by possible departures of evolution in the M_tot versus proxy relations. This issue is discussed in detail below in connection with the dark energy equation-of-state constraint (Section 8.4); here we note only that the systematic uncertainties are approximately 50% of the purely statistical error bars on the dark energy parameters (Ω_Λ, w). Therefore, our cluster data provide a clear independent confirmation for nonzero Ω_Λ.

Comments on the role of geometric information in the cluster mass function test. Cosmological constraints based on fitting the cluster mass function generally use not only information from growth of the structure but also that from the distance–redshift relation because derivation of the high-z mass functions from the data assumes the d(z) and E(z) functions. Quite generally, the estimated mass is a power-law function of these dependences, $\tilde{M} \propto d(z)^\beta \, E(z)^{-\varepsilon }$. Different mass proxies have different β and ε, and thus combine the geometric and growth of structure information in different ways and lead to different degeneracies in the derived cosmological parameters. We find that strongly distance-dependent proxies (such as M_gas, see Paper II) are intrinsically more powerful in constraining the dark energy parameters (Ω_Λ, w). By contrast, distance-independent proxies such as T_X result in poor sensitivity to dark energy but instead better constrain Ω_M. This is well illustrated by the results in Figure 5. The M_gas-based estimates for M_tot result (if we ignore the shape of the mass function information) in degeneracy approximately along the line Ω_M + Ω_Λ = 1. In fact, the evolution of the cluster mass functions derived from M_gas can be made broadly consistent with the Ω_M ≈ 1, Ω_Λ ≈ 0 cosmology if one allows for strong deviations from the CDM-type initial power spectra (Nuza & Blanchard 2006). However, the mass functions estimated from the temperatures of the same clusters are grossly inconsistent with such a cosmology, irrespective of the assumptions on the initial power spectrum (Ω_M = 1 is 8.3σ away from the best fit to the temperature-based mass function, Figure 5). It is encouraging that the 68% CL regions for all three mass proxies overlap near the "concordance" point at Ω_M = 0.25–0.3 and Ω_Λ = 0.7–0.75.

Next, we study constraints on a constant dark energy equation of state, w₀ ≡ p_X/ρ_X, in a spatially flat universe. The analysis using cluster data only is equivalent to the Ω_M − Ω_Λ case (Section 7). We compute the likelihood for the cluster mass functions on a grid of parameters: present dark energy density Ω_X (=1 − Ω_M), w₀, h, and σ₈, then add the HST prior on the Hubble constant (Section 4). Marginalization over nonessential parameters, h and σ₈, gives the likelihood as a function of Ω_M and w₀. We also obtain the equation-of-state constraints combining our cluster data with the three external cosmological data sets (following the reasoning of Dunkley et al. 2009, for the choice of these data sets).

8.1. External Cosmological Data Sets

SN Ia. We use the distance moduli estimated for the Type Ia SNe from the HST sample of Riess et al. (2007), SNLS survey (Astier et al. 2006), and ESSENCE survey (Wood-Vasey et al. 2007), combined with the nearby SN sample (we used a combination of all these samples compiled by Davis et al. 2007). Calculation of the SN Ia component of the likelihood function for the given cosmological model is standard and can be found in any of the above references.

Baryonic Acoustic Oscillations. Detection of the baryonic acoustic peak in the correlation function for large red galaxies in the SDSS survey leads to a good measurement of the combination

$$\begin{equation} \left[\frac{d_A(z)^2}{(cz)^2\,H(z)}\right]^{1/3}\!\!\sqrt{{\Omega _{\rm M}}H_0^2}\;\; \left[\frac{n}{0.98}\right]^{0.35} = 0.469 \pm 0.017 \end{equation} \tag{ 1 }$$

at z = 0.35 (Eisenstein et al. 2005, "SDSS LRG sample"). This prior mostly constrains Ω_M but has some sensitivity also to the dark energy equation of state.

A more recent measurement of the BAO peaks in the combined SDSS and 2dF survey data is presented in Percival et al. (2007a) who determine the BAO distance measure at two redshifts (z = 0.2 and z = 0.35) instead of one in Eisenstein et al. (2005). These new data are somewhat in tension (∼2σ) with the SN+WMAP results (see, e.g., Figure 11 in Percival et al. 2007a), which may artificially tighten the constraints when the BAO data are combined with SN Ia, WMAP, and clusters. We checked, however, that from the combination of SN Ia, WMAP, and SDSS-LRG BAO, we derive the parameter constraints that are essentially equivalent to those in Komatsu et al. (2009), who used the Percival et al. priors. Therefore, the choice of the BAO data set is unimportant in the combined constraints.

WMAP-5. The likelihood for WMAP 5 year data is computed using a simplified approach described in Section  5.4 of Komatsu et al. (2009). This involves a computation, for a given set of cosmological parameters, of three CMB parameters—angular scale of the first acoustic peak, ℓ_A; the so-called shift parameter, R; and the recombination redshift, z_*. The likelihood for the WMAP-5 data is then computed using the covariance matrix for ℓ_A, R, and z_* provided in Komatsu et al. This method is almost as accurate as direct computation of the WMAP likelihood (Wang & Mukherjee 2007) but is much faster, which allowed us to explore the entire multidimensional grid of the cosmological parameters instead of running Markov chain simulations. One additional note is that to compute the CMB likelihood, we had to add the absolute baryon density, Ω_bh², to our usual set of cosmological parameters and then marginalize over it. The reason is that while the average baryon density has very little impact on the rest of our analysis, the CMB data are very sensitive to Ω_bh², thus any variation of h must be accompanied by the corresponding variation of Ω_b without which the computation of the CMB likelihood would be inadequate.

The method outlined above recovers essentially the entire information from the location and relative amplitudes of the peaks in the CMB power spectrum (Wang & Mukherjee 2007). One additional piece of information is the absolute normalization of the CMB power spectra, reflecting the amplitude of density perturbations at the recombination redshift, z_* ≈ 1090. Contrasted with σ₈ determined from our cluster data at z ≈ 0, it constrains the total growth of density perturbations between the CMB epoch and the present, and thus is a powerful additional dark energy constraint.

WMAP-5 Plus Local σ₈. The WMAP team provides the amplitude of the curvature perturbations at the k = 0.02 Mpc⁻¹ scale,

$$\begin{equation} \Delta_{\frak R}^2 = (2.21\pm 0.09)\times 10^{-9}. \end{equation} \tag{ 2 }$$

Section 5.5 in Komatsu et al. (2009) gives the prescription of how to predict this observable for a given set of cosmological parameters and σ₈. A useful accurate fitting formula can also be found in Hu & Jain (2004):

$$\begin{eqnarray} \tilde{\Delta }_{\frak R} &\approx & \frac{\sigma _8}{1.79\times 10^{4}}\, \left(\frac{\Omega _b h^2}{0.024}\right)^{1/3} \left(\frac{{\Omega _{\rm M}}h^2}{0.14}\right)^{-0.563} \nonumber \\ && \times (7.808\,h)^{(1-n)/2}\, \left(\frac{h}{0.72}\right)^{-0.693} \frac{0.76}{G_0} \end{eqnarray} \tag{ 3 }$$

(we adjusted numerical coefficients to take into account that the Hu & Jain approximation uses the CMB amplitude at k = 0.05 Mpc⁻¹ while the WMAP-5 results are reported for k = 0.02 Mpc⁻¹). In this equation, G₀ is the perturbation growth factor between the CMB redshift and the present, normalized to the growth function in the matter-dominated universe: G(z) ≡ (1 + z) δ(z)/δ(z_CMB). This fitting formula helps to understand the nature of the σ₈ versus CMB amplitude constraint. The relation between σ₈ and $\Delta _{\frak R}$ depends on the absolute matter and baryon densities, Ω_Mh² and Ω_bh² (well measured by the CMB data alone), and on the total growth factor, G₀, and the absolute value of the Hubble constant, h. Both of these quantities provide powerful constraints on any parameterization of the dark energy equation of state (Hu 2005), and their combination does so as well.

Inclusion of this information in the total likelihood is straightforward. Given the usual set of cosmological parameters (Ω_X, w₀, h) plus σ₈, one computes

$$\begin{equation} \chi ^2_{\rm CMBnorm} = \big(\tilde{\Delta }_{\frak R}^2\times 10^9 - 2.21\big)^2/0.09^2, \end{equation} \tag{ 4 }$$

where $\tilde{\Delta }_{\frak R}$ can be obtained either from Equation (3) or as described in Komatsu et al. (2009). The χ²_CMBnorm component is then added to the cluster χ² and the sum marginalized over σ₈.

8.2. w₀ from Cluster Data Only

Constraints on the present dark energy density Ω_X and constant equation of state are presented in Figure 7. For comparison, we show separately the results derived only from the evolution of the M_gas and Y_X-based mass functions, and the effect of including the mass function shape information (Section 7 describes the procedure for removing shape information from the cluster likelihood function). We do not consider here the T_X-based mass estimates because they provide little sensitivity to the dark energy parameters. (Section 7). Just like in the Ω_M − Ω_Λ case, evolution of the M_gas and Y_X-based mass functions constrains different combinations of w₀ and Ω_X. The width of the confidence regions across the degeneracy direction is similar but the gas-based results are less inclined giving a little more sensitivity to w₀ for a fixed dark energy density—Δw₀ = ±0.17 from the M_gas-based functions and Δw₀ = ±0.26 from Y_X.

Zoom In Zoom Out Reset image size

Adding the mass function information combined with the HST prior on h breaks the degeneracy along the Ω_X direction. For example, the ellipse in Figure 7 shows the 68% CL region from fitting both the evolution and shape of the M_gas-based mass function. The one-parameter confidence intervals in this case are Ω_X = 0.75 ± 0.04 and w₀ = −1.14 ± 0.21. These results compare favorably with those from other individual methods—SNe, BAO, WMAP (Figure 8), although the SNe and CMB data provide tighter constraints on w₀ for a fixed Ω_X. The real strength of the cluster data is, however, when they are combined with the CMB and other cosmological data sets. The combined constraints are very similar for the M_gas and Y_X-based cluster mass functions, and therefore we discuss only the former hereafter.

Zoom In Zoom Out Reset image size

8.3. w₀ from the Combination of Clusters with Other Data

First, we consider a combination of the cluster data with the WMAP distance priors (see Section  5.4 in Komatsu et al. 2009). Cluster data bring information on growth of density perturbations and normalized distances in the z ≃ 0.0–0.9 interval, and—weakly—on the Ω_Mh parameter. Adding this information reduces the WMAP-only uncertainties on w₀ and Ω_X approximately by a factor of 2 (dark blue region in Figure 9): w₀ = −1.08 ± 0.15, Ω_X = 0.76 ± 0.04.

Zoom In Zoom Out Reset image size

A much more significant improvement of the constraints arises from the σ₈ determination from low-redshift clusters (dark red region in Figure 9). Comparison of the local determination of σ₈ with the CMB normalization mostly provides a measurement of the total perturbation growth factor between z_CMB and the present. This depends more sensitively on w₀ than the evolution of the cluster mass function because of, first, larger redshift leverage, and second, because the perturbation amplitude at high z is measured more accurately by CMB than by 37 clusters from the 400d survey.

Is it appropriate to use the σ₈ versus CMB normalization information in the dark energy constraints or does it require unreasonable interpolation of the dark energy parameterization to high redshifts? We note in this regard that for any combination of the cosmological parameters in the vicinity of the "concordance" model, w₀ ≃ −1, Ω_X = 0.25–0.3, the universe becomes matter dominated and enters the deceleration stage by z ∼ 1.5 − 2; the growth of perturbations is basically fixed after that at G(z) = 1. In other words, the CMB data can be used to safely predict the amplitude of density perturbations at z = 1.5–2 almost independently of the exact dark energy properties. As long as it is appropriate to use a particular dark energy parameterization in the z = 0–2 interval, it is therefore appropriate to use the same model for the joint clusters+WMAP fit.

By itself, adding the σ₈ information does not significantly improve the w₀ and Ω_X constraints (the total extent of the 1σ confidence regions is similar to the WMAP+evolution case), but the confidence region becomes much more degenerate with Ω_X (see the inner red region in Figure 9), which increases the potential for improvement when we combine these results with other cosmological data sets, BAO and SNe.

The combined constraints from all four cosmological data sets are shown in Figure 10 (inner dark red region). The 68% one-parameter confidence intervals are Ω_X = 0.740 ± 0.012 and w₀ = −0.991 ± 0.045. The importance of adding information from our cluster samples is illustrated by a factor of ∼1.5 reduction of the measurement uncertainties with respect to the WMAP+SN+BAO data alone: we obtain w₀ = −0.995 ± 0.067 without clusters (dark blue region in Figure 10; these results are essentially identical to those reported in Komatsu et al. 2009). Perhaps more importantly, including the cluster data also reduces systematic uncertainties by a similar amount (Section 8.4).

Zoom In Zoom Out Reset image size

The best-fit values of the Hubble constant and σ₈ from the combination of all data sets are h = 0.715 ± 0.012 and σ₈ = 0.786 ± 0.011. These values are within 68% confidence intervals of their determination by direct measurements (HST Key Project results for h and fitting the low-z cluster mass function for σ₈). The best-fit combination of the dark energy parameters is also within the 1σ confidence regions for each individual data set included in the constraints (Figure 10). Therefore, the best-fit cosmological model is a good fit to the data. In particular, Figure 17 from Vikhlinin et al. (2009) shows that the mass function models computed in the ΛCDM cosmology (w₀ = −1) provide a very good description of the data.

8.4. Systematic Uncertainties in the w₀ Measurements

We estimate the effect of known sources of systematics on the cosmological constraints by varying the corresponding individual sets of data or internal relations (e.g., evolution in L_X–M_tot entering the survey volume computations) within the estimated 1σ interval. We assume, optimistically, that the current WMAP and BAO data are free from significant systematics (i.e., that they are smaller than statistical uncertainties), and consider systematic errors only in the SN Ia and cluster data sets. In most cases, a single source clearly dominates the systematic error budget for a particular measurement, so we report on only those dominant sources.

The largest known source of systematic error in the SN Ia analysis is the correction for extinction in host galaxies and uncertainties in intrinsic colors of SN Ia (e.g., Frieman et al. 2008). As a measure of systematic uncertainty in the combined SN sample, we use ±0.13 in w₀ for fixed Ω_X, quoted by Wood-Vasey et al. (2007). We implement these errors by computing the SN likelihood in our experiments for (Ω_X, w₀ + 0.13) and (Ω_X, w₀ − 0.13) instead of (Ω_X, w₀).

8.4.1. Main Sources of Cluster Systematics

The largest sources of systematic errors in the cluster analysis are those in the normalization of the M_tot versus proxy relations. They can be separated into two almost independent components: (1) how accurately is the absolute cluster mass scale established by X-ray hydrostatic M_tot estimates in the low-redshift clusters, and (2) how accurately can we predict evolution in the M_tot versus proxy relations, i.e., the relative mass scale between low- and high-redshift clusters. The first component mainly affects the σ₈ measurements and the associated dark energy constraints, while the second component affects the results derived from using only evolution in the cluster mass function (those in Figure 7). Our estimates of the M_tot systematics are discussed extensively in Vikhlinin et al. (2009). For the absolute mass scale (M_tot for fixed Y_X, T_X, or M_gas) at z ≈ 0, we estimate ΔM_sys/M ≲ 9% mainly from the comparison of the X-ray and weak-lensing mass estimates in representative samples. This source of error is implemented by changing the normalization of the M_tot versus Y_X, M_gas, or T_X relations at z = 0 by ±9%. For uncertainties in the evolution of the M_tot versus proxy relations, we estimate ΔM/M ≈ 5% at z = 0.5, mainly from the comparison of the prediction of different models describing observed small deviations of the cluster scaling relations from self-similar predictions, and from the magnitude of these deviations and corresponding corrections we apply to the data. These uncertainties are implemented by multiplying the standard scaling relations by factors of (1 + z)^±0.12.

Comparable to the evolution in the M_tot versus proxy relation are measurement uncertainties in the evolution factor for the L_X–M_tot relation. We do not use L_X to estimate the cluster masses, but the relation is required to compute the survey volume for the high-z sample. The resulting volume uncertainty depends on the mass scale, and can become comparable to the Poisson error for the comoving cluster number density (see Section 5.1.3 in Vikhlinin et al. 2009). We tested how this influences the cosmological fit by varying the parameters of the L_X–M_tot relation within their measurement errors around the best fit (the evolution of L_X for fixed M_tot in our model is parameterized as E(z)^γ and γ is measured to ±0.33; see Section  5.1.3 in Vikhlinin et al. 2009).

Other sources of systematics in the cluster analysis (summarized in Vikhlinin et al. 2009) are negligible compared with those outlined above. We verified also that uncertainties in the intrinsic scatter in the M_tot-proxy relations are not important. The main reason is that in the dark energy constraints, we use high-quality mass proxies (Y_X and M_gas), which should provide mass estimates with small, 7%–10% scatter. Variations of this scatter by up to ±50% with respect to the nominal values do not significantly change the best-fit cosmological parameters. This conclusion is seemingly different from Lima & Hu (2005) because in that paper, they consider proxies with larger scatter (the effect on the cosmological parameter constraints is proportional to scatter squared), and also they assumed that the normalizations in the M_tot versus proxy relation are obtained from self-calibration while we use direct mass measurements for a well observed subsample.

The variations of the best-fit parameters due to the systematics discussed above are reported in Table 2 along with the dominant source of error for each combination of cosmological data sets. For example, variations in the evolution of the M_tot–M_gas and M_tot–Y_X relations affect the best fit to the cluster data only by Δw₀ = ±0.1, while statistical uncertainties are ±0.2 to ±0.3 for fixed Ω_X (Section 8.2); unless the systematics in this case are a factor of 2 larger than our estimates, they are unimportant.

Table 2. Parameter Constraints from the Combination of Clusters with Other Cosmological Data Sets

Parameter	Value	Data Set	Systematic Errors	Dominant Source of Systematics
Flat (Ωk = 0), constant w (w = w0)
w0	−1.14 ± 0.21	evol+shape+h, Section 8.2	±0.10, ±0.08,	Evolution of Mtot vs. proxy relations,
				evolution in LX–Mtot, respectively
w0	−1.08 ± 0.15	evol+CMB, Section 8.3	±0.025	Evolution of Mtot vs. proxy relations
w0	−0.97 ± 0.12	evol+CMB+σ8+BAO, Section 8.3	±0.038	Absolute cluster mass calibration
ΩX	0.732 ± 0.016	CMB+BAO+SN, Section 8.3
w0	−0.995 ± 0.067	CMB+BAO+SN, Section 8.3	±0.076	SN systematics
w0	$\bf -0.991\pm 0.045$	CMB+BAO+SN+evol+σ8, Section 8.3	$\bf \pm 0.022$, $\bf \pm 0.033$	SN systematics, cluster masses
ΩX	$\bf 0.740\pm 0.012$	CMB+BAO+SN+evol+σ8, Section 8.3
h	$\bf 0.715\pm 0.012$	CMB+BAO+SN+evol+σ8, Section 8.3
σ8	$\bf 0.786\pm 0.011$	CMB+BAO+SN+evol+σ8, Section 8.3
Flat (Ωk = 0), constant w (w = w0), nonzero neutrino mass
w0	−1.02 ± 0.055	σ8+CMBν+CMB+BAO+SN+evol, Section 8.5	±0.064	SN systematics.
∑mν	−11110.1 ± 0.12 eV,	σ8+CMBν+CMB+BAO+SN+evol, Section 8.5	±0.1 eV	SN systematics, cluster masses
	<0.33 eV (95% CL)
Flat (Ωk = 0), evolving w: w = w0 + wa(1 − a)
wa + 3.64(1 + w0)	1.0.05 ± 0.17	CMB+BAO+SN+evol+σ8, Section 9.1
Constant w (w = w0), nonflat (Ωk ≠ 0)
w0	−1.03 ± 0.06	CMB+BAO+SN+evol+σ8, Section 9.2

Notes. Codes used in column 3: evol, evolution of the cluster mass function; h, HST prior on Hubble constant; shape, shape of the cluster mass function; CMB, WMAP-5 distance priors; σ₈, comparison of the cluster-derived σ₈ with the CMB power-spectrum normalization (reflecting growth of perturbations between z_CMB and z = 0); BAO, BAO distance prior; SN, SN Ia luminosity distances; CMB_ν, WMAP-5+BAO+SN constraints on the neutrino mass (Section 8.5).

Download table as:  ASCIITypeset image

8.4.2. Systematics in the Combined Constraints

The most interesting case to consider is the reduction in the systematic errors from combining both SN and cluster data with the WMAP and BAO priors. In the SN+CMB+BAO case, the SNe systematics cause variations in the best-fit w₀ by ±0.076 (reduced from ±0.13 for the SN-only case mainly by including WMAP priors). Cluster systematics affects the w₀ constraints from the clusters+WMAP+BAO combination by ±0.04 (dominated by the ±9% uncertainties in the absolute mass scale). The influence of both sources of error is significantly reduced in the combined constraints. We find that the best fit w₀ from SN+clusters+WMAP+BAO is affected by ±0.022 by SN systematics, and by ±0.033 by cluster systematics. The total systematic error in the combined constraint is thus Δw₀ = ±0.04, almost a factor of 2 reduction from ±0.076 achievable without clusters.

We also note that if we significantly underestimate the cluster systematics, the most likely direction is that the cluster total masses are underestimated.¹⁴ If the cluster M_tot are revised high, this would lead to an increase in the derived σ₈, and decrease in w₀ when cluster data are combined with the CMB priors. Dark energy models predicting the equation-of-state parameter significantly above w₀ = −1 will be even less consistent with observations in this case.

8.4.3. Prospects for Further Reduction of Systematic Errors

It is reassuring that all sources of systematic errors we considered affect the dark energy equation-of-state constraints within the statistical measurement errors. This implies that while systematic errors are important, they do not yet dominate the current error budget. The situation will reverse in the future as the data sets expand. More effort will be needed then to reduce the systematics still further. We briefly outline the prospects for reducing the cluster-related systematics. Some of this will happen automatically as the high-z surveys become deeper and cover a larger area. For example, the V(M) uncertainties for our range of redshifts can be eliminated simply by decreasing the flux threshold by a factor of ∼4 compared with the 400d limit, making the sample volume-limited; such an extension will provide also a more accurate measurement of the L_X–M_tot relation. The absolute calibration of M_tot in low-z clusters can be improved by constraining sources of nonthermal pressure (e.g., if turbulence is of any importance for the M_tot estimates, it is easily detectable with an X-ray microcalorimeter), or through stacked weak-lensing analysis (e.g., measuring average lensing shear profiles for a large set of clusters with the same Y_X). To improve limits on nonstandard evolution in the M_tot versus proxy relations, we cannot use direct mass measurements of the high-z objects because they will be degenerate with the assumed distance–redshift relation. Instead, we should improve reliability of numerical models for cluster evolution. The biggest uncertainties in these models at present are related to the processes of gas cooling and star formation, and also to energy feedback from the central AGN. The strategy for future progress can be based on the fact that these processes most strongly affect cluster cores, which we do not use for the mass estimates. We can, therefore, use the data from the central regions to bracket a likely range of uncertainty in the model predictions for the cluster outer regions, where we derive the M_tot proxies. However, even with the current estimated uncertainties, the samples can grow by a factor of ∼4 before the systematics start to dominate. Ultimately, as the cluster surveys detect ∼10⁴ clusters with accurately measured X-ray parameters, the so-called self-calibration techniques (Majumdar & Mohr 2004; Lima & Hu 2004) can be employed to further constrain the evolution in the M_tot versus proxy relations.

8.5. Effects of Nonzero Neutrino Mass

If light neutrinos have masses in the range of a few 0.1 eV, they become nonrelativistic between z_CMB and z = 0, and this transition produces distortions in the matter perturbations power spectrum relative to the prediction of the pure CDM+baryons model. Using approximations of the transfer function from Eisenstein & Hu (1999), it is easy to verify that the effect is approximately proportional to the total mass of neutrinos (more exactly, to ∑m_ν/Ω_M), and the rms fluctuations at cluster scales today are suppressed by approximately 20% if ∑m_ν = 0.5 eV and Ω_M = 0.26. This effect is far outside the measurement uncertainties in σ₈ from clusters (we quote systematic errors of 3% from uncertainties in the M_tot calibration and statistical uncertainties are even smaller, see Table 1). Therefore, neutrino masses in this range (1) may affect the dark energy constraints when cluster data are combined with WMAP (because they will effectively change the relation between σ₈ and the CMB normalization, Equation (3)), and (2) can be tightly constrained by our cluster data.

To test the effect of neutrinos, we ran an additional set of models in which the total neutrino mass was allowed to vary between 0 and 1 eV. For simplicity, we assumed that there are three neutrino species with the same mass, but the final results are not very sensitive to this assumption. The only component of our procedure, which is significantly affected by the nonzero neutrino mass is contrasting the cluster-derived σ₈ with the WMAP normalization of the CMB power spectrum. We can no longer rely on Equation (3) and should instead use the full procedure described in Section  5.5 of Komatsu et al. (2009). Otherwise, the analysis is equivalent to the ∑m_ν = 0 case. The likelihood for all cosmological data sets was computed on our usual grid plus ∑m_ν as an additional free parameter, and then marginalized over Ω_X, h, and σ₈. Finally, we took into account that a combination of WMAP, BAO, and SN data provides some sensitivity to the neutrino mass through the so-called early integrated Sachs-Wolfe effect (see discussion in Section  6.1.3 of Komatsu et al. 2009, and references therein). From this analysis, Komatsu et al. derive a 95% upper limit of ∑m_ν < 0.66 eV. Since our procedure of using WMAP priors (Section 8.1) ignores this additional information, we included it approximately by adding a Gaussian prior ∑m_ν = 0 ± 0.33 eV to the final marginalized likelihood.

The derived constraints on ∑m_ν and w₀ are shown in Figure 11. As expected, when the σ₈ versus CMB normalization constraint is added, there is a degeneracy between the best-fit w₀ and the total neutrino mass. If we were using only clusters and WMAP, the degeneracy would approximately follow the line w₀ + 1 = −0.4∑m_ν and would extend to ∑m_ν ≈ 1.3 eV (the WMAP-only bound on the neutrino mass, Dunkley et al. 2009). This degeneracy is broken, however, when we add the BAO and SN information: low values of w₀ required by clusters+CMB for high values of the neutrino mass are inconsistent with these two data sets. Therefore, a combination of all four data sets can be used to constrain both w₀ and neutrino mass. The best-fit value is ∑m_ν = 0.10 ± 0.12 eV, with a 95% CL upper limit of ∑m_ν < 0.33 eV. This limit is significantly tighter than that is achievable without clusters (<0.66 eV at 95% CL). The constraint on w₀ degrades somewhat compared with the m_ν = 0 case: w₀ = −1.02 ± 0.055 (compared with ±0.045 for m_ν = 0), but is still better than ±0.067 without clusters (see Table 2). To conclude, adding the cluster information allows us to set tight limits on the neutrino mass while still improving the w₀ measurements with respect to the SN+WMAP+BAO case.

Zoom In Zoom Out Reset image size

Our constraints on the neutrino mass are still weaker than the published results from Lyα forest data, ∑m_ν < 0.17 eV (Seljak et al. 2006). Both the cluster and Lyα-based constraints use the same effect—suppression of the power spectrum at small scales by neutrinos—but they have completely different systematics. The main unknown in the Lyα analysis is the thermal state of the low-density intergalactic medium (IGM), usually estimated from numerical simulations; it has been suggested that the thermal state may be more complex than that assumed in previous work thus significantly weakening the m_ν bounds (Bolton et al. 2008). For clusters, the main uncertainty is the absolute mass calibration for low-redshift objects which affects the measurement of σ₈ (Section 6). The 9% systematic uncertainties on ΔM/M that we quote would translate into approximately ±0.075 eV for ∑m_ν, negligible compared with the current statistical uncertainties. We note that if the X-ray cluster mass measurements are wrong by more than 9%, it is almost certainly in the sense that they are underestimated (see footnote 14); the true value of σ₈ will then be higher than our measurement and the bound on the neutrino mass will be even tighter. Therefore, our 95% CL bound of ∑m_ν < 0.33 eV can be considered as a conservative upper limit.

Finally, we demonstrate how our cluster data improves parameter constraints for more general dark energy models. We consider two cases—evolving equation of state, w = w(z), and constant equation of state in a nonflat universe. The results are presented less completely than for the case of constant w in a flat universe. We also do not discuss systematic uncertainties separately for these cases; we checked that the importance of different sources of systematics and their fraction of statistical uncertainties is approximately the same as reported in Section 8.4 for the constant w, a flat universe case.

9.1. w(z) in a Flat Universe

We consider an often used parameterization of the equation-of-state evolution in which w changes linearly with the expansion factor, w(a) = w₀ + w_a(1 − a), or equivalently, w(z) = w₀ + w_az/(1 + z). We do not consider more complex parameterizations because constraints on the evolution term are still weak, and because neither parameterization has a clear physical motivation.

The likelihood function is computed on the Ω_M, w₀, w_a, h, σ₈ grid and then marginalized over Ω_M, h, and σ₈, leading to constraints in the w₀ − w_a plane shown in Figure 12. Constraints on w_a are weak with or without clusters. For example, the model with w₀ = −1.2 and w_a = 1 (leading to w = −0.7 by z = 1) is perfectly consistent with the data. However, clusters make the confidence region substantially narrower (improve w_a constraints for a fixed w₀). A cosmological constant model (w₀ = −1, w_a = 0) is still consistent with the data.

Zoom In Zoom Out Reset image size

Finally, we note that in either case, the degeneracy between w₀ and w_a is almost linear, w_a = A + B w₀. For such degeneracies, constraints on constant w are equivalent to those for evolving w at the pivot redshift, a_p = (1 + z_p)⁻¹ = 1 + 1/B (Hu & Jain 2004). From the slopes of degeneracies in Figure 12, we find z_p ≈ 0.29 without clusters and z_p ≈ 0.38 when cluster information is included. Therefore, our combined constraints on constant w (Section 8) can also be interpreted as those for evolving w at this pivot redshift.

9.2. w₀ in a Nonflat Universe

The final case we consider is constant w in a nonflat universe. The cosmological grid in this case is (Ω_k, Ω_M, w₀, h, σ₈) with the requirement that the dark energy density is Ω_X = 1 − Ω_M − Ω_k. The likelihood is marginalized over Ω_M, h, σ₈, and the constraints on Ω_k and w₀ are shown in Figure 13. Including clusters does not noticeably improve the measurement of Ω_k; by far the most significant contribution to the Ω_k constraint is from the combination of WMAP and BAO data (Komatsu et al. 2009). However, clusters do substantially improve the equation-of-state measurement: w₀ = −1.03 ±  0.06 to be compared with ±0.085 without clusters. A flat ΛCDM model (Ω_k = 0, w₀ = −1) is still consistent with the data within 68% CL.

Zoom In Zoom Out Reset image size

We presented constraints on the cosmological parameters from a new measurement on the galaxy cluster mass function in the redshift range z = 0–0.9. All major sources of information contained in the cluster mass function—its overall normalization and slope at z = 0, and evolution at high redshifts—are determined with our new data with a higher statistical accuracy and smaller systematic errors than before. This leads to much improved and more reliable constraints on the cosmological parameters.

From the normalization of the mass function estimated at low redshifts, we derive the σ₈ parameter degenerate with Ω_M: σ₈(Ω_M/0.25)^0.47 = 0.813 ± 0.013 (stat) ±0.024 (sys). The slope of the low-z mass function is a measure of Ω_Mh: Ω_Mh = 0.184 ± 0.037; combined with the HST prior on h, this is an independent measurement of Ω_M = 0.255 ± 0.043. The matter density can be independently measured with our cluster data using the evolution of the temperature function, yielding consistent results, Ω_M = 0.30 ± 0.05 in a flat ΛCDM model and 0.34 ± 0.08 in a general cosmology.

Evolution of the mass functions between z = 0 and 0.5 (median redshift for our high-z sample) constrains Ω_Λ = 0.83 ± 0.15 in nonflat ΛCDM cosmology, or the dark energy equation-of-state parameter, w₀ = −1.14 ± 0.21, in a spatially flat universe. Inclusion of the information provided by our cluster data also significantly improves the equation-of-state constraints obtained from the combination of multiple cosmological data sets. For example, by combining the 5 year WMAP, most recent SNe measurements, and detection of BAOs in the SDSS with our cluster data, we obtain w₀ = −0.991 ± 0.045 (stat) ±0.040 (sys); both the statistical and systematic errors in the combined constraint are a factor of 1.5–2 smaller than those without clusters. Including cluster information also improves results for an evolving equation-of-state parameter and for constant w in a nonflat universe. A spatially flat ΛCDM model is within the 68% CL interval from the best fit in all cases that we tested.

A good agreement between the geometric and growth of structure-based measurements of w in principle can be used to place limits on modified gravity theories which attempt to explain cosmic acceleration without dark energy (e.g., Wang et al. 2007). When self-consistent models of nonlinear collapse in such theories become available, it should be straightforward to use our cluster data in such tests also.

Comparison of the power-spectrum normalization at z = 0 obtained from clusters with the amplitude of the CMB fluctuations is a sensitive measure of the masses of light neutrinos. We constrain ∑m_ν < 0.33 eV at 95% CL, at the expense of slightly weakening the measurement of w₀ obtained assuming that the neutrino masses are negligibly small.

To facilitate the use of our cluster results in our cosmological studies, we provide at the project Web site¹⁵ machine-readable tables of the likelihood function computed on several cosmological grids.

We thank Jeremy Tinker for providing his mass function model results prior to publication. We also thank W. Hu, O. Gnedin, M. Markevitch, and H. Tananbaum for many useful discussions, and D. Spergel and A. Loeb for comments on the manuscript. Financial support was provided by NASA grants and contracts NAG5-9217, GO5-6120A, NAS8-39073 (A.V., W.R.F., C.J., S.S.M.), GO5-6120C (H.E.); NSF grants AST-0239759 and AST-0507666, NASA grant NAG5-13274 and the Kavli Institute for Cosmological Physics at the University of Chicago (A.K.); Sherman Fairchild Foundation (D.N.); FONDAP Centro de Astrofisica (H.Q.); Russian Foundation for Basic Research grants RFFI 05-02-16540 and RFFI 08-02-00974 and the RAS program OFN-17 (R.B. and A.V.).