Non-parametric reconstruction of the cosmological jerk parameter

The cosmological jerk parameter j is reconstructed in a non-parametric way from observational data independent of a fiducial cosmological model. The Cosmic Chronometer data as well as the Supernovae data (the Pantheon compilation) are used for the purpose. The reconstructed values are found to be consistent with the standard Λ\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\Lambda $$\end{document}CDM model within the 2σ\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$2\sigma $$\end{document} confidence level. The model dependent sets like Baryon Acoustic Oscillation and the CMB Shift data are also included thereafter, which does not significantly help in improving or de-proving the confidence level in favour of Λ\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\Lambda $$\end{document}CDM. The deceleration parameter q is also reconstructed from the same data sets. This is used to find the effective equation of state parameter for the model independent datasets only. Λ\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\Lambda $$\end{document}CDM model is excluded for some part of the evolution in 1σ\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$1\sigma $$\end{document}, but is definitely included in 2σ\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$2\sigma $$\end{document} in the domain (0≤z≤2.36\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$0 \le z \le 2.36$$\end{document}) of all the reconstructions.


Introduction
Even after more than a couple of decades of its discovery [1,2], the accelerated expansion of the universe is yet to be attributed to a universally accepted form of dark energy, responsible for the alleged acceleration. Therefore, the quest for dark energy has been alive along all possible ways. A "reverse engineering", where one makes an attempt to find the characteristics of the matter responsible for a particular evolution history, is amongst the prominent ways for quite a long time. Normally this "reconstruction" is related to figure out a physical characteristics of the matter component, such as the equation of state parameter of the dark energy w DE , or even the potential V (φ) if the dark energy is taken as a scalar field.
Another direction of reconstruction is through the kinematical parameters, such as the deceleration parameter q = a e-mail: pm14ip011@iiserkol.ac.in b e-mail: narayan@iiserkol.ac.in (corresponding author) 2 where a is the scale factor, and H = 1 a da dt , the fractional rate of increase in the linear size of the universe called the Hubble parameter. For a long time, H had been the only cosmological parameter which could be estimated from observational data. As H was found to be evolving, the next higher order derivative of a, namely q was the quantity of interest, Now that q can be measured and is found to be evolving, the third order derivative of a finds a natural importance. Expressed in a dimensionless way, this quantity called the "jerk" is defined as There has already been some work in the reconstruction of a cosmological model through these kinematical parameters. Reconstruction of the deceleration parameter q can be found in the work of Gong and Wang [3,4], Lobo et al. [5], Mamon and Das [6], Cardenas and Motta [7], Gómez-Valent [8] and references therein. Very recently Yang and Gong [9] reconstructed the cosmic acceleration and estimated the redshift at which the transition from a decelerated to an accelerated expansion took place using Gaussian Process. Reconstruction through the jerk parameter has been carried out by Luongo [10], Rapetti et al. [11], Zhai et al. [12], Mukherjee and Banerjee [13,14], Mamon and Bamba [15], Mukherjee et al. [16]. Density perturbations also have been investigated for models reconstructed through the jerk parameter [17]. Very recently a CDM model has been recovered from a reconstruction of cosmographic parameters like q and j [18]. Although the possible importance of the jerk parameter in the game of reconstruction was pointed out long back [19], not much work has been done to utilize its full potential. Also, the work already done is an estimation of parameters with a functional form of j being used as an ansatz. This is necessarily restrictive, as the functional form for j is already chosen.
A more unbiased way is to attempt a non-parametric reconstruction, where the evolution of the relevant quantity is determined directly from observational data without any ansatz a priori. Such attempts normally involve the reconstruction of w DE [20][21][22][23][24][25][26][27]. There are some efforts towards a non-parametric reconstruction of q. Some recent examples can be found in works of Nunes et al. [28], Arjona and Nesseris [29], Velten et al. [30] and Haridasu et al. [31]. However, there is hardly any attempt to model the dark energy through a reconstruction of the jerk parameter in a nonparametric way. Although there is no convincing reason that a reconstruction of kinematic parameters like q or j is more useful than that of a physical quantity like the dark energy equation of state parameter, this indeed provides an alternative route towards the understanding of dark energy in the absence of a convincing physical theory.
In the present work, the jerk parameter j is reconstructed for the first time from the observational data in a non-parametric way. We have utilized various combinations of the Supernova distance modulus data, the Cosmic Chronometer (CC) measurements of the Hubble parameter, the Baryon Acoustic Oscillation (BAO) data and also the Cosmic Microwave Background (CMB) Shift parameter data to examine their effect on the reconstruction.
The reconstruction yields the result that for any combination, the CDM model is well allowed within a 2σ confidence level, and for lower z values, within a 1σ confidence level.
Indeed there are apprehensions that the CMB Shift parameter data depends crucially on a fiducial cosmological model [32] and so does the BAO data [33]. However, we do not ignore them. Our reconstruction is based on the combinations both including and excluding the CMB Shift and the BAO datasets. The final result, when we extract the physical information, that of the effective equation of state parameter w e f f , looks qualitatively very much similar for various combinations of the datasets.
In Sect. 2, the methodology is discussed in brief. The details of the observational data used is discussed in Sect. 3. Section 4 contains the actual reconstruction. The last section includes a discussion of the results obtained.

