Does $\Lambda$CDM really be in tension with the Hubble diagram data?

In this article, we elaborate further on the $\Lambda$CDM"tension", suggested recently by the authors \cite{Lusso:2019akb,Risaliti:2018reu}. We combine Supernovae type Ia (SNIa) with quasars (QSO) and Gamma Ray Bursts (GRB) data in order to reconstruct a model independent Hubble diagram to as high redshifts as possible. Specifically, in the case of either SNIa or SNIa/QSO data, we find that the cosmokinetic parameters extracted from the Gaussian process are consistent with those of $\Lambda$CDM. Including GRBs, in the analysis, we find a tension, which however is not as significant as that mentioned in \cite{Lusso:2019akb, Risaliti:2018reu}. Finally, we argue that the choice of the kernel function used in extracting the luminosity distance might affect the amount of tension.


I. INTRODUCTION
Since the discovery of the accelerated expansion of the Universe from the Supernovae type Ia (SNIa) data [3,4], the combined analysis of various cosmological probes, including those of Cosmic Microwave Background (CMB) [5][6][7], Baryon Acoustic Oscillation (BAO) [8][9][10][11][12][13][14] and cosmic chronometers [15] confirms the aforementioned dynamical result, namely that currently the Universe accelerates. However, the physics of cosmic acceleration is still a mystery, hence the aim in these kind of studies is to provide an explanation regarding the underlying mechanism which triggers such a phenomenon.
In the framework of homogeneous and isotropic Universe, the accelerated expansion can be described by considering either an exotic matter with negative pressure [16][17][18][19][20][21][22] or a modification of gravity [f (R) theories and the like, [23][24][25][26][27]. Among the large family of dark energy and modified gravity models, the simplest case is the spatially flat ΛCDM model for which cold dark matter (CDM) and baryonic matter coexist with the cosmological constant. From the theoretical viewpoint, the ΛCDM model suffers from the well known problems, namely the coincidence and the expected value of the vacuum energy density [28][29][30][31].
On the other hand, despite the fact that the ΛCDM model is found to be in a very good agreement with the majority of cosmological data [7], nonetheless the model seems to be currently in tension with some recent measurements [32,33], related with the Hubble constant H 0 and the present value of the mass variance at 8h −1 Mpc, namely σ 8 . Moreover, Lusso et al. [1] using a combined Hubble diagram of SNIa, Quasars, and gamma-ray bursts (GRBs) found a ∼ 4σ tension between the best fit cosmographic parameters with respect to those of ΛCDM (see also [2]). In the light of the latter results, a heated debate is taking place in the literature and the aim of the present article is to contribute to this debate.
Here, we focus on a model-independent parametrization of the Hubble diagram using the Gaussian process, and investigate its performance against the latest Hubble diagram data. Notice that in this case we need to introduce a kernel function with some hyperparameters which can be optimized in order to fit the data. For more details concerning model-independent methods we refer to [34][35][36][37]. The structure of the paper is as follows. In section II, we introduce the concept of the Gaussian process and we present the corresponding kernel functions that we shall use in the current work. In section III, we discuss the observational data and the procedure of our analysis, while in section IV we provide our results. Finally, in V, we summarize our results and we draw our conclusions.

II. MODEL INDEPENDENT METHOD-GAUSSIAN PROCESS
We consider that the universe is a self-gravitating fluid, endowed with a spatially flat homogeneous and isotropic geometry. In this context, there are two main approaches in order to investigate cosmological data e.g. the luminosity distance. In the first case we impose a cosmological model, hence we estimate the form of the luminosity distance. Then we fit the model to data in order to place constraints on the corresponding parameter space. This is a model-depended method is a sense that different models provide different forms of luminosity distance. Another avenue is to utilize a model independent method in reconstructing the Hubble diagram through the observational data [34][35][36][37]. In this approach we do not need to know apriori the underlying cosmological model. One of the most popular model independent method is the Gaussian process (GP), hence in the present article we test the performance of GP against the available Hubble diagram data.
Briefly, the main steps of the method are the following. Having a data set D our aim is to reconstruct in a model independent way a function f (x) which describes the data. In this case at any point x, the value f (x) is a Gaussian random variable with mean µ(x) and variance V ar(x). Moreover, the function values at any two different points are not independent from each other, hence the covariance function cov(f (x), f (x)) = k(x,x) describes the corresponding correlations. Therefore, having an observational data set (x i , y i ) and considering a kernel function k(x,x), it is straightforward to compute the value of function and its covariance (for more detail see [38]). Concerning the functional form of the kernel, there is a wide range of possibilities. In the current work we restrict our analysis to the following parametrizations: and Notice that (3) and (4) are the so called Matren (ν = 7/2 and ν = 9/2) formulas respectively. Also σ f and l are two hyperparameters which can be constrained from the observational data. Here we use the GAPP code [38] in order to reconstruct f (x) and its derivatives. Specifically, f (x) and its derivatives are given by where GP stands for Gaussian process.

III. OBSERVATIONAL DATA AND METHOD
The luminosity distance is the ideal tool to investigate the Hubble diagram. Our aim is to extend the Hubble relation to as high redshifts as possible, hence in addition to SNIa, we also consider QSOs and GRBs. In particular, bellow we briefly present the type of standard candles, used in the statistical analysis.
• Quasars (QSOs): Furthermore, we use the sample of 1598 QSOs as collected by [2,40]. The redshift interval of the current data is 0.04 < z < 5.1. Notice that, in our analysis we use bin-averaged version of QSOs data.
• In addition to the above data, we use a compilation of 162 GRBs [41][42][43] in the range of 0.03 < z < 9.3. Unlike SNIa, QSOs and GRBs are observed up to very high redshifts (z > 3)at which the distance modulus is more sensitive to the cosmological parameters [44].
The evolution of the distance modulus is given by µ(z) = 5logD L (z) + 25, hence where D L (z) is the luminosity distance from which the normalized comoving distance 1 is written as Notice that H 0 is the Hubble constant and c is the speed of light. Based on the above, we compute the normalized comoving distance data points and then we use them in order to reconstruct the form of D(z) as well as its derivatives. As a matter of fact knowing D(z) and its derivatives, it is straightforward to compute the Hubble function H(z) as well as its first and second derivatives, namely Moreover using the error propagation we obtain Notice that in above formula, we use the same H 0 which has been used to obtain the normalized distance in Eq.(9) and for those quantities with more than one term in uncertainty, we use square root of all terms. For example, for δX = δa + δb + δc + ..., the total uncertainty is δX = (δa) 2 + (δb) 2 + (δc) 2 + .... Following the same notations we compute the deceleration and jerk parameters as well as the corresponding uncertainties. As a function of D(z), these parameters are: and as a function of H(z), In the case of ΛCDM model, namely where Ω(z) = Ω 0 (1 + z) 3 /E(z) 2 and Ω m0 + Ω Λ0 = 1. Lastly, we remind the reader the basic steps of our method (see section II). First the normalized distances D(z) data are given as input to the GAPP code [38]. Second we reconstruct the functional form of D(z) and finally we compute the rest of the cosmological quantities. During the process we consider that the aforementioned data-sets can be treated as statistically independent measurements. This assumption is a rather strong statement given that for example the SNIa, QSO and GRB data are sensitive to luminosity distances and there might be spatial overlap between the various probes, hence this could lead to correlations that might affect the statistical analysis. While this is an important point, unfortunately at the moment there is no standard way to account for it given the lack of the full correlation matrix among the different samples. Therefore, following standard lines we have assumed that the different data-sets are uncorrelated. Within this framework, the corresponding parameter space is given by (H 0 , q 0 , j 0 ).

IV. RESULTS AND DISCUSSION
In this section, we discuss the main results of our analysis. Specifically, in Table (I) we provide an overall presentation of the cosmographic parameters at the present epoch. In the left panel of Fig.(1) we present the evolution of the reconstructed D(z) and its derivatives when using the Gaussian kernel and SNIa data. As expected D (z) decreases as a function of z, hence due to Eq.(10) the Hubble parameter is an increasing function. In the right panel of Fig.(1) we plot the cosmokinetic parameters H(z)/(1 + z), q(z) and j(z) as a function of redshift. Moreover in the case of Marten kernels ν = 7/2 and ν = 9/2 the aforementioned parameters are shown in Figs.(2) and (2) respectively. We observe that the evolution of the kinetic parameters are almost the same with those of Gaussian kernel.
However, when combined SNIa with other probes, such as GRBs, the situation becomes different. Indeed, for SNIa/GRBs we plot in Fig.(3) D(z) and its derivatives versus redshift using the Gaussian (left panel) and Matren ν = 7/2 (right panel) kernels. For both cases we observe that in the evolution of the corresponding derivatives appears oscillations. It is easy to check that the first derivative of D(z) crosses the zero line several times, hence the cosmokinetic parameters diverge at these points. Notice that utilizing the Matren ν = 9/2 kernel the results remain unaltered. We argue that although GRBs may help to reconstruct the cosmic expansion up to z ∼ 10, however there are practical difficulties in achieving this goal in the case of Gaussian process.
Moreover, the results of SNIa/QSO combination are presented in Fig.(4) and Tab.(I). In this case, we observe that D(z) slowly decreases prior to z ∼ 3, hence a small oscillation appears at that redshift. Furthermore, we find that both Gaussian and Matren kernels provide similar results and in contrast to SNIa/GRB case, here the first derivative of the D(z) does not cross the zero line at 1 σ level.
Lastly, we combine SNIa, GRBs and QSOs in order to compute the reconstructed comoving distance for all kernels. As an example in Fig.(5) we plot the evolution of D(z) and the corresponding derivatives in the case of Matren ν = 7/2 kernel. Again we verify that there are epochs which are located at large redshifts and for which D (z) crosses the zero line (similar behavior is found for the other kernels). Now we focus on Tab. (I) which shows the cosmokinetic parameters at the present time for various data and kernels explored in this study. Considering only the traditional standard candles (SNIa), we find that the Hubble constant is close to 70Km/sec/Mpc regardless the form of kernel, while q 0 and j 0 are consistent (within 1σ) with those of ΛCDM. Combining SNIa and GRB data, we find that the value of H 0 does not change significantly and it remains close to 70Km/sec/Mpc. In the case of Gaussian kernel, the current value of the deceleration parameter is in agreement with that of ΛCDM at 1 σ level. For the Matren's kernels the extracted value of q 0 is marginally consistent with ΛCDM with q 0 < q Λ,0 . Concerning j 0 , our results are similar to those of [1], however the corresponding uncertainties are larger (by a factor of 2.5-4) than those of [1], implying that the extracted jerk parameters are consistent with the the predictions of ΛCDM at 2σ level. Combining SNIa and QSO datasets, we find that for all kernels the cosmokinetic parameters are in a very good agreement (with 1σ) with those of ΛCDM model.
Finally, in the case of the Gaussian kernel the combination SNIa/QSOs/GRBs indicates that the extracted values of q 0 and j 0 are ∼ 3σ away from those of ΛCDM. However, the opposite situation holds in the case of Matren's kernels, namely both q 0 and j 0 are consistent (due to large uncertainties) with the predictions of ΛCDM. In a nutshell, for the usual standard candles (SNIa data) and for the combination SNIa/QSOs we find that the cosmokinetic parameters extracted from the Gaussian process are consistent with ΛCDM. However, including GRBs in the analysis we find a relatively week tension of the ΛCDM model which lies between 2 and 3 σ levels respectively. The fact that the moderate tension appears due to the presence of the GRBs data in the analysis indicates that GRBs appear to be anything but standard candles. Moreover, the combined SNIa/QSO/GRB analysis shows that the choice of the kernel function might affect the amount of tension significantly. Indeed in the case of Matren's kernels we produce cosmokinetic parameters which are consistent with those of ΛCDM, while using the Gaussian kernel it seems that the ΛCDM model is in tension with the measurements (q 0 , j 0 ).

V. CONCLUSION
It is well known that the concordance ΛCDM model fits accurately the current cosmological data [7], nonetheless it has been proposed that the model is not without its problems. Indeed there are indications that the ΛCDM model is in tension with some important measurements [32,33], namely the Hubble constant H 0 and the present value of the mass variance at 8h −1 Mpc, namely σ 8 . In this context, Lusso et al. [1] using a combined Hubble diagram of SNIa, Quasars, and Gamma-Ray Bursts (GRBs) found a ∼ 4σ tension between the best fit cosmokinetic parameters with respect to those of ΛCDM (see also [2]). Whether the above tensions are the result of yet unknown systematic errors or indicate some underlying new Physics is still an open issue. Therefore, on this subject an intense debate is taking place in the literature and the aim of the present work is to contribute to this debate. The same plot using Matren (ν = 7/2) kernel and the same data set.
In particular, we combined the traditional standard candles (SNIa data) with other extragalactic sources (Quasars and GRBs) to reconstruct, in a model independent way, the Hubble diagram to as high redshifts as possible and to compute the corresponding cosmokinetic parameters, namely deceleration and jerk parameters. Using only the SNIa data we found that the cosmokinetic parameters extracted from the Gaussian process are consistent with those of ΛCDM. Also in the case of SNIa/QSO combination, we found that for all kernels the cosmokinetic parameters are in a very good agreement (with 1σ) with those of ΛCDM model.
On the other hand combining SNIa with Quasars and GRBs we revealed some tension, which however is not as significant as that mentioned in [1,2]. Finally, we found that the choice of the kernel function used in extracting the comoving distance might affect the amount of tension. If this is the case, then the ΛCDM cosmokinetic tension is not the result of some underlying new Physics.