Cluster counts. II. Tensions, massive neutrinos, and modified gravity

TL;DR

This paper investigates whether extensions to the standard cosmological model, such as massive neutrinos and modified gravity, can resolve tensions between CMB and galaxy cluster observations, finding that neutrinos do not help but modified gravity shows potential.

Contribution

The study introduces a simple modified gravity model with a non-standard growth index and analyzes its impact on cluster counts and tensions in cosmological data.

Findings

Massive neutrinos do not alleviate the cluster-CMB tension.

A modified gravity model can help break degeneracies in cluster data.

Calibration close to standard value suggests need for new physics.

Abstract

The $\Lambda$CDM concordance model is very successful at describing our Universe with high accuracy and few parameters. Despite its successes, a few tensions persist; most notably, the best-fit $\Lambda$CDM model, as derived from the Planck CMB data, largely overpredicts the abundance of SZ clusters when using their standard mass calibration. Whether this is a sign of an incorrect calibration or the need for new physics remains a matter of debate. Here we examined two simple extensions of the standard model and their ability to release this tension: massive neutrinos and a simple modified gravity model via a non-standard growth index $\gamma$. We used both the Planck CMB and SZ cluster counts as datasets, with or without local X-ray clusters. In the case of massive neutrinos, the SZ calibration $(1-b)$ is constrained to $0.59^{+0.03}_{-0.04}$ (68\%), more than 5$\sigma$ away from its standard value $\sim0.8$. We found little correlation between $\sum m_\nu$ and $(1-b)$, corroborating previous conclusions derived from X-ray clusters; massive neutrinos do not alleviate the cluster-CMB tension. With our simple $\gamma$ model, we found a large correlation between calibration and growth index but contrary to local X-ray clusters, SZ clusters are able to break the degeneracy between the two thanks to their extended $z$ range. The calibration $(1-b)$ was then constrained to $0.60^{+0.05}_{-0.07}$, leading to an interesting constraint on $\gamma=0.60\pm 0.13$. When both massive neutrinos and modified gravity were allowed, preferred values remained centred on standard $\Lambda$CDM values, but $(1-b)\sim0.8$ was allowed (though only at the $2\sigma$ level) provided $\sum m_\nu\sim0.34 $ eV and $\gamma\sim0.8$. We conclude that massive neutrinos do not relieve the cluster-CMB tension and that a calibration close to the standard value $0.8$ would call for new physics in the gravitational sector.

Equations

$$f(a)={\Omega }_m(a{)}^{\gamma },$$

$$T=A_{T-M}(h{M}_{\Delta }{)}^{2/3}{\left(\frac{{\Omega }_m\Delta (z)}{178}\right)}^{1/3}(1+z),$$

$$E^{-\beta }(z)\left[\frac{D_A^2(z){\bar{Y}}_{500}}{10^{-4}\mathrm{Mp}{\mathrm{c}}^2}\right]=Y_{\ast }{\left[\frac{h}{0.7}\right]}^{-2+\alpha }{\left[\frac{(1-b)M_{500}}{6\times 10^{14}{M}_{\odot }}\right]}^{\alpha },$$

Full text

1 Introduction

The Lambda cold dark matter ($\Lambda$CDM) model is remarkably successful at describing the majority of observations relevant to cosmology. These include probes of the background evolution of our Universe, its early times, and the dynamics and evolution of matter perturbations. Perhaps the most striking success of the $\Lambda$CDM model is its ability to explain the observed fluctuations of the cosmic microwave background (CMB), as measured most notably by the Planck satellite, which led to a determination of the cosmological parameters in the $\Lambda$CDM model with unprecedented accuracy.

However, some tension is also noted, particularly between the $\Lambda$CDM cosmological parameters favoured by Sunyaev-Zel’dovich (SZ) cluster counts, and those derived from the angular power spectra of CMB temperature and polarisation fluctuations (Planck Collaboration et al., 2014, 2016b). This particular tension appears only when the standard mass calibration $(1-b)=0.8$ of SZ clusters (defined as the ratio of the hydrostatic to true mass) is used, leading to a significant discrepancy in the $\sigma_{8}$ (the present-day amplitude of matter fluctuations) versus $\Omega_{m}$ (the dark matter density) plane, or equivalently a discrepancy on the calibration parameter when left free, compared to its standard value (see Fig. 3 in Sakr et al., 2018, and Fig. 1 in the present article).

