Towards a possible solution to the Hubble tension with Horndeski gravity

The Hubble tension refers to the discrepancy in the value of the Hubble constant $H_0$ inferred from the cosmic microwave background observations, assuming the concordance $\Lambda$CDM model of the Universe, and that from the distance ladder and other direct measurements. In order to alleviate this tension, we construct a plausible dark energy scenario, within the framework of Horndeski gravity which is one of the most general scalar-tensor theories yielding second-order equations. In our set-up, we include the self-interactions and nonminimal coupling of the dynamical dark energy scalar field which enable very interesting dynamics leading to a phantom behaviour at low redshifts along with negative dark energy densities at high redshifts. These two features together make this model a promising scenario to alleviate the Hubble tension for appropriate choices of the model parameters. Towards a consistent model building, we show that this set-up is also free from both the gradient and ghost instabilities. Finally, we confront the predictions of the model with low redshift observations from Pantheon, SH0ES, cosmic chronometers and BAO, to obtain best fit constraints on model parameters.


I. INTRODUCTION
One of the most significant discoveries in modern cosmology is the observational evidence of the accelerated expansion of the present Universe, which has been independently inferred from the observations of distant type Ia supernovae [1,2].Such an accelerated expansion is thought to be caused by an unknown yet dominant component of the Universe called dark energy.Within the regime of general relativity (GR), the Standard Model of Cosmology, also known as the Λ-Cold Dark Matter (ΛCDM) model, can explain this phenomenon with the help of a cosmological constant Λ.Despite being the simplest model describing the Universe, the ΛCDM framework is plagued by fundamental issues like the cosmological constant problem and lack of physical understanding regarding the late time acceleration of the Universe, which might be arising due to some unknown dynamics at large scales.There are mainly two approaches that have been followed in order to address these issues -either by considering GR to be the accepted theory of gravity and modifying the dark energy components [3] or by modifying the gravitational theory itself in a way that leads to an accelerated expansion of the Universe, while reducing to GR in certain limits [4].
In the former case, one disregards the cosmological constant Λ and instead introduces a scalar field that mimics the dynamical dark energy, thereby leading to the late time accelerated expansion.On the other hand, a modification of the underlying gravitational theory is usually done by introducing some extra degrees of freedom (DoF) in addition to Einstein's gravity or GR.One of the ways to manifest this is via the scalar-tensor theories of gravity [5] wherein the extra DoF are mixed with curvature, and on subjecting the action to a conformal transformation, they get nonminimally coupled with the matter.The idea is that one should be able to make use of such modified gravity theories to generate extra DoF that mimic the late time acceleration at large scales and hence resolve the problems that the ΛCDM model faces.However, with respect to the local measurements, these DoF should be suppressed so that Einstein's gravity, which agrees very well with observations, becomes the underlying theory of gravity.The required suppression of such extra DoF so as to be consistent with local measurements is usually achieved via various screening mechanisms such as the chameleon [6] and the Vainshtein mechanisms [7].
Over the years, several discrepancies in the standard ΛCDM model have also emerged with the advancement of precision cosmology.These are related to the mismatch in the values of cosmological parameters obtained from direct observations as compared to those inferred from the concordance ΛCDM model.One of the most interesting and statistically significant discrepancies is in the value of Hubble constant H 0 inferred from the cosmic microwave background (CMB) data [8] which assumes a ΛCDM model and those obtained from distance ladder measurements conducted by the SH0ES collaboration [9].Although one can attribute this mismatch to the systematic errors in direct measurement of H 0 , multiple observations from alternative methods have also indicated a tension with H 0 inferred from ΛCDM [10][11][12][13][14][15][16].At present, the observations from SH0ES give a high value of H 0 = (73.3± 1.04) km/s/Mpc, with a 5σ tension with that inferred from Planck.The Hubble tension thus raises a direct question about our understanding of the expansion of the Universe, thereby emphasizing the need to look for possible modifications in the standard ΛCDM picture.Apart from the H 0 tension, the σ 8 (the cosmological parameter related to the clustering of matter in an 8h −1 Mpc radius sphere) tension refers to a discrepancy between the estimation from the CMB and the measurement from surveys like BOSS [17,18].Some other discrepancies are associated with the Ly-α baryon acoustic oscillations (BAO), which provide an independent measurement of the expansion history at higher redshifts.The Ly-α measurement from eBOSS at z ∼ 2.33 prefers a lower value of H(z) than ΛCDM model [19], although with a mild statistical significance of 1.5σ.
Recently, a large number of attempts have been made to alleviate the Hubble tension using different methods.A comprehensive summary of these approaches can be found in [10,[20][21][22][23][24][25][26][27], most of which can be categorised either as an early-time or a late-time modification.For a particular cosmological model, one can obtain the value of H 0 via a measurement of the angular scale of sound horizon θ s = r s (z rec )/d A (z rec ), where z rec is the redshift at recombination, r s (z rec ) is the sound horizon, and d A (z rec ) is the corresponding angular diameter distance.Since θ s is constrained to high precision by CMB data [8], various solutions to the H 0 tension must necessarily respect this constraint.The early-time solutions usually employ a reduction in r s (z rec ), which would require increasing d A (z rec ) in order to keep θ s fixed, automatically implying an increase in the value of H 0 .Some of these approaches include combining the BAO and the Big Bang Nucleosynthesis (BBN) data [28,29], imposing constraints on early dark energy [30][31][32][33][34][35][36][37], hot new early dark energy [38,39], cascading dark energy [40], early dark energy with thermal friction [41,42], and using massive neutrinos [43][44][45].However, these approaches encounter a myriad of problems, and it is understood that early-time solutions are not sufficiently capable of alleviating the H 0 tension [46].
On the other hand, late-time solutions incorporate mechanisms to alter the rate of expansion at redshifts close to the present, which inherently implies z ≪ z rec , thus resulting in a higher value of H 0 without disturbing either r s (z rec ) or d A (z rec ), and hence θ s .The usual prototype to bring about this is to work with a phantom dark energy model, where the equation of state behaves as w < −1 instead of the ΛCDM model (with w = −1).To this end, a few of the numerous attempts include employing "quintessence" potentials to solve the H 0 and σ 8 tensions simultaneously [47], considering a quintessence field that tracks the equation of state during recombination and then boosts the late-time expansion [48], tuning the late-time acceleration of the Universe to cause a monotonous increase in the expansion rate [49], employing a phenomenological two-parameter family of dark energy models [50], interacting dark energy models [51][52][53], using new gravitational scalar-tensor theories that are free of ghosts by definition [54], and many more [55][56][57][58].
In this paper, we introduce a plausible dark energy scenario to provide a possible resolution of the Hubble tension, which is motivated within the framework of Horndeski gravity [59] -one of the most general scalar-tensor theories involving second-order field equations.Horndeski theory provides a rich phenomenology to modify the underlying gravitational theory by introducing an additional degree of freedom, a scalar field, leading to interesting implications in inflationary [60][61][62][63] well as dark energy physics [64][65][66].Recently Horndeski theories have been exploited to address Hubble tension as well [54,67,68].We consider a special case of the Horndeski setup with a dynamical dark energy scalar field ϕ involving self-interaction terms and nonminimal coupling.This leads to interesting features in the dynamics of the late universe and paves the way for a resolution of the Hubble tension.We emphasize two important characteristics of the model, the first being phantom behavior (implying that the dark energy equation of state w ϕ < −1) at present while the other is negative dark energy density (ρ ϕ < 0) at high redshifts.The phantom dark energy models are well motivated in literature to provide a faster expansion phase leading to large H 0 , thereby resolving the Hubble tension [69][70][71][72][73][74], but at the same time being disfavoured by BAO, σ 8 (or S 8 ) measurements, unless accompanied by a phantom crossing behaviour at some epoch [25,75,76].On the other hand, a negative dark energy density at high redshifts has been preferred by various data and observational reconstructions [43,[77][78][79][80][81][82][83], and has been found to be helpful in resolving many anomalies and cosmological tensions like BAO Ly-α measurements, H 0 and S 8 tensions [84][85][86][87][88][89][90].
In this work, we show that with the combined features mentioned above, our model promises to effectively address the H 0 tension.We also comment upon the predictions of the model towards simultaneous alleviation of H 0 and σ 8 tensions based on the discussion by Heisenberg et al. [75].The dark energy scenario proposed in this work is a late-time modification in the dynamics of the universe, thereby not affecting the early universe physics before recombination.Towards a consistent, well-behaved theory, we show that the model is free of gradient and ghost instabilities, which may otherwise arise in more general scalar-tensor theories.Finally, we compare the model predictions against low redshift data from Pantheon, SH0ES, BAO, and cosmic chronometers (CC) [9,[91][92][93] and employ a Markov Chain Monte Carlo (MCMC) analysis to obtain best fit constraints on model parameters.
Our paper is organised as follows: in Section II, we briefly introduce the Horndeski theory and its dynamics and state the conditions for the avoidance of any gradient and ghost instabilities.In Section III, we present our model as a special case of Horndeski Lagrangian with a particular choice of interaction terms.In this section, we further discuss in detail the dynamics and features of the model and confront the model predictions with available data.Finally, in Section IV, we summarize our results and conclude with some discussions and future directions.

II. HORNDESKI THEORY AND ITS DYNAMICS
We work in the framework of the Horndeski theory which is a generalised scalar-tensor theory in four dimensions.The Lagrangian is constructed out of the metric tensor and a scalar field and leads to second order equations of motion, thereby being free of Ostrogradsky instability.The complete Lagrangian for the Horndeski theory (or equivalently, for the generalized Galileons) is given by where where R is the Ricci scalar, G i are four independent arbitrary functions of ϕ and X, and Thus, in the Horndeski gravity, the complete action can be given as, where, L M accounts for the matter (and the radiation) component.Such a theory can well explain the dark energy phase of Universe modelled by a scalar field ϕ coupled to gravity (minimal or non-minimal coupling) depending on the choice of G i (ϕ, X).In a flat FRW background with ds 2 = −dt 2 + a(t) 2 dx 2 , two Friedmann equations are given by [94], In fact, the Eqs. 3 and 4 can be written in a simplified manner as [95], where ρ ϕ and p ϕ are the energy density and pressure of the scalar field given as, In the above equations, κ 2 = 1/M 2 Pl , and from now on, we set κ 2 = 1.Finally, the equation for the evolution of the scalar field is given by [94], where, In the above equations ρ M and p M are the energy density and pressure of the matter field (and radiation), which satisfy the continuity equation, where w M = p M /ρ M is the equation of state for the fluid.The complete background evolution of the Universe can be obtained by solving the above set of equations.However, due to the non-trivial and non-canonical structure of the Horndeski Lagrangian, sometimes these models can have instability issues in the evolution of the perturbations, rendering the background evolution inappropriate.It is therefore necessary to keep a check on the parameters of the theory, even if one is not concerned with the evolution of perturbations for a given dark energy model.The Laplacian or gradient instability is related to the propagation speed of the scalar and tensor perturbations which arise in the regime where the square of sound speed of perturbations becomes negative.This leads to an unstable growth of the perturbation modes on small scales.Another crucial one is the ghost instability which arises when the sign of the kinetic term of the perturbations goes negative.We shall work in the parameter space wherein we avoid all these instabilities so as to have a well-behaved theory.Following the standard technique of linear cosmological perturbation theory, one can obtain the second-order action for the scalar and tensor perturbations in Horndeski theory as [94], where R is the scalar curvature perturbation while h ij are the tensor modes or gravitational waves (GWs), and Q S and Q T are given by Similarly, c S and c T are the propagation speeds of the scalar and tensor modes, respectively, and are given as, and where, For consistent dynamics, free of any gradient and ghost instabilities for both scalar and tensor modes, one must satisfy Another non-trivial characteristic of Horndeski theories is that the speed of GWs can evolve with time for models containing G 4 (ϕ, X) and G 5 (ϕ, X), as can be seen from Eq. 23 and in fact using Eqs.18 and 21, can be stated as The recent observations of GWs from the LIGO-VIRGO collaboration and their electromagnetic counterparts put a very stringent bound on speed of GWs, such that [96,97] which implies that a consistent dark energy model should not violate Eq. 24 in any regime.

III. OUR SCENARIO WITH SELF-INTERACTIONS AND NON-MINIMAL COUPLING
In this section, we present a model constructed within the framework of the Horndeski Lagrangian, comprising of non-zero G 2 , G 3 , and G 4 terms of Eqs. 1 such that, Here we choose, G 4 as a function of the dark energy scalar field ϕ i.e.G 4 (ϕ, X) = G 4 (ϕ) and G 5 (ϕ, X) = 0.This choice ensures that the speed of gravitational waves is always luminal, i.e. c T = 1, as can be verified easily by Eq. 23.Our model as specified by Eq. 25 can also be seen as an extension to the quintessence dark energy scenario [98], where ϕ is the scalar field that drives the dark energy phase, X = φ2 /2 is the kinetic term, and V = V (ϕ) is the potential of the field.Clearly, G 3 = 0, G 4 = 1/2 corresponds to the quintessence case, while ΛCDM scenario can be obtained by choosing For the model given by Eq. 25 the energy density and pressure of the dynamical scalar field can be written using Eqs.7 and 8 as, Therefore, using the above expressions for ρ ϕ and p ϕ in Eqs. 3 and 4 one can obtain the two Friedman equations.Further, the evolution of the dark energy scalar field for our setup of Eq. 25 can be written using Eq. 9 as, Finally, Eqs.5,6 and 28 form the set of three equations that provide the evolution of the cosmological background in the presence of the dynamical dark energy scalar field with energy density and pressure given by Eqs. 26 and 27, respectively.These equations are solved numerically to obtain the evolution of relevant background quantities such as the Hubble parameter H(z), dark energy equation of state w ϕ (z), etc.Throughout the analysis, we fix the potential of the dark energy scalar field to a linear one, i.e., V (ϕ) = V 0 ϕ.The initial conditions on H, ϕ, and φ are set at the time of recombination, i.e. at z rec = 1100.We set H(z rec ) = H ΛCDM (z rec ) in H 0 units, while we choose ϕ(z rec ) = 2 and φ(z rec ) = 10 −3 in Planck units.Also, the parameter V 0 (of dimension [M Pl ] 3 ) of the dark energy potential is kept fixed throughout the analysis, with V0 = 1.1 in H0 units.
The three crucial model parameters in our analysis are c1, c2, and c3, which control the strength of couplings as can be seen from Eq. 25.Specifically, c1 and c2 contained inside G3 control the self-interaction terms of the scalar field and its derivative, while c3 in G4 controls the strength of nonminimal coupling to gravity.These two interactions arising from G3 and G4, respectively, lead to distinct effects in the dynamics of the expansion of the universe.Specifically, the G3 term in Eq. 25, with an appropriate choice of c1 and c2 can give a phantom behavior (w ϕ < −1) at low redshifts leading to a faster expansion compared to ΛCDM, while the nonminimal coupling to gravity due to the choice of G4 in Eq. 25, leads to negative dark energy density at high redshifts (around z ≥ 3).This can be qualitatively seen from Eq. 26, where the appropriate choice of c3(> 0) can lead to ρ ϕ < 0 in regimes where V (ϕ) is subdominant, given that φ is really small.In Sec.III B, we shall discuss these effects in detail.As mentioned in Sec.I, negative dark energy at high redshifts has been found to be effective in resolving anomalies and cosmological tensions.In the next section, we will demonstrate these effects for different choices of c1, c2, and c3.Also, the dimensions of the three model parameters, c1, c2, c3 in terms of Planck mass are [M Pl ] 0 , [M Pl ] −3 and [M Pl ] 1 respectively.The case with c1 = c2 = c3 = 0 gives the usual quintessence scenario, with a background similar to the ΛCDM.
Another interesting feature of the Horndeski theory involving a non-trivial G4 term is that it can lead to an evolution of the Planck mass on cosmological scales.This leads to a modification in the luminosity distance compared to the standard GR due to additional friction arising from the running Planck mass [99], referred to as the GW distance.Future GW distance observations from the Einstein Telescope will constrain the deviations from GR, arising in modified gravity theories like Horndeski gravity [100].Furthermore, the presence of Galileon-type interaction in our model, i.e.G3(X), leads to the Vainshtein mechanism, thereby reconciling with GR inside the Vainshtein radius [101].

A. Alleviating the Hubble tension -our model results
In this section, we present the results of the numerical analysis of the cosmological background for our model given by Eq. 25.As mentioned earlier, we fix the initial conditions on ϕ, φ, and H at recombination zrec, along with a fixed value of V0, to solve the set of background equations numerically.The top left panel of Fig. 1 shows the evolution of the normalized Hubble parameter H(z) for different values of c1, c2 and c3 while the top right panel shows the corresponding evolution of the dark energy equation of state w ϕ (z).The results of the standard ΛCDM model are also displayed for comparison.In order to visualize the salient features of the model, we also confront the corresponding evolution of the Hubble parameter with observational data from BAO (bottom left) and cosmic chronometers (bottom right), as shown in Fig. 1.The model can favour larger values of the Hubble parameter at present in comparison to the ΛCDM case for specific choices of model parameters.This can be seen in the top left panel of Fig. 1 that the model gives a comparatively faster expansion at low redshifts than ΛCDM, thereby large H0.This can clearly be attributed to the phantom behavior of the dark energy scalar (w ϕ < −1) field around the present epoch in the top right panel of Fig. 1.
However, at higher redshifts, the H(z) for the model is less than that of ΛCDM, which is due to the negative energy density of the scalar field around those redshifts.This feature of ρ ϕ < 0 at high redshifts is traced in the evolution of w ϕ in the right panel of Fig. 1, which crosses a singularity around the epoch (z ∼ 2 − 5) where ρ ϕ switches sign.Around this point, the dark energy density transits to positive values, giving appropriate evolution at present.This type of singularity in w ϕ (z), arising due to negative dark energy in some regime, has been studied in literature to address anomalies in Ly-α BAO measurements In the bottom panel, the model is confronted with the measurement of H(z) from SN [9], BAO, and CC data [92,93].
of H(z), as well as supported by the observational reconstructions [78,84,85].Recently it has been pointed out ( [80,83] and references therein) that the negative energy density of dark energy scalar field ρ ϕ at high redshifts, wherein the effects of dark energy are really subdominant or insignificant, can lead to resolution of cosmological tensions like H0, S8 etc. Indeed, the total energy density of the universe is always positive, even when ρ ϕ is negative at higher redshifts.Since the contribution of ρ ϕ to the total energy density remains really insignificant at earlier epochs, the singularity in the equation of state does not affect the evolution of relevant background entities like H, Ḣ, etc. Thereby the two crucial characteristics of the dark energy component in our model i.e. the phantom behavior at low redshifts and negative energy density at higher redshifts, collectively are capable of explaining the observed measurement of Hubble parameter from independent observations like Supernovae, BAO and CC as demonstrated in the lower panel of Fig. 1.To illustrate further, these two characteristics of the model lead to a sign switching in δH(z) (= H(z) − HΛCDM(z)), i.e. δH(z) < 0 at high redshifts while δH(z) > 0 around the present, which is seen as a necessary condition to be satisfied by any late dark energy model trying to address the Hubble tension [75,90].These conditions are obtained to ensure that any late-time modifications do not disturb the precisely measured angular size of the sound horizon around recombination by Planck.In Sec.III B, we shall further illustrate the distinct effects arising from the individual interaction terms, i.e. from G3 and G4 of Eq. 25, one at a time.
For a consistent model building, it is necessary to ensure that there are no instabilities, as discussed in Sec.II.The necessary conditions to avoid such instabilities in a Horndeski setup are given by Eq. 22.As shown in Fig. 2, from the evolution of the speed of perturbations of scalar modes cS and the parameter QS, it is evident that our model is free of gradient and ghost instabilities.Moreover, as discussed earlier in Sec.III, for our setup given by Eq. 25, the speed of tensor perturbations is always luminal thereby, we do not worry about tensor perturbations.
Finally, we comment upon the implications of such a model towards resolution of σ8 (or S8) tension as well, in reference to the points made by Heisenberg et al. in [75].The former reference discusses the necessary conditions in the background evolution of a dark energy model for simultaneous alleviation of H0 and σ8 tensions in a model-independent way.Our model meets the criteria suggested by [75], i.e., it exhibits a phantom crossing behaviour.Thereby, without going to the evolution of perturbations, we highlight that our model indeed shows the possibility of resolving the σ8 tension as well.Due to the complicated structure of perturbations in a Horndeski setup, we leave the evolution of perturbations for future work.B. Understanding the effects of the choice of G3 and G4 in our model As discussed earlier, the proposed dark energy model given by Eq. 25 involves non-zero G3 and G4 terms.From Eq. 1, one can see that a non-trivial G3 invokes self-interactions of dark energy field while the choice of G4 introduces non-minimal coupling to gravity.In this section, we will try to understand the individual effects of G3 and G4 on the dynamics of our model.As mentioned in Sec III, the role of G3 is to lead to a phantom behaviour at low redshifts.The condition for a phantom behaviour is given by For the case when G4 is turned off, putting c3 = 0 in Eqs. 26 and 27, the above condition resorts to, φ2 (1 − 2c1 − c2 φ + 3c2H φ) < 0 which can be achieved with appropriate choices of c1 (predominantly) and c2.A non-zero c2 is necessary to ensure the stability condition Qs > 0 is not violated in any regime.It is to be noted that with c3 = 0, one cannot go into the ρ ϕ < 0 regime as can be seen from Eq. 26, since the terms O( φ2 ) < V (ϕ).Now, for the case when only the non-minimal coupling term G4 is switched on, putting c1 = c2 = 0, one gets which can give rise to ρ ϕ < 0 for a suitable choice of c3 whenever V (ϕ) < 6c3ϕH 2 in some regimes, where we ignore the terms O( φ) being very small.Since we are working with a linear potential, the above condition can be given as V0 < 6c3H 2 , or V0 < 2c3ρc, where ρc = 3H 2 .This can easily be satisfied at high redshifts where dark energy density is really insignificant.We now show these results obtained from numerical solutions of the background equations Eqs. 5, 6 and 28 for these two cases discussed above.The left panel in Fig. 3 shows the evolution of H(z)/(1 + z) for cases where either G3 or G4 is turned on.Note that we do not vary the scalar field mass (quantified by V0) throughout the analysis.For the case with G3 turned on (i.e. with c3 = 0), one can obtain a higher value of H0 arising due to the phantom behaviour of dark energy around the present epoch.This can be visualized from the green and brown plots in the evolution of H(z) and w ϕ (z) in Fig. 3.But at the same time, the equation of state does not cross the phantom divide or become singular in any regime; thereby, such a setup may not address the BAO Ly-α measurements of H(z) at high redshifts.In fact, it will not satisfy the condition of δH(z) < 0 at high redshifts, as required for successful alleviation of the Hubble tension [75,90].Such type of scenarios where the dark energy equation of state makes a sudden transition to phantom behaviour at very low redshifts, so-called "hockey stick" models, are not considered suitable (see for example [102]).In the other case, when only G4 term is turned on (by putting c1 = c2 = 0), we do not obtain H0 significantly higher than ΛCDM.This is also understood from the evolution of w ϕ (z) in the red and blue curves in Fig. 3, which rather behaves like a quintessence field at present redshift with w ϕ > −1.Interestingly, such a setup gives rise to a singularity in w ϕ (z), which leads to a phantom crossing behaviour, the position of which can be controlled depending on the choice of parameter c3.Because of this feature in w ϕ , arising due to negative dark energy at high redshifts, such setup can explain the H(z) observations from BAO Ly-α measurements.In conclusion, the phantom behaviour of the model can be attributed to the G3 term (controlled by c1 and c2) while the negative dark energy at high redshifts arises due to the choice of G4 (controlled by c3) in our model given by Eq. 25.The former effect gives rise to a large H0 at present, explaining the observations from Supernovae [9] while the latter gives the necessary phantom crossing condition to successfully resolve H0 tension, in addition to explaining the measurements of H(z) from BAO and SDSS at higher redshifts [92,93].Hence, the combined characteristics drive the model towards a successful alleviation of the Hubble tension.FIG.3: The evolution of H(z)/(1 + z) confronted with measurement of Hubble parameter from [9], BAO, and CC data [92,93] (left) and the corresponding equation of state of dark energy scalar field (right) for the cases where either of G 3 and G 4 are turned on to see the effects of these individual terms on the background cosmological dynamics.

C. Best fit constraints on model parameters
In this section, we assess our model against the low redshift data to obtain constraints on the model parameters.This is achieved by employing a Markov Chain Monte Carlo (MCMC) analysis, which samples the parameter space of c1, c2 and c3, to obtain best fit constraints on these parameters and thereby on the derived background quantities of interest (here H0, ΩM ) for different likelihoods.The data used in our analysis are as follows: • The SH0ES measurement of H0 = 73.3± 1.04 km/s/Mpc [9], modelled with a Gaussian likelihood.
• The expansion rate measurements from BAO, compiled in Table I.
• The cosmic chronometer (CC) data on H(z) as compiled in Table III of [103] and the covariance matrix obtained [104] following the method discussed in [105].
For the sampling process, the three free model parameters are provided with uniform priors as summarised in Table II.In Fig 4 the results of the analysis are presented as posterior distributions of the model parameters and the derived background quantities for individual data (top panel) and for the combined data (bottom panel) obtained using the GetDist Python package [109].We also tabulate the mean values and 1σ constraints on these parameters obtained from the MCMC sampling in Table III.
For the case of individual likelihoods, the posteriors of parameters c1 and c2 overlap to a great extent for all the cases.Interestingly for the parameter c3, which controls the strength of nonminimal coupling and causes the dark energy density to possess negative values at high redshifts, the constraints from BAO differ largely from those obtained for the other three data, namely CC, Pantheon and SH0ES.The BAO constraints seem to prefer a large value of c3, thereby inclined to strongly favour the presence of negative dark energy at high redshifts and, thus, explain the anomalous Ly-α measurement of expansion rate, as also discussed in [84,107].But unlike BAO, the other data (CC, Pantheon and SH0ES) seem to be less sensitive to the nonminimal coupling parameter c3, as evident from Table III, only providing an upper bound on it.This can also be understood, given that these data mostly constrain the expansion history at low redshifts (z < 2), while the negative dark energy density feature is effective at somewhat larger redshifts.
A similar behaviour is seen in the joint likelihood analysis (with all data combined), where the parameter c3 favours very small values with an upper bound as shown in Table III, indicating that the combined data does not constrain well the nonminimal coupling parameter.But any non-zero value of c3, irrespective of its amplitude, implies the presence of negative dark energy.As discussed earlier, this incorporates interesting features like phantom crossing in the model, which is essential in order to successfully address H0 and σ8 tension without disturbing the CMB peaks [75].The other two parameters c1 and c2 are well-constrained with nearly Gaussian posteriors and similar mean values for different data.This indicates the preference for a phantom behaviour of dark energy at low redshifts (due to the non-zero c1 and c2), which leads to a large value of H0 ∼ 72 km/sec/Mpc, thereby reducing the tension with local measurements.However, in order to obtain the complete picture, one needs to include the CMB data as well to obtain parameter constraints, which is beyond the scope of the present work.Indeed, including the CMB data would help us impose more stringent constraints, especially on the nonminimal coupling parameter c3, indicating whether or not the data favours the presence of negative dark energy at high redshifts.The mean values of all the parameters along with their 1-σ constraints are listed here.

IV. CONCLUSIONS AND DISCUSSIONS
The Hubble tension is certainly one of the most intriguing problems in modern cosmology and has recently received a lot of attention.Numerous attempts are ongoing as well as have been made in order to bridge the gap between the ΛCDM-predicted value of H0 and the values obtained through the late-time (low redshift) observations.While the mismatch can certainly be attributed to the systematic errors associated with different direct or indirect (or both) measurements, the tension may also be a true reflection of reality, requiring exotic new physics and, perhaps, a dramatic revision to our current understanding of cosmic evolution on large scales.
In this work, we have proposed a new dark energy model within the framework of Horndeski gravity as a plausible scenario towards the resolution of H0 tension.The dark energy expansion is governed by a dynamical scalar field involving non-trivial self-interactions and non-minimal coupling, motivated by the Horndeski Lagrangian.The proposed model resorts to late-time modification approaches to resolve H0 tension.The two crucial characteristics of our model, which are negative dark energy density at high redshifts (giving rise to phantom crossing) and a phantom behaviour around the present epoch, make it an appealing candidate among various solutions to alleviate the Hubble tension.As discussed in detail, both these conditions suffice for explaining the observed measurement of the Hubble parameter from independent observation and at the same time, ensuring that the CMB measurements of the angular size of the sound horizon from Planck remain unaffected.We have further shown that our scenario remains free of various issues, such as gradient instability or superluminal propagation and, thus, is consistent from a physical model-building perspective.We have also briefly mentioned that our model shows the possibility of resolving the σ8 tension as well due to the presence of a phantom crossing scenario [75].To check the feasibility of our model, we use low redshift data from Pantheon, SH0ES, BAO and CC, and employing an MCMC analysis, obtain the best fit constraints on model parameters.The analysis shows that the proposed dark energy model can support larger values of H0(∼ 72 km/sec/Mpc) due to the phantom behaviour of the dark energy field at low redshifts (due to self-coupling parameters c1 and c2), thereby resolving the H0 tension.
In fact, the next relevant step in this direction would be to study the evolution of linear density perturbations in such Horndeski models to see the effects on CMB and the matter power spectrum.We consider this as a future prospect to perform complete data analysis, including full CMB data, and explore the parameter space of such late Universe modifications arising within the Horndeski theory, along with examining the possibility of alleviating the σ8 tension simultaneously.Finally, it would be very interesting to construct other consistent models within this framework, which can simultaneously alleviate the H0 and σ8 tensions without affecting the CMB measurements.We leave these fascinating possibilities for future work.

c 1 c 1 c 1 1 ΛCDMFIG. 1 :
FIG.1:The evolution of the normalized Hubble parameter (top left) and equation of state of dark energy scalar field (top right) for different combinations of c 1 , c 2 , and c 3 .In the bottom panel, the model is confronted with the measurement of H(z) from SN[9], BAO, and CC data[92,93].

c 1 1 FIG. 2 :
FIG. 2:The evolution of the speed of scalar perturbations c s (left) and stability parameter Q s (right) for a specific choice of model parameters c 1 , c 2 and c 3 .The plots indicate that the model is free of gradient and ghost instabilities as specified by conditions in Eq. 22.

TABLE I :
Baryon acoustic oscillation (BAO) expansion rate measurements from SDSS collaboration.The references of these measurements are provided in the table.

TABLE II :
The priors on model parameters used in our analysis.The prior range is the same for all the likelihoods used.

TABLE III :
We present the best fit constraints on various model parameters after MCMC analysis with different low redshift data.