The methodology
At the outset, we do not assume any fiducial model for the universe except that it is given by a spatially flat, isotropic and homogeneous metric given by ds 2 = −c 2 dt 2 + a 2 (t)(dr 2 + r 2 dθ 2 + r 2 sin 2 θ dφ 2 ). (2) We pretend that we do not even know the Einstein equations and pick up only the kinematical quantities. We define the reduced Hubble parameter as, h(z) = H (z) H 0 . A subscript 0 indicates the value of the quantity at the present epoch and z is the redshift given as 1 + z = a 0 a . The luminosity distance of any object (such as a Supernova), can be obtained as For convenience, we define a dimensionless comoving distance, Combining Eqs. (3) and (4) and taking derivative with respect to z, we obtain the relation between Hubble parameter and the comoving distance as, where a prime denotes the derivative with respect to z. In terms of the dimensionless quantities h, D and their derivatives, the cosmological deceleration parameter and the jerk parameter can be written as The uncertainty in q(z), σ q is estimated by the standard technique of error propagation from Eq. (7), Similarly, the uncertainty associated with j (z), σ j is obtained from Eq. (8) as As D (z) is connected to h(z) through Eq. (6), the uncertainty σ D is related to σ h by In order to implement the reconstruction, the widely used Gaussian processes (GP) [34][35][36], which is a powerful model-independent technique, is adopted. This is a distribution over functions which generalizes the idea of a Gaussian distribution for a finite number of quantities to the continuum. Given a set of data points one can use Gaussian processes to reconstruct the most probable underlying continuous function describing the data, and also obtain the associated confidence levels, without assuming a concrete parametrization of the aforesaid function. It requires only a probability on the target function f (z). For a detailed overview one can refer to the Gaussian Process website. 1 In cosmology, GP has attracted a wide application in reconstructing or testing models without an a priori fiducial model [9,24,[37][38][39][40][41][42][43][44][45][46][47]. For a pedagogical introduction to GP, one can refer to Seikel et al. [37]. The publicly available GaPP 2 (Gaussian Processes in Python) code has been used in this work.
Assuming that the observational data, such as the distance data D, or Hubble data H , obeys a Gaussian distribution with a mean and variance, the posterior distribution of reconstructed function f (z) can be expressed as a joint Gaussian distribution of different data-sets involving D or H . In this process, the key ingredient is the covariance function k(z,z) which correlates the values of different D(z) and H (z) at redshift points z andz. The covariance function k(z,z) depends on a set of hyperparameters (e.g. the characteristic length scale l and the signal variance σ f ). This approach also provides a robust way to estimate derivatives of the function in a stable manner. The hyperparameter l corresponds roughly to the distance one needs to move in the input space before the function value changes significantly, while σ f describes the typical change in the function value.
The choice of covariance function, given in (12) affects the reconstruction to some extent. Here we have used the Matérn (ν = 9 2 , p = 4) covariance [34] between two redshift points separated by |z −z| distance units, as in equation (13). This leads to the most reliable and stable results amongst the other significant choices [48].

Observational datasets
Supernova distance modulus data, Cosmic Chronometer data, radial Baryon Acoustic Oscillation data and CMB Shift Parameter data have been utilized in this work.

Cosmic chronometer data
Cosmic Chronometer (CC) H (z) data points are measured by calculating the differential ages of galaxies [49][50][51][52][53][54], as a  [54] function of the redshift z and is given by A complete compilation of the CC data used in this work is enlisted in Table 1.

BAO data
An alternative compilation of the H (z) data can be deduced from the radial BAO peaks in the galaxy power spectrum, or from the BAO peak using the Ly-α forest of QSOs, which are based on the clustering of galaxies or quasars [55][56][57][58][59][60][61][62][63][64][65][66][67]. Table 2 includes almost all data reported in various surveys so far. One may find that some of the H (z) data points from clus- tering measurements are correlated since they either belong to the same analysis or there is an overlap between the galaxy samples. Here in this paper, we mainly consider the central value and standard deviation of the data into consideration. Therefore, we assume that they are independent measurements as in [68].

Reconstructed H 0
After the preparation of CC or/and BAO data, we utilize the GP method to reconstruct the Hubble parameter H (z) and the results are shown Fig. 1. The value of the Hubble constant H 0 obtained from this model independent reconstruction is shown in Tables 3 and 4. Further, we normalize the datasets to obtain the dimensionless or reduced Hubble parameter h(z) = H (z)/H 0 . Considering the error of Hubble constant, we calculate the uncertainty in h(z) as, where σ H 0 is the error associated with H 0 .

SN-Ia data
For the supernova data, we use the recent Pantheon compilation [69] consisting of 1048 SNIa, which is the largest spectroscopically confirmed SN-Ia sample by now. It consists of different supernovae surveys, including SDSS, SNLS, various low-z and also some high-z samples from HST. We include the covariance matrix along with systematic errors in our calculation. The numerical data of the full Pantheon SNIa catalogue and a detailed description is publicly available. 3,4 The distance modulus of each supernova can be estimated as where d L is the luminosity distance as in Eq. (3). The distance modulus of SN-Ia can be derived from the observation of light curves through the empirical relation (Tripp formula [70]) where X 1 and C are the stretch and colour correction parameters, m * B is the observed apparent magnitude and M B is the absolute magnitude in the B-band for SN-Ia while α and β are two nuisance parameters characterizing the luminositystretch, and luminosity-colour relations respectively. M is a distance correction based on the host-galaxy mass of the SN-Ia and B is a distance correction based on predicted Plots for H (z) reconstructed from CC data (left), BAO data (middle), and combined CC + BAO data (right). The solid black line is the "best fit" and the associated 1σ , 2σ and 3σ confidence regions are shown in lighter shades biases from simulations. Usually, the nuisance parameters α and β are optimized simultaneously with the cosmological model parameters or are marginalized over. However, this method is model dependent.
We shall adopt the following strategy. Based on the new approach called BEAMS with Bias Corrections (BBC) [71] the nuisance parameters in the Tripp formula [70] were retrieved and the observed distance modulus is reduced to the difference between the corrected apparent magnitude m B and the absolute magnitude M B i.e., μ = m B − M B . In the Pantheon sample by Scolnic et al. [69], the corrected apparent magnitude m B = m * B + α X 1 − βC along with M and B corrections are reported. We shall avoid marginalizing the over nuisance parameters α, β but marginalize over the Pantheon data for M B in combination with the H (z) data, for the reconstructed value H 0 . The values of α and β thus become irrelevant for the present method of estimation.
With the smooth function H (z) reconstructed from the Hubble data, we use a simple trapezoidal rule [72] for the calculation of the comoving distance, The uncertainty in d c is obtained by error propagation formula, where, Thus, we obtain the smooth function of the comoving distance d c (z) and its error σ d c (z) from the Hubble data. This d c (z) can be rewritten in a dimensionless form by scaling with H 0 c , such that The subscript GP indicates Gaussian Process. Further, we convert the distance modulus of SN-Ia to the normalized comoving distance through the relation (4) For this we perform another Gaussian Process on the apparent magnitudes m B of the SN-Ia data and reconstruct them at the same redshift z as that of the Hubble data. The total uncertainty or error propagation μ and D in μ and D respectively are estimated following the standard practice. The total uncertainty matrix of distance modulus is given by, where C stat and C sys are the statistical and systematic uncertainties respectively. The uncertainty of D(z) is propagated from that of μ and H 0 using the standard error propagation formula, where σ H 0 is the uncertainty of Hubble constant, the superscript 'T ' denotes the transpose of any matrix, D 1 and D 2 are the Jacobian matrices, where D is a vector whose components are the normalized comoving distances of all the SN-Ia. The absolute magnitude of SN-Ia is degenerate with the Hubble parameter H 0 . We get the constraints on M B by minimizing the χ 2 function, considering a uniform prior M B ∈ [−35, −5]. The χ 2 function is given by, where D = (D − D GP ) and = D + σ 2 D GP respectively. The Markov Chain Monte Carlo (MCMC) analysis is performed for maximizing the likelihood function. We adopt a python implementation of the ensemble sampler for MCMC, the emcee, 5 introduced by Foreman-Mackey et al. [73]. The results are plotted using the GetDist 6 module of python, developed by Lewis [74]. Plots for the marginalized M B constraints are shown in the first three columns ([a], [b], [c]) of Fig. 2. The best fit results are given in Table 3.

CMB shift parameter data
The so-called shift parameter is related to the position of the first acoustic peak in the power spectrum anisotropies of the cosmic microwave background (CMB). However the shift parameter R is not directly measurable from the cosmic microwave background, and its value is usually derived from data assuming a spatially flat cosmology with dark matter and cosmological constant, where z c = 1089 is the redshift of recombination. We use the CMB shift parameter R = 1.7488 ± 0.0074 [75] and matter density parameter is marginalized assuming a fiducial CDM model. The χ 2 for CMB Shift parameter data is given by 1

H 0 data
In view of the known tussle between the value of H 0 as given by the Planck data [76], and that from HST observations of 70 long-period Cepheids in the Large Magellanic Clouds by the SH0ES team [77], reconstruction using both of them have been carried out separately. The recent global and local measurements of H 0 = 67.27 ± 0.60 km s −1 Mpc −1 for TT+TE+EE+lowE (P18) [76] and H 0 = 74.03 ± 1.42 km s −1 Mpc −1 (R19) [77] with a 4.4σ tension between them, are considered for the purpose.

The reconstruction
The reconstructed functions h(z), D(z), and their respective derivatives are plotted against z for different sets of the data, and shown in Figs. 3, 4, 5, 6, 7 and 8. The black solid line is the best fit curve, and the black dashed line represents the CDM model with m0 = 0.3. The shaded regions correspond to the 68, 95 and 99.7% confidence levels (CL). The specific points marked with error bars represent the observational data used in reconstruction. For the Pantheon data, Eqs. (23) and (25) are used to estimate the D data points and the uncertainty D from the observed μ and μ respectively. For the CC and BAO data, we consider Eq. (15) and convert the H -σ H data to h-σ h data set. From (6) we can clearly see D (z) is related to h(z). So, we can take into account the h data points, the uncertainty associated σ h , and represent it using Eqs. (6) and (11). Thus, given a set of observational data points we have used the Gaussian Process to construct the most probable underlying continuous function h(z) or D(z) describing the data, along with their derivatives, and have also obtained the associated confidence levels.

Reconstruction of q
Using the reconstructed values of h(z), D(z) and their derivatives in Eq. (7), the deceleration parameter q is now reconstructed. As already mentioned, the Gaussian Process is employed for this. The plots are shown in Fig. 9 for various combination of the datasets. The shaded regions correspond to the 68, 95 and 99.7% confidence levels (CL). The black solid lines show the "best fit" values of the recon-  Although the best fit curve appears to deviate from the a monotonic behaviour, the deceleration parameter as given by the CDM is included generally in 1σ , and at most in the 2σ confidence level.

Reconstruction of j
We now reconstruct the cosmological jerk parameter j using the Gaussian Process from the reconstructed function h(z), D(z) and their higher order derivatives using Eq. (8). If one use the standard Einstein equations with a cold dark matter and a Cosmological constant, j is a constant whose value is unity. Results for the reconstructed jerk is given in Fig. 10. The shaded regions correspond to the 68, 95 and 99.7% con-fidence levels (CL). The black solid lines show the "best fit" values of the reconstructed function, and the black dashed lines corresponds to the CDM model with m0 = 0.3 and j = 1. The plot shows that the CDM model, is allowed within a 2σ error bar. Plots for the "best fit value" of the jerk parameter clearly indicate that j has an evolution, and also, this evolution may well be non-monotonic. The reconstructed best fit q 0 and j 0 at the present epoch is shown in Table 5.
In Fig. 10, we used the value of H 0 , generated by the reconstruction, as described in Sect. 3.3. We further examine if the two different strategies for determining value of H 0 affect our reconstruction differently. We proceed with the analysis similar to that above except we add the P18 or R19 data to the H (z) data tables. For comparison one can refer to Tables 3 and 4 to get an insight as to how including the H 0 measurement affects our reconstruction. Plots for the   reconstructed j (z) along with a comparison for the model independent CC + Pantheon data is shown in Fig. 11.

Reconstructing the effective EOS
We now relax our pretension of not knowing Einstein equations. We use the definition of deceleration parameteṙ in Einstein equations, where ρ and p are the total energy density and pressure contribution from all components constituting the Universe. Therefore, the effective equation of state parameter is Using the reconstructed q(z) for the different data combinations, we arrive at the effective EoS parameter w eff nonparametrically. The evolution of w eff is shown in Fig. 12. The black solid line represent the best fit curve. The shaded regions show the uncertainty associated with w eff corresponding to the 1σ , 2σ and 3σ confidence level. The reconstructed values of w eff at z = 0 is shown in Table 6.
For the CDM model, the dark matter contributes only to the energy density while the cosmological constant con-tributes to both the energy density and pressure. The effective EOS (33) thus takes the following form.
Considering the value of m0 = 0.299 ± 0.013 from the CMB Shift parameter marginalization we can calculate the value of the effective EOS for the CDM model to be −0.701 with ±0.013 uncertainty at z = 0 using the standard error propagation method. For higher redshift (z > 1.5), the reconstructed w eff in the present work indicates a nonmonotonic behaviour. However, the corresponding w eff for the CDM model is included definitely in the 2σ confidence level.

Fitting function for j (z)
In this section we attempt to write an approximate fitting function for the reconstructed jerk parameter in the low redshift range 0 < z < 1 for the model independent CC+Pantheon, CC+Pantheon+R19 and CC+Pantheon+P18 combination. We consider a polynomial function for j (z) with respect to redshift z as,   Fig. 13 Plots showing a comparison between the reconstructed jerk j (z) and the estimated j fit (z) using the combined CC + Pantheon data (left), CC + Pantheon + P18 data (middle) and the CC + Pantheon + R19 data (right). The black solid line is the best fit curve from GP reconstruction. The 1σ and 2σ CL are shown in dashed and dotted lines respectively. The bold dot-dashed line represents the best fit function from χ 2 -minimization. The associated 1σ uncertainty is shown by the shaded region This polynomial is non-linear in z, but it is linear in j i 's. Thus estimating the above equation by the method of least squares or χ 2 minimization holds.
We define the χ 2 function as Here in this work we perform the fitting using a trial and error estimation for different orders of i in Eq. (35). To check the goodness of the fit we calculate the minimized χ 2 for every ith order fitting. The value of the reduced χ 2 ν = χ 2 ν , where ν signifies the degrees of freedom, is estimated. This procedure entails to go from order to order in the polynomial and getting the best-fitting χ 2 , and truncating once the reduced χ 2 falls below one. We start with i = 1 followed by i = 2 and so on and check for which order i the value of χ 2 ν < 1. The measure of χ 2 ν obtained for the three cases studied are mentioned below. Again, the estimated j i 's along with their 1σ uncertainties are given. A comparison between the reconstructed j (z) and estimated j fit are shown in Fig. 13. We also plot the correlations between the parameters j i 's in Fig. 14. We note that the process fails for z > 1, but we can do a reasonable estimate for z < 1. So we show the plots only in the domain 0 ≤ z ≤ 1.
For CC + Pantheon data, The and j 2 = −0.188 +0.403 −0.403 . χ 2 ν = 0.817. In all the three cases, the coefficients in the polynomial are estimated by the χ 2 minimization technique, and the polynomial is truncated once the reduced χ 2 falls below one to prevent over-fitting. If we proceed on fitting with any higher order polynomial, the 1σ uncertainty for the fitted function will not be contained within the 1σ error margin of j (z) reconstructed by GP.

Discussion
The major aim of the present work is a reconstruction of the cosmological jerk parameter j from the observational datasets. The reconstruction is non-parametric, so j is unbiased to any particular functional form to start with. Also, it does not depend on the theory of gravity, only except the assumption that the universe is described by a 4 dimensional spacetime geometry and it is spatially flat, homogeneous and isotropic. It deserves mention that although a non-parametric reconstruction is there in the literature for quite some time now for reconstructing physical quantities like the equation of state parameter or the quintessence potential and the cosmographic quantity like the deceleration parameter q, it has hardly been used to reconstruct the jerk parameter. We have reconstructed q from the datasets, in order to self-consistently reconstruct the j-parameter and the effective equation of state.
Kinematical quantities, that can be defined with the metric (namely the scale factor a) alone, form the starting quantities of interest in the present case. As the deceleration parameter q is now an observed quantity and is found to evolve, the next higher order derivative, the jerk parameter is the focus of attention. Surely the parameters made out of even higher derivatives like snap (4th order derivative of a), crack (5th order derivative) etc. could well be evolving [78]. But we focus on j which is the evolution of q, the highest order derivative that is an observationally measured quantity. For a parametric reconstruction of j, one can still start from the higher order derivatives [79,80] and integrate back to j, and estimate the parameters, coming in as constants of integration, with the help of data. But this does not form the content of the present work as already mentioned. Instead we have directly reconstructed q and j non-parametrically from the datasets, without assuming any kind of parametrizations to start with.
It is found that for various combinations of datasets, the jerk parameter corresponding to the CDM model is included in the 2σ confidence level.
The effective equation of state parameter w eff is linear in q, so the plots for both of them will look similar. We use the reconstruction of q to plot w eff against the redshift z in Fig. 12. The plots indicates that w eff has an evolution which is not necessarily monotonic. The plots also indicate that the universe might have another stint of accelerated expansion in the recent past before entering into a decelerated phase and finally giving way to the present accelerated expansion. For the reconstruction of w eff , the model dependent data sets like BAO and CMB data are not included.
For the reconstructed j, this non-monotonicity is not apparent when only CC or Supernova (Pantheon) data are individually employed, and j decreases slowly for the former and increases rapidly for the latter (Fig. 10). When these two data sets are combined, the non-monotonicity appears and this nature is preserved even when the model dependent datasets like CMB Shift and BAO data are included. It should also be noted that the exclusion of the CMB Shift and BAO does not (i.e., when CC and Pantheon data are used) seriously affect the agreement with CDM within the 2σ error bar.
It may be noted that we obtained the marginalized constraints on H 0 , M B and m0 , by keeping the nuisance parameters α and β fixed using the BBC framework. As there may be correlations between these parameters keeping them fixed may adversely affect the model independent nature of the reconstruction to an extent. However, we have employed the BEAMS with Bias Corrections (BBC) [71] so that the model dependence could be minimized.
There is a very recent work by Bengaly [81] which shows, in a model independent way, that the accelerated expansion of the universe is correct even in 7σ ! So the observational constraints on the kinematic parameters find even more importance. The present work is an attempt to reconstruct the jerk parameter in a non-parametric way. Some of the recent parametric reconstructions of j show that the present value of j indicates that the CDM model is not included in the 1σ level [15,16]. The present work also shows that the evolution of j may exclude the corresponding value in the CDM for a part of the evolution in 1σ , but at a 2σ level j = 1 is indeed included. It deserves mention that a recent model independent study [82] shows that if the GRB (Gamma Ray Burst) data are included, the current value of j is quite different from the standard CDM value of j = 1.
It deserves mention that a very recent work by Steinhardt, Sneppen and Sen [83] points out some errors in the quoted values of the redshift z in the Pantheon dataset. We have worked out the reconstruction of j with the corrected values of z given in [83]. There is hardly any qualitative difference in the plots. The only noticeable difference found is in the lower middle panel of Fig. 10 for the best fit curve where Pantheon + BAO + CMB data are combined. However, even in 1σ , there is no change in the conclusion. The result that CDM is always included in 2σ always remains valid.
Acknowledgements The authors would like to thank the anonymous referee for his constructive suggestions and valuable comments that led to a definite improvement of the work. Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecomm ons.org/licenses/by/4.0/. Funded by SCOAP 3 .