A critical aspect of this tension is that the mass of a cluster (including its dark matter halo) is not directly observable. The above tension is therefore present only when additional priors on cluster masses are used, and it more generally relies on the use of scaling relations between halo mass, redshift, and the observable cluster properties (known as mass proxies). Based on the assumption of hydrostatic equilibrium, the Planck Collaboration used X-ray observations from XMM-Newton in their analysis to derive masses of the intra-cluster medium (Arnaud et al., 2010) and corrected the value by 20% to account for the incomplete thermalisation of the gas. This correction is expressed as a ‘bias’ $b$ in mass, defined as $(1-b)=0.8$ in the fiducial case.

Should the aforementioned tension be confirmed, we would be forced to consider extensions or alternatives to the standard $\Lambda$CDM model of cosmology. Since the CMB itself is mostly sensitive to the physics of the early Universe, one can reconcile it with cluster observations by introducing a modification to the cosmological model that has a significant impact preferentially at late times. More specifically, a new theory leading to a lower growth rate of structures would be required in order to predict a lower abundance of clusters, which would better reflect the data.

Such growth can be achieved in several ways; one possibility is to add mass to neutrinos in the standard cosmology (beyond their minimal experimental mass, 0.06 eV). Among other effects, massive neutrinos slow down the growth of matter perturbations during the matter and dark energy-dominated eras on scales smaller than their free-streaming length (see Lesgourgues & Pastor, 2012, for a review on the effect of neutrinos in cosmology). Alternatively, one could consider replacing the standard cosmological constant with a different form of dark energy. A fluid with a varying equation of state (the so-called quintessence models) or the addition of a scalar field can achieve the desired result, but a more radical approach consists of modifying the current theory of general relativity itself (see Clifton et al., 2012; Ishak, 2019, for extensive reviews).

In Sakr et al. (2018), we pointed out a similar discrepancy between the abundance of local X-ray clusters and CMB cosmology. We found that massive neutrinos do not solve the issue, while a simple modified gravity model parametrised by a modified growth rate, the $\gamma$ model, can reconcile the Planck cluster mass calibration with the Planck CMB cosmology at the expense of a relatively high $\gamma\sim 0.9$ (a value potentially already excluded by current data; see L’Huillier et al., 2018; Zarrouk et al., 2018). In our study, we combined constraints from the redshift-mass distribution of SZ clusters with the latest CMB measurements from the Planck satellite, performing a Bayesian analysis through a Markov chain Monte Carlo (MCMC) analysis. We varied the parameters of the standard model, the cluster mass calibration $(1-b)$, and the additional parameters introduced in the context of the extended models we considered. They included the total mass of three massive degenerate neutrinos and the phenomenological $\gamma$ parameter of our simple modified gravity model.

In Sect. 2, we describe the formalism used for predicting cluster abundances, as well as the extensions to the standard model we examined. In the same section, we also detail the datasets used in this work, and the implementation of the MCMC analysis which samples the posterior probability distribution function. We present and discuss our results in Sect. 3, while we summarise and discuss our conclusions in Sect. 4.

2 Methodology and data

2.1 Cluster abundances as a cosmological probe

Clusters of galaxies provide tight constraints on cosmological models, as they (mostly) probe the growth of matter fluctuations directly throughout the history of the Universe. Their use as cosmological probes has been extensively studied and exploited in the past. The methodology that we adopted here is closely related to the approach followed by Ilić et al. (2015) and Sakr et al. (2018). We refer the interested reader to those two articles for more details, as we only provide the main steps here.

We computed the expected numbers of clusters (as a function of mass and redshift) in a given cosmological model via integrals over the mass function, with the latter being taken from Despali et al. (2016) (thereafter D16). This more recent mass function was chosen over available alternatives (notably the popular Tinker et al., 2008, mass function), although this choice has only a minimal impact on our results and conclusions, as can be seen in Fig. 1 when comparing the filled (D16) and dash-dotted (Tinker) purple contours. The D16 mass function was optimised to work with the virial radius as a definition for clusters, which is the convention we used for X-ray clusters, although this choice had little impact on observables when modelling was done self-consistently (see Sakr et al., 2018). We did, however, use the ‘critical $M_{500}$’ convention (i.e. the density contrast of clusters is defined as 500 times the critical density of the Universe) for SZ clusters, as the available Planck SZ dataset adopts this particular convention. We used the analytical formulae provided by the D16 authors in order to convert their mass function between their fiducial (virial) cluster mass definition and any other definition (critical $M_{500}$ here).

In the case of massive neutrinos, we used the so-called CDM prescription to amend the mass function (Costanzi et al., 2013; Castorina et al., 2014) unless otherwise stated. The filled and dashed orange contours of Fig. 1 illustrate the effect of the CDM prescription in the analysis; as discussed in Sakr et al. (2018), we include its effects for the sake of completeness, even though it led to only marginal differences. As seen in Fig. 1, it led to slightly lower values of $\sigma_{8},$ and concurrently high values of $\Omega_{m},$ but with slightly increased error bars on both.

When considering a simple modified gravity model, as originally proposed by Linder (2005), we modelled the linear growth rate of structures like so

[TABLE]

where $\gamma$, the growth index, is here a free parameter generalising the approximation $\gamma\sim 0.545$ in standard gravity (Peebles, 1980). The growth rate $f$ is defined as $f(a)\equiv d\ln G(a)/d\ln a$, where $G$ is the (scale-independent) linear growth factor of matter density perturbations $\delta_{m}$ defined as $\delta_{m}(a)=\delta_{m,0}G(a)/G(0)$.

In all extensions to the standard model considered here, we assume the non-linear regime would follow the same threshold $\delta_{c},$ and the same virial density, $\Delta_{v},$ as in the standard picture. Scaling laws are then used to relate theoretical cluster quantities (i.e. mass) to observable quantities. For the X-ray temperature $T$ versus cluster mass $M$, we use the standard scaling relation

[TABLE]

where $A_{T-M}$ is the normalisation parameter, $\Delta$ is the density contrast chosen for the definition of a cluster, expressed with respect to the total background matter density111Eq. (2) is different when the contrast density $\Delta_{c}$ is expressed with respect to the critical density at redshift $z$. of the Universe at redshift $z$, and $M_{\Delta}$ is the mass of the cluster according to the same definition (as mentioned before, in our case for X-ray clusters, we used the virial radius and mass definition). Similarly, we used the following scaling relation between the expected mean SZ signal (denoted by the letter $\bar{Y}$) and mass $M_{500}$

[TABLE]

where the quantity $D_{A}$ refers to the angular diameter distance and $E(z)\equiv H(z)/H_{0}$. We note here the aforementioned hydrostatic mass bias $(1-b)$, while $\alpha$, $\beta$, and $Y_{\ast}$ are additional free parameters in the SZ scaling law. A fourth free parameter exists, namely $\sigma_{\ln Y}$; it corresponds to the scatter of the log-normal distribution that the measured $Y_{500}$ is assumed to follow (and of which the mean $\bar{Y}_{500}$ is given by Eq. (3))222Planck Collaboration et al. (2016b) provides best-fit values and uncertainties for all four of those parameters, derived from a fit based on X-ray observations. In our analysis, we follow the conventions used by the Planck Collaboration by fixing the $\beta$ parameter to its best-fit value ($=2/3$), letting $\sigma_{\ln Y}$ free with Gaussian priors derived from its fit ($=0.075\pm 0.01$), and leaving $\alpha$ completely free. We do, however, fix $Y_{\ast}$ to its best-fit value ($\log Y_{\ast}=-0.186$) instead of leaving it free, noting that it is degenerate with the bias parameter $(1-b)$..

2.2 Datasets and numerical tools

We followed a standard MCMC analysis using two publicly available codes: the CosmoMC package (Lewis & Bridle, 2002; Lewis, 2013) and the Monte Python code (Audren et al., 2013) that respectively use the CAMB and CLASS Boltzmann codes for the computation of relevant cosmological quantities: the matter power spectrum and the power spectra of CMB anisotropies. Those two codes were used to explore our full-parameter space under the constraint of our data sets extracted from CMB measurements (Planck Collaboration et al., 2016c), baryon acoustic oscillations (BAO) data (Anderson et al., 2014; Font-Ribera et al., 2014), X-ray-selected clusters (Ilić et al., 2015) and an SZ clusters sample (Planck Collaboration et al., 2016b), separated or combined.

The 2015 Planck CMB (Planck Collaboration et al., 2016a) and BAO (Anderson et al., 2014; Font-Ribera et al., 2014) likelihoods were already interfaced with the two aforementioned MCMC codes. The SZ likelihood (Planck Collaboration et al., 2014) is present in the CosmoMC package, while we implemented it ourselves in the Monte Python code. We also used our dedicated module for computing the likelihood associated with the X-ray cluster sample (Ilić et al., 2015). While massive neutrinos are already implemented in both the CLASS and CAMB codes, we implemented there the other $\Lambda$CDM extension we considered, namely, the $\gamma$ model. Additionally, we modified our X-ray and SZ clusters likelihood modules accordingly to support the massive neutrinos case and $\gamma$ parametrisation, all according to the recipes described earlier.

2D confidence regions and 1D constraints on parameters were then derived (as is usual in a Bayesian formalism) from the posterior distribution resulting from the MCMC chains. All parameters limits mentioned in the text are to be read as 68% confidence limits.

3 Results

The long-standing, so-called discrepancy between SZ cluster counts and the Planck CMB measurements in the $\Lambda$CDM model is illustrated in Fig. 1 in the $\Omega_{m}-\sigma_{8}$ plane, with contours from the two probes: in blue for the 2015 Planck CMB data, and purple for SZ cluster counts (the calibration parameter $(1-b)$ being fixed to the Planck Collaboration fiducial value of $0.8$). One can see that the separation between the two contours is significantly large and remains essentially identical when considering the new Planck 2018 data release (green contours). Allowing massive neutrinos with varying mass has the effect of pushing both the CMB contours (light blue and light green contours respectively for the 2015 and 2018 Planck CMB data) and the SZ cluster contours (orange) to higher values of $\Omega_{m}$ and lower values of $\sigma_{8}$, elongating the two in the direction of the $\Omega_{m}-\sigma_{8}$ degeneracy, with the CDM prescription shifting the contours slightly to higher $\Omega_{m}$ and lower $\sigma_{8}$.

The SZ clusters-CMB tension is therefore not relieved by the possible presence of massive neutrinos, leading to essentially identical conclusions in the case where only a sample of local X-ray clusters are used (Sakr et al., 2018). It should be noted that in Fig. 1 we added constraints from Big Bang nucleosynthesis (BBN) and BAO data to restrict $\Omega_{b}$ and $H_{0}$ to realistic values for cluster contours. We checked that considering wider priors results in more elongated contours in the direction of the $\Omega_{m}-\sigma_{8}$ degeneracy, but without relieving the tension.

3.1 The $\mathit{\Lambda}$CDM case

In order to examine this tension in a quantitative way, we combine CMB data and the SZ cluster sample, adding or not adding the constraints obtained with the X-ray cluster sample, and freeing the SZ and X-ray calibration factors, as well as the slope of the mass-SZ relation and its dispersion. In the standard $\Lambda$CDM model (as defined by the Planck Collaboration, which includes one massive neutrino with a minimal mass of 0.06 eV) one can see from Fig. 2 that, in the cases of CMB data with SZ and/or X-ray clusters, the constraints on all cosmological parameters are essentially identical. We obtain $(1-b)=0.58\pm 0.03$ and note a clear correlation between the $(1-b)$ calibration factor for SZ clusters and the $A_{T-M}$ calibration for X-ray clusters (at virial mass), indicating the consistency of the constraints obtained from the two types of cluster samples.

3.2 Are massive neutrinos a solution?

The possible role of massive neutrinos in alleviating the cluster-CMB tension is discussed in past literature, with somewhat different conclusions, (Rozo et al., 2013; Dvorkin et al., 2014; Mantz et al., 2015; Roncarelli et al., 2015; Bull et al., 2016; Salvati et al., 2018) often without clear indication as to whether conclusions depend on a calibration choice. In Fig. 3, we present the constraints on the SZ cluster calibration parameter $(1-b)$ from the combination of CMB data sets with the SZ cluster sample in the presence of massive neutrinos (blue contours). We obtain a constraint on the calibration $(1-b)=0.585^{+0.031}_{-0.037}$ essentially identical to the minimal neutrino case. Furthermore, no strong correlation appears in the calibration-neutrino masses plane, similarly to what is obtained from the X-ray sample alone (Sakr et al., 2018). We also notice that the constraint on neutrino masses is quite noticeably improved by combining clusters with the CMB (with 68% limits being reduced by about $\sim 28\%$) compared to the CMB-only constraint (black curve). The standard Planck mass calibration ($(1-b)\sim 0.8$) is formally rejected at more than 5 $\sigma$. However, for illustration purposes, we show the consequence of putting a strong (Gaussian) prior on $(1-b)$ consistent with the standard Planck calibration, namely $(1-b)=0.8\pm 0.01$ (green contours). In this case, the favoured neutrino masses are found to be non-zero: $\sum m_{\nu}=0.70\pm 0.15$ eV, but the corresponding best-fit model has a $\Delta\chi^{2}\approx 20$ compared to the case without a prior on $(1-b)$. This result is consistent with previous investigations obtaining non-zero neutrino masses using such a prior. It is therefore important to realise that Planck CMB and Planck SZ cluster counts are perfectly consistent with one another, but as mentioned before, this combination of data requires an SZ mass calibration of $(1-b)\sim 0.6$, significantly lower than the standard calibration $(1-b)\sim 0.8$.

3.3 Introducing modified gravity

If the standard calibration $(1-b)\sim 0.8$ is consolidated by observations, the conclusion of the previous sections leads us to examine alternative theories to $\Lambda$CDM. One option is to consider a variable dark energy equation of state $w$ instead of a simple cosmological constant. However, Salvati et al. (2018) found that a constant $w$ does not lead to a significant change in the required calibration. Any proposed alternative explanation should lead somehow to a lower amplitude of matter fluctuations at low redshift without altering the remarkable agreement with the CMB seen in the standard cosmological picture. Such a lower $\sigma_{8}$ would naturally come from a late-time deviation of the growth of matter fluctuations with respect to $\Lambda$CDM, a possibility that could naturally arise from a modification of the theory of general relativity for gravity. In the following, we will consider a phenomenological approach through the introduction of the growth index $\gamma$ following Eq. (1). We follow the same approach as in previous sections and perform an MCMC analysis to explore the space of the standard cosmological parameters, adding the calibration $(1-b)$ and the index $\gamma$ as additional free parameters under the constraints of our data sets. As previously found in Sakr et al. (2018), we expect increasing $\gamma$ to yield a lower $\sigma_{8},$ and thereby a higher calibration $(1-b)$ - corresponding to clusters being less massive. Similarly to the previous section, we use the combination of CMB and SZ cluster data, with or without adding our X-ray cluster sample to constrain the parameters of the model.

Our results are summarised in Fig. 4; the SZ calibration parameter $(1-b)$ shows a clear correlation with $\gamma$, as anticipated. However, a significant difference appears compared to the case where only the X-ray sample is used: as the SZ sample spans a wide range of (comparatively) higher redshifts, it provides additional constraints on the growth rate. As a consequence, more extreme values of $\gamma$ are rejected by the SZ cluster data, resulting in closed contours in the $\gamma$ versus $(1-b)$ plane (blue contours). The situation remains essentially unchanged when we add the X-ray sample data (orange contours). The constraint on the calibration reads $(1-b)=0.602^{+0.053}_{-0.065}$, and the standard calibration $(1-b)=0.8$ is rejected at more than 99%, while $\gamma=0.60\pm 0.13$ is consistent with the growth rate in the standard $\Lambda$CDM model ($\sim 0.545$). This illustrates the fact that even a simple modified version of the growth rate cannot easily accommodate the calibration $(1-b)=0.8$. Again for illustration purposes, we show the consequences of putting a strong prior on $(1-b),$ consistent with standard Planck calibration ($=0.8\pm 0.01$, green contours). The value of $\gamma$ is then found to be constrained to $\gamma=0.917\pm 0.088,$ consistent with conclusions derived from the X-ray sample alone (Sakr et al., 2018).

3.4 Allowing massive neutrinos and modified gravity

Here, we repeat the previous analysis, allowing both massive neutrinos and a modified gravity through the growth index $\gamma$. As can be seen from Fig. 5 and 6, when combining CMB and SZ clusters, a degeneracy is present between the effects of massive neutrinos and modified gravity. The preferred model remains, however, close to the standard picture: $(1-b)=0.651^{+0.043}_{-0.12}$, $\gamma=0.65^{+0.13}_{-0.16}$ and $\sum m_{\nu}\leq 0.197$, although an SZ calibration of $(1-b)\sim 0.8$ is allowed within the two $\sigma$ posterior contours. The addition of X-ray clusters reduces the contours, although in a marginal way. When constraining $(1-b)\sim 0.8$ with a strong prior, the preferred parameters are still quite far away from the standard model: $\gamma=0.834^{+0.080}_{-0.094}$ and $\sum m_{\nu}=0.34^{+0.13}_{-0.15}$ eV, but we note that a model with $\gamma\sim$0.65 and $\sum m_{\nu}\sim 0.1$ is within the two $\sigma$ contours.

3.5 Freeing the amplitude of the late-time matter power spectrum

As a final illustration of what our current data sets actually tell us, we performed the same analysis once again using a combination of CMB and cluster data in the $\Lambda$CDM picture, but allowed the present-day amplitude of matter fluctuations $\sigma_{8}$ to be an additional free parameter, no longer derived from an extrapolation at late times of the primordial amplitude $A_{s}$. More precisely, we assume that the growth rate between redshift zero and $\sim 1$ (the redshift range containing both the SZ and X-ray cluster samples) follows the standard $\Lambda$CDM growth rate, but modified at earlier times through a process about which we do not assume anything. In practice, we rescale the $\Lambda$CDM matter power spectrum for $z\in[0,1]$ by a constant factor, so that we obtain the desired input $\sigma_{8}$ while keeping the same value of the primordial $A_{s}$. Such an approach allows for (comparatively) more freedom than the $\gamma$ parametrisation, as the latter only produces a very specific redshift evolution for the growth rate.

In Fig. 7, we show confidence contours in calibration parameters plane, for example, $(1-b)$ versus $A_{T-M}$, constrained by the combination of CMB data and our two cluster samples: X-ray and SZ. The various contours shown refer to the cases we examined in the previous sections, with the addition of the new model described in this section ($\Lambda$CDM $+$ free $\sigma_{8}$ at present epoch). We observe that all contours lay along the correspondence line333We derived this line using a subset of the Planck SZ clusters sample, for which both temperature and SZ measurements are available. (dashed black) between $(1-b)$ and $A_{T-M,}$ and that in all cases the preferred SZ mass calibration $(1-b)$ remains around 0.6 (see 1D posteriors in Fig. 8) even in the case where $\sigma_{8}$ is treated as a free parameter. In the latter case, a value of $(1-b)\sim 0.8$ is still acceptable (meaning within reasonable confidence limits), but would thus require a significant modification of the growth rate at earlier times in order to reach the required current value of $\sigma_{8}$.

4 Discussion and conclusions

In the present work, we have examined the constraints on the calibration constants in the scaling relations of the Planck SZ and a local X-ray cluster sample in the context of various cosmological models and assumptions. In our approach, we treat these calibrations (namely $(1-b)$ for SZ clusters, and $A_{T-M}$ for X-ray clusters) as free parameters without any priors. Neither massive neutrinos nor a modified gravity model (parametrised by the $\gamma$ growth index) allow us to reconcile the Planck data (CMB and SZ clusters) with the standard SZ mass calibration value of $(1-b)\sim 0.8$. Consequently, this implies that there is no simple solution to the so-called clusters-CMB tension, which may be more accurately described as a tension between Planck data and the empirical calibration of the mass-SZ observable (which yields the $(1-b)\sim 0.8$ value). Despite this observation, we explore the possibility of more extended models allowing massive neutrinos, and a modification of gravity, as well as the approach of freeing the amplitude of matter fluctuations at low redshift (as measured by $\sigma_{8}$). Nevertheless, the SZ cluster counts combined with CMB data still prefer a low calibration ($(1-b)\sim 0.6$) in all considered models. We therefore find it striking that no cosmological model appears to prefer the (comparatively) high calibration of $(1-b)\sim 0.8$.

On a final note, this work was finished within a week of the Planck 2018 CMB likelihood being published, hence the use of the older 2015 data instead. We did check, however, that for the $\Lambda$CDM case notably, using this new CMB data yields a slight difference on the preferred values of $(1-b)=0.622\pm 0.028,$ which still leaves the level of tension above $6\sigma$ (see Fig. 8). Therefore, we do not expect significant changes to the other results of our study and to our conclusions.

Acknowledgements.

SI was supported by the European Structural and Investment Fund and the Czech Ministry of Education, Youth and Sports (Project CoGraDS - CZ.02.1.01/0.0/0.0/15_003/0000437). ZS was supported by a grant of excellence from the Agence des Universités Francophones (AFU). This work has been carried out thanks to the support of the OCEVU Labex (ANR-11-LABX-0060) and the A*MIDEX project (ANR-11-IDEX-0001-02) funded by the “Investissements d’Avenir” French government program managed by the ANR. It was also conducted using resources from the IN2P3 Lyon computing center. We thank the referee for their relevant comments and recommendations.