# A model-independent reconstruction of dark sector interactions

Rodrigo von Marttens  Observatório Nacional, 20921-400, Rio de Janeiro, RJ,
Brasil  Javier E. Gonzalez  Facultad de Ciencias e Ingeniería, Universidad de
Manizales, 170002, Manizales, Colombia  Departamento de Física, Universidade
Federal do Rio Grande do Norte, 59072-970, Natal, RN, Brasil  Jailson Alcaniz
Observatório Nacional, 20921-400, Rio de Janeiro, RJ, Brasil  Valerio Marra
Núcleo de Astrofísica e Cosmologia (Cosmo-ufes) & Departamento de Física,
Universidade Federal do Espírito Santo, 29075-910, ES, Brasil  INAF -
Osservatorio Astronomico di Trieste, via Tiepolo 11, 34131 Trieste, Italy
IFPU – Institute for Fundamental Physics of the Universe, via Beirut 2, 34151,
Trieste, Italy  Luciano Casarini  Departamento de Física, Universidade Federal
de Sergipe, 49000-000 São Cristóvão, SE, Brasil  Institute of Theoretical
Astrophysics, University of Oslo, P.O. Box 1029 Blindern, N-0315 Oslo, Norway

(August 29, 2024)

###### Abstract

Relaxing the conventional assumption of a minimal coupling between the dark
matter (DM) and dark energy (DE) fields introduces significant changes in the
predicted evolution of the Universe. Therefore, testing such a possibility
constitutes an essential task not only for cosmology but also for fundamental
physics. In a previous communication [Phys. Rev. D99, 043521, 2019], we
proposed a new null test for the $\Lambda$CDM model based on the time
dependence of the ratio between the DM and DE energy densities which is also
able to detect potential signatures of interaction between the dark
components. In this work, we extend that analysis avoiding the $\Lambda$CDM
assumption and reconstruct the interaction in the dark sector in a fully
model-independent way using data from type Ia supernovae, cosmic chronometers
and baryonic acoustic oscillations. According to our analysis, the
$\Lambda$CDM model is consistent with our model-independent approach at least
at $3\sigma$ CL over the entire range of redshift studied. On the other hand,
our analysis shows that the current background data do not allow us to rule
out the existence of an interaction in the dark sector. Finally, we present a
forecast for next-generation LSS surveys. In particular, we show that Euclid
and SKA will be able to distinguish interacting models with about 4% of
precision at $z\approx 1$.

Cosmology: theory – dark energy – dark matter – large-scale structure

##  I Introduction

The standard $\Lambda$-cold dark matter ($\Lambda$CDM) model provides a
remarkably successful description of the universe on large scales Aghanim  _et
al._ (2020); Abbott  _et al._ (2019a, b). In this model, the dark energy and
dark matter components are minimally coupled to each other and dominate the
structure and evolution of the universe in late times. Recently, the validity
of the $\Lambda$CDM model has been questioned given the discrepancies and
tensions between early and the late time Universe measurements Verde  _et al._
(2019); Riess (2019); Asgari  _et al._ (2020), and a possible alternative
could be models with a non-minimally coupling between dark energy and dark
matter Zimdahl and Pavon (2001); Alcaniz and Maia (2003).

As it is well known, there is no known fundamental principle that prevents a
coupling term $Q$ between the energy components of the cosmological dark
sector. Physically, $Q$ represents an energy exchange between the dark
components, which necessarily violates adiabaticity at the same time that it
brings about important consequences on the model predictions. Currently, a
number of analysis show that some particular classes of interacting models are
able to provide a good description of the data von Marttens  _et al._ (2019a);
Benetti  _et al._ (2019); Xia and Wang (2016); Kumar and Nunes (2017); Yang
_et al._ (2019a); Martinelli  _et al._ (2019); Cid  _et al._ (2019); Yang  _et
al._ (2017); vom Marttens  _et al._ (2017a, 2014); da Silva  _et al._ (2020);
Gonzalez  _et al._ (2018), sometimes alleviating some of the tensions of the
standard cosmology Di Valentino  _et al._ (2017); Yang  _et al._ (2018a);
Kumar  _et al._ (2019); Pan  _et al._ (2019a); Di Valentino  _et al._ (2020);
Yang  _et al._ (2018b); Pan  _et al._ (2019b, 2020). From the theoretical
point of view, however, critiques to these models do exist, being mainly
related to the absence of a natural guidance from fundamental physics on the
coupling term, which leads to a phenomenological choice of $Q$ Valiviita  _et
al._ (2008); Costa and Alcaniz (2010); Cai and Su (2010); Carneiro  _et al._
(2019).

In this paper, instead of assuming a parameterization of $Q$ a priori, we take
a different route and reconstruct physical quantities directly related to the
coupling term from observational data. Different reconstruction methods have
been widely used in Cosmology in different contexts Zhao  _et al._ (2017);
Joudaki  _et al._ (2018); L’Huillier  _et al._ (2018); Shafieloo  _et al._
(2018); L’Huillier  _et al._ (2019); Keeley  _et al._ (2019, 2020); Liao  _et
al._ (2020); Keeley  _et al._ (2021). The method we adopt in this work is the
Gaussian Process Seikel  _et al._ (2012).

Extending previous results von Marttens  _et al._ (2019b), we employ a
reconstruction method to map the evolution of the interaction in the dark
sector over a large range of redshift in a model-independent way using data
from type Ia supernovae, cosmic chronometers and baryonic acoustic
oscillations. Our results show a good agreement with the standard cosmology
($Q=0$), although the existence of an interaction in the dark sector cannot be
ruled out. For completeness, we also perform a forecast analysis for some
next-generation galaxy surveys and discuss their ability to constrain the
possibility of a non-minimally coupling between the dark components.

This work is outlined as follows: In Sec. II we present a general framework
for describing the cosmological background with an unified dark sector. In
particular we discuss how the two important features affect the dynamics of
the dark sector: the dynamical character of the DE component and the
possibility of the dark components interact with each other. In Sec. III we
develop the formalism for depicting a general interacting dark sector. A
slight assumption that the dark sector interaction depends on the energy of
the involved components is made, but no model is specified. In Sec. IV we
finally examine the main point of the paper: how to perform a model-
independent reconstruction of a dark sector interaction. We develop the
relevant equations for the reconstruction, analyze the data employed and
present the results obtained from $H\left(z\right)$ measurements (Cosmic
Chronometers and BAO) and type Ia SNe (Pantheon). Sec. V is focused to shed a
light on how future surveys will be able to improve the results of the
proposed analysis. In this analysis we use forecasts for measurements of
$H\left(z\right)$ from J-PAS Benitez  _et al._ (2014); Bonoli  _et al._
(2020), DESI Aghamousa  _et al._ (2016), Euclid Amendola  _et al._ (2013) and
SKA Bacon  _et al._ (2020). Lastly, Sec. VI is dedicated to the concluding
remarks.

##  II Unified dark sector: background description

Let us consider the FLRW line element ($c=1$)

$ds^{2}=dt^{2}-a^{2}\left(t\right)\left[dr^{2}+r^{2}\left(d\theta^{2}+\sin^{2}\theta
d\phi\right)\right]\,,$ (1)

where $a$ is the scale factor. For the material content, we assume a Universe
composed of four components: radiation (denoted by the subindex $r$), baryons
($b$), cold dark matter ($c$) and dark energy ($x$). At the background level,
each one of these components is described as a fluid with an equation of state
(EoS) $p_{i}=w_{i}\rho_{i}$, and while the radiation and baryonic components
are separately conserved, DM is allowed to interact with DE. Initially, we do
not assume any particular parameterization for the DE component, which means
that $w_{x}$ is allowed to be a function of the scale factor. As will be
discussed later, we will also consider the specific case $w_{x}=-1$, which, in
the absence of interaction, is equivalent to a cosmological constant
$\Lambda$.

Within the general relativistic framework and assuming spatial flatness, the
background dynamics is given by the Friedmann equations,

$\displaystyle\left(\dfrac{\dot{a}}{a}\right)^{2}\equiv H^{2}=\dfrac{8\pi
G}{3}\rho\,,$ (2a) $\displaystyle\dfrac{\ddot{a}}{a}=-\dfrac{4\pi
G}{3}\left(\rho+3p\right)\,,$ (2b)

where $H\equiv\dot{a}/a$ is the Hubble rate and the absence of an index in
$\rho$ and $p$ indicates that both quantities refer to the cosmic substratum.
Radiation and baryons satisfy the usual background energy conservation
equation,

$\displaystyle\dot{\rho}_{r}+4H\rho_{r}=0\qquad$ $\displaystyle\Rightarrow$
$\displaystyle\qquad\rho_{r}=\rho_{r0}\,a^{-4}\,,$ (3a)
$\displaystyle\dot{\rho}_{b}+3H\rho_{b}=0\qquad$ $\displaystyle\Rightarrow$
$\displaystyle\qquad\rho_{b}=\rho_{b0}\,a^{-3}\,,$ (3b)

where the subindex $0$ denotes that the corresponding quantities are evaluated
at $a=a_{0}=1$.

Before treating the dark components separately, let us now combine the DM and
DE components in order to describe the dark sector as an effective unified
dark component. For the unified description, the total dark energy density is
the sum of DM and DE energy densities whereas the total dark pressure is the
DE pressure,

$\rho_{d}=\rho_{c}+\rho_{x}\qquad{\rm and}\qquad p_{d}=p_{c}+p_{x}=p_{x}\,.$
(4)

From the above expression we also find

$p_{d}=\left(\frac{w_{x}}{1+r}\right)\rho_{d}\,,$ (5)

where $r$ is defined as the ratio between CDM and DE energy densities
($r\equiv\rho_{c}/\rho_{x}$). In the FLRW cosmology context, both functions
$w_{x}$ and $r$ must depend on the scale factor, and the term in the
parentheses of Eq. (5) can be seen as an effective EoS parameter of the
unified dark component, i.e.,

$w_{d}\left(a\right)=\frac{w_{x}\left(a\right)}{1+r\left(a\right)}\,.$ (6)

Such a unified dark component must satisfy the energy conservation equation,

$\dot{\rho}_{d}+3H\big{[}1+w_{d}\left(a\right)\big{]}\rho_{d}=0\,.$ (7)

which, for a general time-dependent EoS parameter $w_{d}\left(a\right)$, has
the well-known solution,

$\rho_{d}=\rho_{d0}\exp\left[-3\int\frac{1+w_{d}\left(a\right)}{a}da\right]\,.$
(8)

By splitting the total energy density in radiation, baryons and dark sector,
Eq. (2b) can be rewritten as,

$3H^{2}=8\pi G\left(\rho_{r}+\rho_{b}+\rho_{d}\right)\,,$ (9)

where $\rho_{r}$, $\rho_{b}$, and $\rho_{d}$ are given respectively by Eqs.
(3a), (3b) and (8). Eqs. (8) and (9) show that all information about the dark
sector contained in the background expansion (Hubble rate and, consequently,
any distance measurement) comes from $w_{d}(a)$. We refer the reader to Ref.
von Marttens  _et al._ (2020) where it was pointed out as a way to establish
an explicit mapping between dynamical DE models and interacting dark sector
models, so that they have identical Hubble rates. In practice, this mapping
relates models that measure distances identically, and then, they can not be
distinguished by distance-based observables.

The unified EoS parameter $w_{d}$ also contains all the dark sector’s
contribution to the deceleration parameter

$q\equiv-\dfrac{\ddot{a}a}{\dot{a}^{2}}=\dfrac{1}{2}\sum_{i}\Omega_{i}\left(1+3w_{i}\right)\,,$
(10)

where $\Omega_{i}\equiv 8\pi G\rho_{i}(a)/3H^{2}(a)$ is the density parameter
of the $i$-th component ($i=r$, $b$ and $d$). Taking all the components into
account separately, Eq. (10) reduces to

$\displaystyle q$ $\displaystyle=$
$\displaystyle\dfrac{1}{E^{2}}\bigg{[}\Omega_{r0}a^{-4}+\dfrac{1}{2}\Omega_{b0}a^{-3}+\frac{1}{2}\Big{(}E^{2}-\Omega_{r0}a^{-4}$
(11) $\displaystyle-\Omega_{b0}a^{-3}\Big{)}\left(1+3w_{d}\right)\bigg{]}\,,$

where $E\equiv H/H_{0}$. From Eq. (11), it is straightforward to obtain that
the acceleration condition ($q<0$) can be also formulated in terms of $w_{d}$
as,

$\displaystyle
w_{d}<-\frac{1}{3}\dfrac{E^{2}+\Omega_{r0}a^{-4}}{E^{2}-\Omega_{b0}a^{-3}-\Omega_{r0}a^{-4}}\,.$
(12)

From Eq. (6), one can identify that the unified dark EoS parameter $w_{d}$ has
two time-dependent degrees of freedom: $w_{x}\left(a\right)$ and
$r\left(a\right)$, being each one of them related to a dynamical feature of
the dark sector. The DE EoS parameter is related to the dynamical nature of
the DE component, as will be seen in more detail in Sec. III. Assigning a
dynamic behavior to the DE component is one of the most common alternatives to
the standard cosmological model Ratra and Peebles (1988). In this context,
several parameterizations have already been proposed for $w_{x}$ as a function
of the scale factor Chevallier and Polarski (2001); Linder (2003); Wetterich
(2004); Barboza and Alcaniz (2008).

On the other hand, the ratio between DM and DE energy densities
$r\left(a\right)$ is associated to the existence (or not) of an interaction
between the dark components. In order to understand how it is associated to a
dark sector interaction, it is convenient to introduce its derivative with
respect to cosmic time,

$\dot{r}=r\left(\frac{\dot{\rho}_{c}}{\rho_{c}}-\frac{\dot{\rho}_{x}}{\rho_{x}}\right)\,.$
(13)

Eq. (13) can be combined with the conservation equations for the dark
components in order to write a differential equation for the dynamics of
$r\left(a\right)$. If the dark sector interacts, a source function appears in
the term in the parenthesis of Eq. (13), and it will directly affect the time
evolution of the ratio between DM and DE energy density. The interacting case
will be discussed in more detail in Sec. III.

For the $\Lambda$CDM model, the dark components are assumed to be independent,
and the DE is characterized by the cosmological constant ($w_{x}=-1$) so that
the background energy conservation for the dark components are given by

$\displaystyle\dot{\rho}_{c}+3H\rho_{c}=0\,,$ (14)
$\displaystyle\dot{\rho}_{x}=0\,.$ (15)

Substituting Eqs. (14) and (15) in Eq. (13), one find,

$\dot{r}+3Hr=0\qquad\Rightarrow\qquad r=r_{0}\,a^{-3}\,,$ (16)

with the Hubble rate written as

$H^{2}=\Omega_{r0}\,a^{-4}+\Omega_{b0}\,a^{-3}+\Omega_{d0}\frac{1+r_{0}\,a^{-3}}{1+r_{0}}\,,$
(17)

where $r_{0}=\rho_{c0}/\rho_{x0}=\Omega_{c0}/\Omega_{x0}$. If $w_{x}=-1$, Eq.
(16) means that any deviation from $r\propto a^{-3}$ indicates the existence
of a dark sector interaction. In general, interacting models are
phenomenologically proposed by an ansatz for a source function in the energy
conservation equation, but they can also be equivalently proposed by an ansatz
for the function $r\left(a\right)$ Fuño  _et al._ (2016); vom Marttens  _et
al._ (2017b). In Ref. von Marttens  _et al._ (2019b), Eq. (17) was used to
introduce a new null test sensitive to the existence of interaction in the
dark sector.

##  III Interacting dark sector

As mentioned earlier, cosmologies with energy exchange between the dark
sector’s components constitute a viable alternative to the standard model. The
observational viability of specific classes of interacting models has been
investigated through the usual observational tests Amendola (2000); Barrow and
Clifton (2006); Clemson  _et al._ (2012); Yang and Xu (2014); Pan  _et al._
(2015); Caprini and Tamanini (2016); Yang  _et al._ (2018c, d, 2019b, 2019c),
as well as through model-independent analyses Zhou  _et al._ (2019).

At the background level, being DM and DE described by perfect fluids, this
non-gravitational coupling between dark components can be characterized by a
scalar source term $Q$, which is the time component of the covariant
derivative of the energy-momentum tensor $T^{\mu\nu}_{c\ ;\mu}=-T^{\mu\nu}_{x\
;\mu}=\left(Q,\vec{0}\right)$ or, equivalently,

$\displaystyle\dot{\rho}_{c}+3H\rho_{c}=Q\,,$ (18a)
$\displaystyle\dot{\rho}_{x}+3H\rho_{x}\left(1+w_{x}\right)=-Q\,.$ (18b)

Clearly, from the above equations, the direction of the energy transfer
depends on the sing of the source term: if $Q$ is positive, one finds DE
decaying into CDM whereas the opposite occurs if $Q$ is negative.

In this paper, we assume that the interaction term has the form $Q=3H\gamma
R\left(\rho_{c},\rho_{x}\right)$, where $\gamma$ is a free constant parameter
and $R$ is a general function that depends on the energy densities of the
components involved in the interaction. The sign of the parameter $\gamma$
gives the direction of the interaction, while its absolute value gives the
strength. Note that, irrespective of the function
$R\left(\rho_{c},\rho_{x}\right)$, the $\Lambda$CDM limit is always recovered
for $\gamma=0$. The function $R\left(\rho_{c},\rho_{x}\right)$ has unit of
energy density.

Replacing eqs. (18a) and (18b) into eq. (13), we obtain

$\dot{r}-3Hr\left[f\left(r\right)+w_{x}\right]=0\,,$ (19)

where

$f\left(r\right)\equiv\gamma R\left(\dfrac{\rho_{c}+\rho_{x}}{\rho_{c}\
\rho_{x}}\right)\,.$ (20)

As shown in Ref. von Marttens  _et al._ (2019a), this approach for describing
dark sector interactions is particularly interesting because, for any choice
of $f\left(r\right)$, eqs. (18a) and (18b) can always be decoupled. Moreover,
several interacting DE models proposed in the literature can be recovered. For
example, the very general case
$Q=3H\gamma\rho_{c}^{\alpha}\rho_{x}^{\beta}\left(\rho_{c}+\rho_{x}\right)^{\sigma}$,
with $\alpha+\beta+\sigma=1$, corresponds to $f\left(r\right)=\gamma
r^{\alpha-1}\left(1+r\right)^{\sigma+1}$.

Models with interaction in the dark sector are generally proposed to have a
specific parameterization for the source function $Q$. However, from the
definition of the function $f\left(r\right)$, it is direct to see that an
equivalent approach can be followed starting from the choice of
$f\left(r\right)$. These two functions are related by

$Q=3Hf\left(r\right)\left(\dfrac{\rho_{c}\
\rho_{x}}{\rho_{c}+\rho_{x}}\right)\,.$ (21)

From Eq. (19), one can also show that for a given $w_{x}\left(a\right)$, the
correspondence between $Q$ and $f\left(r\right)$ can be extended to
$r\left(a\right)$, in the sense that any choice of $r\left(a\right)$ leads to
an specific solution for $f\left(r\right)$, given by

$f\left(r\right)=\dfrac{a}{3r}\dfrac{dr}{da}-w_{x}\left(a\right)\,,$ (22)

and, consequently, to a specific solution for $Q$ provided by eq. (21).

Tab. 1 shows three examples for the correspondence between $f\left(r\right)$,
$Q$ and $r\left(a\right)$ for the case $w_{x}=-1$. This one-to-one mapping is
crucial to establish that, since there is no difference between choosing
$f\left(r\right)$, $Q$ or $r\left(a\right)$ for specifying a particular
interacting DE model, a model-independent reconstruction of any of these
quantities provides information regarding whether or not such dark sector’s
interaction exists.

$f\left(r\right)$ | $Q$ | $r\left(a\right)$
---|---|---
$\gamma$ | $3H\gamma\dfrac{\rho_{c}\rho_{x}}{\rho_{c}+\rho_{x}}$ | $r_{0}a^{-3\left(1-\gamma\right)}$
$\gamma(1+r)$ | $3H\gamma\rho_{c}$ | $\dfrac{\left(\gamma+1\right)r_{0}}{a^{3\left(\gamma+1\right)}\left(\gamma+\gamma r_{0}+1\right)-\gamma r_{0}}$
$\gamma\left(1+\dfrac{1}{r}\right)$ | $3H\gamma\rho_{x}$ | $\dfrac{a^{-3\left(\gamma+1\right)}\left(-\gamma a^{3\left(\gamma+1\right)}+\gamma+\gamma r_{0}+r_{0}\right)}{\gamma+1}$
Table 1: Some specific cases used to illustrate the correspondence between
$f\left(r\right)$, $Q$ and $r\left(a\right)$, considering $w_{x}=-1$.

##  IV Reconstructing the dark sector interaction

We shall now discuss in more detail how to reconstruct a possible dark sector
interaction in a model-independent way directly from the data. To this end we
address the following points in this section. First, we discuss the formalism
employed, showing the explicit equations used to perform the model-independent
reconstruction. Second, we present the datasets that we adopt in the analysis.
Finally, we present and explore our main results.

###  IV.1 Set of relevant equations

The Friedmann equation (9) can be used to write $w_{d}$ in terms of the Hubble
rate. To do so, it is necessary to replace the solutions (3b) and (8)111From
now on, since we will use only low-$z$ data, the radiation component will be
neglected. in eq. (9), and then solving the resulting equation for $w_{d}$.
This leads to the following expression

$\small{1+w_{d}\left(z\right)=\frac{\left(1+z\right)}{3}\frac{d}{dz}\left\\{\ln\left[\frac{1}{\Omega_{d0}}\left(\frac{H^{2}}{H_{0}^{2}}-\frac{\Omega_{b0}}{\left(1+z\right)^{-3}}\right)\right]\right\\}},$
(23)

which can be reduced to

$w_{d}\left(z\right)=\dfrac{E\left[3E-2\left(1+z\right)E^{\prime}\right]}{3\Omega_{b0}\left(1+z\right)^{3}-3E^{2}}\,,$
(24)

where a prime denotes derivative with respect to redshift. From eq. (6), it is
direct to obtain the $\Lambda$CDM expression for $w_{d}\left(z\right)$, i.e.,

$w_{d}^{\Lambda{\rm
CDM}}\left(z\right)=-\dfrac{\Omega_{x0}}{\Omega_{x0}+\Omega_{c0}\left(1+z\right)^{3}}\,.$
(25)

Eq. (24) shows explicitly how $w_{d}$ can be obtained from a reconstruction of
measurements of the expansion rate. However, it is worth emphasizing that in
this scenario deviations from $\Lambda$CDM result (25) do not necessarily mean
an interacting dark sector since both $w_{x}\left(z\right)$ and
$r\left(z\right)$ have not been specified. In this case, deviations from the
$\Lambda$CDM result may be due to a dynamic DE or an interacting dark sector.
That is the essence of the so-called (background) dark degeneracy Kunz (2009);
Carneiro and Borges (2014).

In order to break such degeneracy, we must first specify the DE EoS
$w_{x}\left(a\right)$. For simplicity, in this work we choose an interacting
vacuum DE, i.e., $w_{x}=-1$. Now, combining eqs. (24) and (6) the ratio
between DM and DE energy densities is given by

$r\left(z\right)=\dfrac{-2\left(z+1\right)EE^{\prime}+6E^{2}-3\Omega_{b0}\left(z+1\right)^{3}}{E\left[3E-2\left(z+1\right)E^{\prime}\right]}\,.$
(26)

Therefore, any deviation from $r\left(z\right)\propto\left(1+z\right)^{3}$
means a dark sector’s interaction. In this sense, as discussed earlier, the
$r\left(z\right)$ reconstruction itself already has all the information on
whether or not the dark sector interacts. Combining the eqs. (6) and (12), it
is possible to write an analogous acceleration condition in terms of
$r\left(z\right)$, i.e.,

$r\left(z\right)<\dfrac{3\left(E^{2}-\Omega_{b0}a^{-3}-\Omega_{r0}a^{-4}\right)}{E^{2}-\Omega_{r0}a^{-4}}-1\,.$
(27)

An equivalent procedure can be followed to obtain the DE EoS parameter given
an expression of $r\left(z\right)$. However, note that a significant advantage
of reconstructing $r\left(z\right)$ instead of $w_{x}\left(z\right)$ is that
no assumption on the density parameters of the dark sector’s components are
required in the former case. On the other hand, when $w_{x}\left(z\right)$ is
reconstructed, the prior information on $r\left(z\right)$ contains implicitly
the value of $\Omega_{c0}$ through the relation
$r_{0}=\Omega_{c0}/\left(1-\Omega_{c0}-\Omega_{b0}\right)$ Escamilla-Rivera
and Capozziello (2019).

Finally, the last quantity we reconstruct is the interacting function
$f\left(r\right)$. Combining eqs. (22) and (26), we obtain

$f\left(r(z)\right)=\frac{2\left\\{\Omega_{b0}\left(z+1\right)^{4}E\left[\left(z+1\right)E^{\prime\prime}-2E^{\prime}\right]+\Omega_{b0}\left(z+1\right)^{5}E^{\prime
2}-\left(z+1\right)E^{3}\left[10E^{\prime}+\left(z+1\right)E^{\prime\prime}\right]+3\left(z+1\right)^{2}E^{2}E^{\prime
2}+9E^{4}\right\\}}{E\left[3E-2\left(z+1\right)E^{\prime}\right]\left[-2\left(z+1\right)EE^{\prime}+6E^{2}-3\Omega_{b0}\left(z+1\right)^{3}\right]}\,.$
(28)

Eqs. (24), (26) and (28) constitute the set of equations that will be
reconstructed in our analysis. We do not choose to reconstruct the source term
$Q$ because, according to its relation with $f\left(r\right)$, it would be
necessary to solve the background energy conservation equations for the dark
components (using the reconstructed solution for $f\left(r\right)$), and the
error propagation would make this approach impracticable.

###  IV.2 Cosmological data

To obtain a fully model-independent reconstruction it is not enough to use
non-parametric statistical methods, but it is also necessary to use model-
independent data. For this reason, we consider the following cosmological data
divided into Hubble expansion rate and SN Ia luminosity distance measurements.

####  IV.2.1 $H(z)$ measurements

We construct an $H(z)$ compilation with independent measurements provided by
the following techniques:

  * •

Hubble constant, $H_{0}$: we adopt the cosmology-independent determination
$H_{0}=(75.35\pm 1.68$) km/s/Mpc, which was obtained from the latest SH0ES
analysis using the SN Ia distance-redshift relation calibrated via Cepheid
variables Camarena and Marra (2020).

  * •

Cosmic Chronometers (CC): it is possible to determine the Hubble rate of the
Universe by computing the age difference $\Delta t$ between passively-evolving
galaxies at close redshifts. The main requirements of the considered galaxy
samples are: they have similar metallicities, low star formation rates and the
average age of their stars far exceeds $\Delta t$. In the differential age
method, the derivative of the cosmic time with respect to the redshift
($dt/dz$) is approximated by the ratio of the variation of galaxies’ age with
redshift ($\Delta t/\Delta z$). This correspondence is plausible by assuming
that the analyzed galaxies were formed at the same time in the past. CC data
are estimated without assuming any cosmological model and currently the 31
available datapoints cover a wide redshift interval ($z\in[0.07,1.965]$) Zhang
_et al._ (2014); Simon  _et al._ (2005); Moresco  _et al._ (2016); Stern  _et
al._ (2010); Moresco (2015); Ratsimbazafy  _et al._ (2017)222All data points
are presented in Tab. I of the Ref. Marra and Sapone (2018)..

  * •

Baryon Acoustic Oscillation (BAO): we also consider the $H(z)$ estimates from
the anisotropic BAO signal detected in the Luminous Red Galaxy (LRG)
clustering and the quasar Ly-$\alpha$ forest. This technique does not measure
directly $H(z)$, but the combination of the Hubble rate and the scale of the
sound horizon at the drag epoch,
$d_{H}/r_{d}\equiv\left[H\left(z\right)r_{d}\right]^{-1}$. In that sense, it
is necessary to calibrate all the BAO data because of their dependence on the
sound horizon scale. Thus, in order to convert the combination
$d_{H}\left(z\right)/r_{d}$ into a model-independent measurement of the Hubble
parameter as presented in Tab. 2, we made use of the model-independent result
for the sound horizon obtained from low-redshift standard rulers
$r_{d}=\left(101.2\pm 1.8\right)h^{-1}$ Mpc Verde  _et al._ (2017)333We use
here the value obtained in Verde  _et al._ (2017) that considers the condition
of a spatially flat Universe.. Combining, this result with the aforementioned
local measurement of $H_{0}$ Camarena and Marra (2020), the sound horizon is
given by

$r_{d}=(134.0\pm 4.0)\ \rm{Mpc}\,.$ (29)

Furthermore, to extract the 3D BAO feature from the large-scale matter
distribution of the universe, it is necessary to use a fiducial cosmological
model and one of the most important discussions about BAO data concerns its
dependence on this fiducial model. Recent works have shown that the
constraints from BAO estimates are model-independent for a wide class of
cosmological models Carter  _et al._ (2020); Bernal  _et al._ (2020).

In this work, in order to avoid double counting with the BAO data used to
obtain the sound horizon, we do not use the first two data points, at
$z=[0.38,0.51]$, provided in the full analysis of the galaxies distribution in
the SDSS DR12 LRG Alam  _et al._ (2017).

From the eBOSS DR16 LRG catalog Alam  _et al._ (2020), we use one data point
at $z=0.698$. We also use the $d_{H}\left(z\right)/r_{d}$ measurement from the
eBOSS DR16 QSO catalog at $z=1.48$. Finally, the analyses of the cross-
correlation and the auto-correlation of the quasar Ly-$\alpha$ absorption from
eBOSS DR16, that provide two uncorrelated estimates of
$d_{H}\left(z\right)/r_{d}$ rate at high-$z$444We use here the Gaussian
approximation of the eBOSS DR16 Ly-$\alpha$ data (auto and cross), which is a
conservative approximation..

BAO data
---
Catalog | $z$ | $H\left(z\right)$ | $\sigma_{H\left(z\right)}$ | Ref.
eBOSS DR16 LRG | $0.698$ | $113.0$ | $4.16$ | Alam  _et al._ (2020)
eBOSS DR16 QSO | $1.48$ | $170.0$ | $7.73$
eBOSS DR16 Ly-$\alpha$ (cross) | $2.33$ | $251.5$ | $10.63$
eBOSS DR16 Ly-$\alpha$ (auto) | $2.33$ | $247.0$ | $11.64$
Table 2: BAO data points.

This calibration using one data of our set produces correlations between
BAO-$H(z)$ measurements and $H_{0}$. We calculate the total correlation matrix
of the $H(z)$ dataset by considering the BAO-$H(z)$, the sound horizon (29)
and $H_{0}$ values and uncertainties. The final matrix correlation of $H_{0}$
and BAO-$H(z)$ data is given by

$\small\mathcal{C}_{H}=\begin{bmatrix}2.83&4.25&6.37&9.38&9.25\\\
4.25&17.32&15.59&23.01&22.66\\\ 6.37&15.59&59.81&34.53&33.96\\\
9.38&23.01&34.53&112.95&50.24\\\ 9.25&22.66&33.96&50.24&135.44\\\
\end{bmatrix},$ (30)

where the first row and column correspond to the $H_{0}$ correlation
coefficients. Note that the CC data is uncorrelated, thus its correspondent
block in the total covariance matrix is diagonal.

####  IV.2.2 Type Ia Supernovae

The second dataset we use is the Pantheon catalog of type Ia Supernovae (SNe
Ia), which contains 1048 data points in the range $0.01<z<2.26$, with their
respective covariance matrix (including statistical and systematic
errors)555The data as well as its covariance matrix can be found at
github.com/dscolnic/Pantheon. The data obtained from SN Ia is not directly the
Hubble rate, but the $B$-band apparent magnitude $m_{B}$, which is related to
the luminosity distance $d_{L}\left(z\right)$ by

$m_{B}=5\log_{10}\left[\dfrac{d_{L}\left(z\right)}{1{\rm
Mpc}}\right]+25+M_{B}\,,$ (31)

where $M_{B}$ is the $B$-band absolute magnitude of the SNe Ia.

We directly relate $d_{L}$ with the comoving distance ($d_{c}$) by considering
a spatially flat universe via:

$d_{c}(z)=\frac{d_{L}(z)}{1+z}=\frac{10^{(m_{B}-M_{B}-25)/5}}{1+z}.$ (32)

Thus, to convert this $m_{B}$ SNe Ia catalog into comoving-distance data it is
just necessary to know the absolute magnitude value. In general, for parameter
selection analyses, $M_{B}$ is considered as a nuisance parameter. Here,
however, we fix it at the value Camarena and Marra (2020)

$M=-19.2334\pm 0.0404\,{\rm mag}\,.$ (33)

The previous approach allows us to directly calculate the Hubble rate function
with the reconstruction of the first derivative of $d_{c}$ through its
definition,

${d_{c}}\left(z\right)=\int_{0}^{z}\dfrac{d\tilde{z}}{H\left(\tilde{z}\right)}\quad\Rightarrow\quad
H\left(z\right)=\dfrac{1}{d_{c}^{\prime}\left(z\right)}\,.$ (34)

It is worth to highlight that although the $M_{B}$ considered is also model
independent, this value is not relevant in our analysis because all the
interaction quantities depend only on the normalized Hubble parameter
($E\equiv H(z)/H_{0}=d_{c}^{\prime}(0)/d_{c}^{\prime}(z)$) and therefore the
multiplicative factor $10^{-M/5}$ in Eq. (32) cancels. In this way, the SNe
data becomes almost model independent.

Another important aspect in our SN Ia analysis is that the SN Ia data
distribution along the redshift is not uniform. In fact, there are many more
SNe Ia at low-$z$ than at high-$z$, as shown in Fig. 1. For this reason, we
show only the data set within the redshift range used in our analysis. The
redshift threshold for the SN Ia reconstruction will be discussed on a case-
by-case basis in Sec. IV.3.

Figure 1: Distribution of the SN Ia data points in terms of the redshift. The
dashed lines indicate $z=0.8$ and $z=1.5$, dividing the redshift range into
three intervals. The first interval contains $92.18\%$ of the Pantheon’s
catalog data (966 SNe Ia), the second interval includes $7.25\%$ (76 SNe Ia),
remaining only $0.57\%$ (6 SNe Ia) in the last interval.

Figure 2: Model-independent reconstruction of the EoS parameter of the unified
dark fluid. The solid line is the best fit of the GP reconstruction, the
dashed line is the $\Lambda$CDM result, and the dotted line is the
acceleration condition, given by Eq. (12). Left panel: Result obtained using
the $H(z)$ data. Right panel: Result obtained using the SN Ia data.

####  IV.2.3 Reconstruction

The model-independent reconstruction is performed by applying the Gaussian
Process (GP) method. We use the GaPP python library Seikel  _et al._ (2012)
with a square exponential covariance function and optimize its hyperparameters
by maximizing the GP’s likelihood to obtain the reconstruction of $H(z)$,
$d_{c}(z)$ and their derivatives. As mentioned earlier, the quantities that we
will reconstruct are $w_{d}(z)$, $r(z)$ and $f(z)$, which are respectively
given by Eqs. (24), (26) and (28). For the first data set, we extract the
Hubble rate directly from the data. The normalization of the Hubble rate is
made considering the value of $H\left(z=0\right)=H_{0}$ obtained with the
reconstruction. On the other hand, for the second data set, the reconstructed
quantity is initially the $d_{c}$, and the (normalized) Hubble rate is
obtained following Eq. (34). Note that, compared with the first data set, the
SN Ia analysis will always need one more derivative, which may affect the
error propagation.

In order to reconstruct $w_{d}(z)$, $r(z)$ and $f(z)$ (Eqs. (24), (26) and
(28)), we perform a Monte Carlo sampling taking into account the mean values
and the complete covariance matrix of $H_{0}$, $H(z)$, $H^{\prime}(z)$ and
$H^{\prime\prime}(z)$ jointly described by a GP (multivariate Gaussian
distribution), for the first data set. In the case of SN data, we perform the
Monte Carlo sampling with $d^{\prime}_{c}(0)$, $d^{\prime}_{c}(z)$,
$d^{\prime\prime}_{c}(z)$ and $d^{\prime\prime\prime}_{c}(z)$ mean values and
their covariance matrix (see Eq. (2.6) of Ref. Seikel  _et al._ (2012)).
Another prior information required by the Monte Carlo sampling of eqs. (24),
(26) and (28) is the baryon density parameter $\Omega_{b0}$. Here, we use a
result obtained in Ref. Pettini and Cooke (2012),
$\omega_{b}\equiv\Omega_{b0}h^{2}=0.0223\pm 0.0009$. When combined with the
local measurement of $H_{0}$ presented in Sec. IV.2.1, we obtain
$\Omega_{b0}=0.0393\pm 0.0024$.

Figure 3: Model-independent reconstruction of the ratio between CDM and DE
energy densities. The solid line is the best fit of the GP reconstruction, the
dashed line is the $\Lambda$CDM result, and the dotted line is the
acceleration condition, given by Eq. (27). Left panel: Result obtained using
the $H(z)$ data. Right panel: Result obtained using the SN Ia data.

###  IV.3 Results

In this section we present and discuss our results for the reconstruction of
the possible interaction in the cosmological dark sector. The first quantity
reconstructed is the unified dark EoS parameter, given by Eq. (24). As already
mentioned, deviations from $w_{d}(z)=w_{d}^{\Lambda{\rm CDM}}\left(z\right)$
do not necessarily mean a dark sector interaction, but can also represent a
dynamical DE or even a combination of both.

Figure 4: Model-independent reconstruction of the dark sector interaction
function in terms of the redshift. The dashed horizontal line in zero,
indicates the $\Lambda$CDM (non-interacting) result. Left panel: Result
obtained using the $H(z)$ data. Right panel: Result obtained using the SN Ia
data.

Fig. 2 shows the $w_{d}\left(z\right)$ reconstruction. Whereas its left panel
shows the results obtained using $H(z)$ measurements, the right panel shows
the results from SN Ia data. In both cases, one can see that $\Lambda$CDM
model666From now on, when $\Lambda$CDM model is compared to our results, it is
implicit that we refer to the best fit of the Planck analysis using data from
TT,TE,EE+lowE+lensing+BAO Aghanim  _et al._ (2020). is consistent at least at
$3\sigma$ CL over the entire redshift range. Furthermore, in both cases, the
best fits of the reconstructions indicate that the universe switched from a
decelerated phase to an accelerated one at around $z\approx 0.6$. We estimate
the CL of $w_{d}$ at $z=0$ up to $5\sigma$ via Monte Carlo sampling for the
two data sets and confirm the present cosmic acceleration with a higher CL.
Our results indicate current acceleration at $\sim 5.9\sigma$ CL with the
$H(z)$ data and, if we assume the gaussianity is maintained at higher levels,
the SNe Ia data confirms cosmic acceleration at $\sim 12\sigma$ CL. It is
worth emphasizing that these results only assume a dark sector composed of a
pressureless matter (whether interacting or not) combined with a general DE
component, which is utterly free from a hypothesis about its nature. In the
SNe Ia analysis, even though all data points have been used, we restrict
ourselves to show only the result for $z<1.5$. For higher values of the
redshift, the SNe Ia analysis can not properly constraint $w_{d}(z)$, and it
has no physical meaning.

Imposing now $w=-1$, we reconstruct the ratio between DM and DE energy
densities, given by Eq. (26). Now any deviation from
$r\left(z\right)\propto\left(1+z\right)^{3}$ indicates the existence of a non-
gravitational interaction between the dark components. The results for the
$H\left(z\right)$ and SN Ia data sets are presented respectively in the left
and right panels of Fig. 3. As well as the first case, both analyses have an
excellent agreement with the $\Lambda$CDM model, so that the reconstruction is
consistent in $1\sigma$ CL with the standard cosmological model across the
whole redshift interval analyzed. Once again, as the first analysis, the best
fits of the $r\left(z\right)$ reconstructions predict a transition from a
decelerated phase to an accelerated phase at $z\approx 0.6$. Using Eq. (27),
the cosmic acceleration can be assessed in terms of $r(z)$. The evidence that
universe now ($z=0$) experiences an accelerated expansion remains at more than
$6\sigma$ for both data sets. Even though this analysis considers only
deviation associated with dark sector interactions, similar conclusions in
comparison to the $w_{d}(z)$ analysis are indeed expected because of the dark
degeneracy. Here, the SN Ia analysis only delivers reasonable constraints for
$z<0.8$.

Figure 5: Parametric combination of the model-independent reconstructions for
the interacting function $f$ and the ratio between CDM and DE energy densities
$r$. In each panel, the result is compared with one of the interacting models
shown in Tab. 1 with $\gamma=\pm 0.1,\pm 0.25,\pm 0.5$ corresponding to bands
in dim gray, standard grey and light grey respectively. Top panel:
$f\left(r\right)=\gamma$. Middle panel: $f\left(r\right)=\gamma(1+r)$. Bottom
panel: $f\left(r\right)=\gamma(1+1/r)$.

The last reconstruction we perform is for the interacting function, defined in
Eq. (20), in terms of the redshift, which is given by Eq. (28). As already
presented in Tab. 1, this quantity identifies uniquely in a simple way the
source term $Q$, and for this reason, its reconstruction can be seen as a
model-independent reconstruction of the dark sector interaction. The results
using $H(z)$ and SN Ia data are shown respectively in the left and right
panels of Fig. 4. In both analysis the $\Lambda$CDM model agrees within
$2\sigma$ CL for the entire redshift range.

Finally, since the relation between the interaction function and the source
term is more evident when $f$ is written in terms of the ratio between CDM and
DE energy densities, we take advantage of the $r(z)$ and $f(z)$ results to
perform a parametric investigation of $f(r)$. In Fig. 5, we show the $1\sigma$
CL results for the interacting function in terms of the ratio between CDM and
DE energy densities, being the red diamonds the results obtained with the
$H(z)$ measurements and the black squares obtained with the SN Ia data. In
each of the three panels we compare our results with the models highlighted in
Tab. 1, where the regions between $\gamma=\pm 0.1$, $\gamma=\pm 0.25$ and
$\gamma=\pm 0.5$ are respectively associated to dim gray, standard grey and
light grey. The top, middle and bottom panels of Fig. 5 show the result for
the cases $f(r)=\gamma$, $f(r)=\gamma(1+r)$ and $f(r)=\gamma(1+1/r)$. Our
results confirm that the $\Lambda$CDM model (dashed line in $f(r)=0$) is
compatible with our results in about $1\sigma$ CL, but also indicate that
current data has no strength to discard the interacting models.

##  V Future perspectives

With the advent of the next-generation LSS surveys, the subsequent few years
promise to be particularly fruitful for observational cosmology. Among others,
we can highlight J-PAS Benitez  _et al._ (2014); Bonoli  _et al._ (2020), DESI
Aghamousa  _et al._ (2016), EUCLID Amendola  _et al._ (2013) and SKA Bacon
_et al._ (2020). These surveys will provide measurements of $H\left(z\right)$
from BAO with good precision in a redshift interval where there are not many
SNe Ia ($z\gtrsim 0.8$). In that sense, this data will play a crucial role as
a complementary data for unveiling the nature of the dark sector of the
Universe. In this work, we use the forecast for the (normalized) Hubble rate
from J-PAS, DESI, EUCLID and SKA (bands 1 and 2) presented in Benitez  _et
al._ (2014); Bacon  _et al._ (2020). The specifications of the catalogs we use
are the following:

$\displaystyle\qquad\bullet\ {\rm J-PAS:}$ $\displaystyle 0.31<$
$\displaystyle\ z<3.91$ $\displaystyle,\qquad$ $\displaystyle N=54$ (35)
$\displaystyle\qquad\bullet\ {\rm DESI:}$ $\displaystyle 0.15<$
$\displaystyle\ z<1.691$ $\displaystyle,\qquad$ $\displaystyle N=30$ (36)
$\displaystyle\qquad\bullet\ {\rm Euclid:}$ $\displaystyle 0.7<$
$\displaystyle\ z<2.0$ $\displaystyle,\qquad$ $\displaystyle N=14$ (37)
$\displaystyle\qquad\bullet\ {\rm SKA\ B1:}$ $\displaystyle 0.5<$
$\displaystyle\ z<2.9$ $\displaystyle,\qquad$ $\displaystyle N=12$ (38)
$\displaystyle\qquad\bullet\ {\rm SKA\ B2:}$ $\displaystyle 0.1<$
$\displaystyle\ z<0.4$ $\displaystyle,\qquad$ $\displaystyle N=4$ (39)

where $N$ is the number of points in each survey. The results of the relative
error $\sigma_{E/E}$ are shown in Fig. 6. As can be seen, the surveys will be
able to provide a very important complementary data to investigate the DE
component and the late-time cosmic expansion. In particular, we find that
Euclid and SKA (band 1) will deliver sub-percent measurements of
$H\left(z\right)$ at $0.5<z<2$. We also confirm the very good performance of
J-PAS at low-$z$ ($z\lesssim 0.7$), as reported in Aparicio Resco  _et al._
(2020).

Our forecast analysis is similar to one presented in Bengaly  _et al._
(2020a); Bengaly (2019); Bengaly  _et al._ (2020b). Using the data presented
in Fig. 6, we build a $\Lambda$CDM ($\Omega_{m}=0.3$777Note that, since only
the normalized Hubble parameter is used, no assumption for $H_{0}$ is
required.) mock catalog for each survey, and use GP to reconstruct them.
First, we apply the reconstructions of the (normalized) Hubble rate in Eq.
(24). In this case, we forecast the EoS parameter for the dark sector as a
whole, i.e., no assumption on the time-dependence of the DE component, or
whether it interacts or not with CDM is made. The relative error on
$\sigma_{w_{d}\left(z\right)}/w_{d}\left(z\right)$ is presented in the top
panel of Fig. 7. Finally, assuming that DE is described by $w_{x}=-1$, but
admitting the possibility of an interaction in the dark sector, we apply the
reconstructed (normalized) Hubble rate in Eq. (26). In this case, we forecast
the ratio between CDM and DE energy densities, which, for a given
$w_{x}\left(a\right)$, determines uniquely an interaction between the dark
components. In that sense, this analysis shows the strength of the next-
generation LSS surveys for constraining interacting models in a model-
independent way. The result for the relative error on
$\sigma_{r\left(z\right)}/r\left(z\right)$ is presented in the bottom panel of
Fig. 7. In both forecast analyses, we highlight that Euclid and SKA (band 1)
will be able to constraint $w_{d}\left(z\right)$ and $r\left(z\right)$ with
$\sim 4\%$ of precision at $z\approx 1$.

Figure 6: Forecasted relative error on $\sigma_{E/E}$.

Figure 7: Results of the forecast analysis. Top panel: Forecasted relative
error on $\sigma_{w_{d}\left(z\right)}/w_{d}\left(z\right)$. Bottom panel:
Forecasted relative error on $\sigma_{r\left(z\right)}/r\left(z\right)$.

##  VI Conclusions

A number of recent analysis have shown that models with interaction in the
dark sector are able to provide a good description of the current
observational data, which makes them a promising alternative to the standard
cosmology. However, the absence of guidance from a fundamental theory on the
coupling term $Q$ makes its choice arbitrary.

In this paper we tested the possibility of a non-minimal coupling between the
CDM and DE components by reconstructing the functions $w_{d}(z)$, $r(z)$ and
$f(r)$, as defined in Secs. II, and III. Using currently available data of
type Ia Supernova, Cosmic Chronometers and Baryonic Acoustic Oscillations, we
followed the model-independent approach discussed in Sec. IV and showed that
the current observations show a good agreement (within 3$\sigma$ level) with
the $\Lambda$CDM hypothesis of uncoupled dark components ($Q=0$), although an
interaction between them cannot be completely ruled out. For completeness, we
also performed a forecast analysis for the next-generation of large-scale
structure surveys considering their predictions for the normalized expansion
rate $E(z)$. We found that J-PAS will have a better performance at low-$z$
($\lesssim 0.7$) when compared with DESI and that Euclid and SKA (band 1) will
be able to constrain the functions $w_{d}\left(z\right)$ and $r\left(z\right)$
with $\sim 4\%$ of precision at $z\approx 1$. Such results clearly show that
the upcoming data from the next-generation surveys will play a crucial role in
probing the possibility of a non-minimal coupling between the components of
the cosmological dark sector.

## Acknowledgements

It is a pleasure to thank Carlos A. P. Bengaly for sharing the data used in
the forecast analysis. RvM acknowledges support from the Programa de
Capacitação Institucional PCI/ON/MCTI. JA acknowledges support from CNPq
(Grants no. 310790/2014-0 and 400471/2014-0) and FAPERJ (grant no. 233906).
JEG, VM and LC thanks CNPq (Brazil) for partial financial support. This
project has received funding from the European Union’s Horizon 2020 research
and innovation programme under the Marie Skłodowska-Curie grant agreement No
888258.

## References

  * Aghanim  _et al._ (2020) N. Aghanim _et al._ (Planck), Astron. Astrophys. 641, A6 (2020), arXiv:1807.06209 [astro-ph.CO] .
  * Abbott  _et al._ (2019a) T. Abbott _et al._ (DES), Phys. Rev. Lett. 122, 171301 (2019a), arXiv:1811.02375 [astro-ph.CO] .
  * Abbott  _et al._ (2019b) T. Abbott _et al._ (DES), Phys. Rev. D 99, 123505 (2019b), arXiv:1810.02499 [astro-ph.CO] .
  * Verde  _et al._ (2019) L. Verde, T. Treu, and A. Riess (2019) arXiv:1907.10625 [astro-ph.CO] .
  * Riess (2019) A. G. Riess, Nature Rev. Phys. 2, 10 (2019), arXiv:2001.03624 [astro-ph.CO] .
  * Asgari  _et al._ (2020) M. Asgari _et al._ (KiDS) (2020) arXiv:2007.15633 [astro-ph.CO] .
  * Zimdahl and Pavon (2001) W. Zimdahl and D. Pavon, Phys. Lett. B 521, 133 (2001), arXiv:astro-ph/0105479 .
  * Alcaniz and Maia (2003) J. Alcaniz and J. Maia, Phys. Rev. D 67, 043502 (2003), arXiv:astro-ph/0212510 .
  * von Marttens  _et al._ (2019a) R. von Marttens, L. Casarini, D. Mota, and W. Zimdahl, Phys. Dark Univ. 23, 100248 (2019a), arXiv:1807.11380 [astro-ph.CO] .
  * Benetti  _et al._ (2019) M. Benetti, W. Miranda, H. A. Borges, C. Pigozzo, S. Carneiro, and J. S. Alcaniz, JCAP 12, 023 (2019), arXiv:1908.07213 [astro-ph.CO] .
  * Xia and Wang (2016) D.-M. Xia and S. Wang, Mon. Not. Roy. Astron. Soc. 463, 952 (2016), arXiv:1608.04545 [astro-ph.CO] .
  * Kumar and Nunes (2017) S. Kumar and R. C. Nunes, Phys. Rev. D 96, 103511 (2017), arXiv:1702.02143 [astro-ph.CO] .
  * Yang  _et al._ (2019a) W. Yang, O. Mena, S. Pan, and E. Di Valentino, Phys. Rev. D 100, 083509 (2019a), arXiv:1906.11697 [astro-ph.CO] .
  * Martinelli  _et al._ (2019) M. Martinelli, N. B. Hogg, S. Peirone, M. Bruni, and D. Wands, Mon. Not. Roy. Astron. Soc. 488, 3423 (2019), arXiv:1902.10694 [astro-ph.CO] .
  * Cid  _et al._ (2019) A. Cid, B. Santos, C. Pigozzo, T. Ferreira, and J. Alcaniz, JCAP 03, 030 (2019), arXiv:1805.02107 [astro-ph.CO] .
  * Yang  _et al._ (2017) W. Yang, S. Pan, and D. F. Mota, Phys. Rev. D 96, 123508 (2017), arXiv:1709.00006 [astro-ph.CO] .
  * vom Marttens  _et al._ (2017a) R. vom Marttens, L. Casarini, W. Zimdahl, W. Hipólito-Ricaldi, and D. Mota, Phys. Dark Univ. 15, 114 (2017a), arXiv:1702.00651 [astro-ph.CO] .
  * vom Marttens  _et al._ (2014) R. vom Marttens, W. Hipólito-Ricaldi, and W. Zimdahl, JCAP 08, 004 (2014), arXiv:1403.0427 [astro-ph.CO] .
  * da Silva  _et al._ (2020) W. da Silva, J. Gonzalez, R. Silva, and J. Alcaniz, (2020), arXiv:2009.10216 [astro-ph.CO] .
  * Gonzalez  _et al._ (2018) J. Gonzalez, H. Silva, R. Silva, and J. Alcaniz, Eur. Phys. J. C 78, 730 (2018), arXiv:1809.00439 [astro-ph.CO] .
  * Di Valentino  _et al._ (2017) E. Di Valentino, A. Melchiorri, and O. Mena, Phys. Rev. D 96, 043503 (2017), arXiv:1704.08342 [astro-ph.CO] .
  * Yang  _et al._ (2018a) W. Yang, A. Mukherjee, E. Di Valentino, and S. Pan, Phys. Rev. D 98, 123527 (2018a), arXiv:1809.06883 [astro-ph.CO] .
  * Kumar  _et al._ (2019) S. Kumar, R. C. Nunes, and S. K. Yadav, Eur. Phys. J. C 79, 576 (2019), arXiv:1903.04865 [astro-ph.CO] .
  * Pan  _et al._ (2019a) S. Pan, W. Yang, C. Singha, and E. N. Saridakis, Phys. Rev. D 100, 083539 (2019a), arXiv:1903.10969 [astro-ph.CO] .
  * Di Valentino  _et al._ (2020) E. Di Valentino, A. Melchiorri, O. Mena, and S. Vagnozzi, Phys. Rev. D 101, 063502 (2020), arXiv:1910.09853 [astro-ph.CO] .
  * Yang  _et al._ (2018b) W. Yang, S. Pan, E. Di Valentino, R. C. Nunes, S. Vagnozzi, and D. F. Mota, JCAP 09, 019 (2018b), arXiv:1805.08252 [astro-ph.CO] .
  * Pan  _et al._ (2019b) S. Pan, W. Yang, E. Di Valentino, E. N. Saridakis, and S. Chakraborty, Phys. Rev. D 100, 103520 (2019b), arXiv:1907.07540 [astro-ph.CO] .
  * Pan  _et al._ (2020) S. Pan, W. Yang, and A. Paliathanasis, Mon. Not. Roy. Astron. Soc. 493, 3114 (2020), arXiv:2002.03408 [astro-ph.CO] .
  * Valiviita  _et al._ (2008) J. Valiviita, E. Majerotto, and R. Maartens, JCAP 07, 020 (2008), arXiv:0804.0232 [astro-ph] .
  * Costa and Alcaniz (2010) F. Costa and J. Alcaniz, Phys. Rev. D 81, 043506 (2010), arXiv:0908.4251 [astro-ph.CO] .
  * Cai and Su (2010) R.-G. Cai and Q. Su, Phys. Rev. D 81, 103514 (2010), arXiv:0912.1943 [astro-ph.CO] .
  * Carneiro  _et al._ (2019) S. Carneiro, H. Borges, R. von Marttens, J. Alcaniz, and W. Zimdahl, (2019), arXiv:1909.10336 [gr-qc] .
  * Zhao  _et al._ (2017) G.-B. Zhao _et al._ , Nature Astron. 1, 627 (2017), arXiv:1701.08165 [astro-ph.CO] .
  * Joudaki  _et al._ (2018) S. Joudaki, M. Kaplinghat, R. Keeley, and D. Kirkby, Phys. Rev. D 97, 123501 (2018), arXiv:1710.04236 [astro-ph.CO] .
  * L’Huillier  _et al._ (2018) B. L’Huillier, A. Shafieloo, and H. Kim, Mon. Not. Roy. Astron. Soc. 476, 3263 (2018), arXiv:1712.04865 [astro-ph.CO] .
  * Shafieloo  _et al._ (2018) A. Shafieloo, B. L’Huillier, and A. A. Starobinsky, Phys. Rev. D 98, 083526 (2018), arXiv:1804.04320 [astro-ph.CO] .
  * L’Huillier  _et al._ (2019) B. L’Huillier, A. Shafieloo, E. V. Linder, and A. G. Kim, Mon. Not. Roy. Astron. Soc. 485, 2783 (2019), arXiv:1812.03623 [astro-ph.CO] .
  * Keeley  _et al._ (2019) R. E. Keeley, S. Joudaki, M. Kaplinghat, and D. Kirkby, JCAP 12, 035 (2019), arXiv:1905.10198 [astro-ph.CO] .
  * Keeley  _et al._ (2020) R. E. Keeley, A. Shafieloo, B. L’Huillier, and E. V. Linder, Mon. Not. Roy. Astron. Soc. 491, 3983 (2020), arXiv:1905.10216 [astro-ph.CO] .
  * Liao  _et al._ (2020) K. Liao, A. Shafieloo, R. E. Keeley, and E. V. Linder, Astrophys. J. Lett. 895, L29 (2020), arXiv:2002.10605 [astro-ph.CO] .
  * Keeley  _et al._ (2021) R. E. Keeley, A. Shafieloo, G.-B. Zhao, J. A. Vazquez, and H. Koo, Astron. J. 161, 151 (2021), arXiv:2010.03234 [astro-ph.CO] .
  * Seikel  _et al._ (2012) M. Seikel, C. Clarkson, and M. Smith, JCAP 06, 036 (2012), arXiv:1204.2832 [astro-ph.CO] .
  * von Marttens  _et al._ (2019b) R. von Marttens, V. Marra, L. Casarini, J. E. Gonzalez, and J. Alcaniz, Phys. Rev. D99, 043521 (2019b), arXiv:1812.02333 [astro-ph.CO] .
  * Benitez  _et al._ (2014) N. Benitez _et al._ (J-PAS), (2014), arXiv:1403.5237 [astro-ph.CO] .
  * Bonoli  _et al._ (2020) S. Bonoli _et al._ , (2020), arXiv:2007.01910 [astro-ph.CO] .
  * Aghamousa  _et al._ (2016) A. Aghamousa _et al._ (DESI), (2016), arXiv:1611.00036 [astro-ph.IM] .
  * Amendola  _et al._ (2013) L. Amendola _et al._ (Euclid Theory Working Group), Living Rev. Rel. 16, 6 (2013), arXiv:1206.1225 [astro-ph.CO] .
  * Bacon  _et al._ (2020) D. J. Bacon _et al._ (SKA), Publ. Astron. Soc. Austral. 37, e007 (2020), arXiv:1811.02743 [astro-ph.CO] .
  * von Marttens  _et al._ (2020) R. von Marttens, L. Lombriser, M. Kunz, V. Marra, L. Casarini, and J. Alcaniz, Phys. Dark Univ. 28, 100490 (2020), arXiv:1911.02618 [astro-ph.CO] .
  * Ratra and Peebles (1988) B. Ratra and P. J. E. Peebles, Phys. Rev. D37, 3406 (1988).
  * Chevallier and Polarski (2001) M. Chevallier and D. Polarski, Int. J. Mod. Phys. D 10, 213 (2001), arXiv:gr-qc/0009008 .
  * Linder (2003) E. V. Linder, Phys. Rev. Lett. 90, 091301 (2003), arXiv:astro-ph/0208512 .
  * Wetterich (2004) C. Wetterich, Phys. Lett. B 594, 17 (2004), arXiv:astro-ph/0403289 .
  * Barboza and Alcaniz (2008) J. Barboza, E.M. and J. Alcaniz, Phys. Lett. B 666, 415 (2008), arXiv:0805.1713 [astro-ph] .
  * Fuño  _et al._ (2016) A. R. Fuño, W. Hipólito-Ricaldi, and W. Zimdahl, Mon. Not. Roy. Astron. Soc. 457, 2958 (2016), arXiv:1409.7706 [astro-ph.CO] .
  * vom Marttens  _et al._ (2017b) R. F. vom Marttens, L. Casarini, W. Hipólito-Ricaldi, and W. Zimdahl, JCAP 01, 050 (2017b), arXiv:1610.01665 [astro-ph.CO] .
  * Amendola (2000) L. Amendola, Phys. Rev. D62, 043511 (2000), arXiv:astro-ph/9908023 [astro-ph] .
  * Barrow and Clifton (2006) J. D. Barrow and T. Clifton, Phys. Rev. D73, 103520 (2006), arXiv:gr-qc/0604063 [gr-qc] .
  * Clemson  _et al._ (2012) T. Clemson, K. Koyama, G.-B. Zhao, R. Maartens, and J. Valiviita, Phys. Rev. D85, 043007 (2012), arXiv:1109.6234 [astro-ph.CO] .
  * Yang and Xu (2014) W. Yang and L. Xu, Phys. Rev. D89, 083517 (2014), arXiv:1401.1286 [astro-ph.CO] .
  * Pan  _et al._ (2015) S. Pan, S. Bhattacharya, and S. Chakraborty, Mon. Not. Roy. Astron. Soc. 452, 3038 (2015), arXiv:1210.0396 [gr-qc] .
  * Caprini and Tamanini (2016) C. Caprini and N. Tamanini, JCAP 1610, 006 (2016), arXiv:1607.08755 [astro-ph.CO] .
  * Yang  _et al._ (2018c) W. Yang, S. Pan, and J. D. Barrow, Phys. Rev. D97, 043529 (2018c), arXiv:1706.04953 [astro-ph.CO] .
  * Yang  _et al._ (2018d) W. Yang, S. Pan, R. Herrera, and S. Chakraborty, Phys. Rev. D98, 043517 (2018d), arXiv:1808.01669 [gr-qc] .
  * Yang  _et al._ (2019b) W. Yang, S. Pan, L. Xu, and D. F. Mota, Mon. Not. Roy. Astron. Soc. 482, 1858 (2019b), arXiv:1804.08455 [astro-ph.CO] .
  * Yang  _et al._ (2019c) W. Yang, S. Vagnozzi, E. Di Valentino, R. C. Nunes, S. Pan, and D. F. Mota, JCAP 1907, 037 (2019c), arXiv:1905.08286 [astro-ph.CO] .
  * Zhou  _et al._ (2019) Z. Zhou, T. Zhang, and T. Li, Eur. Phys. J. C 79, 527 (2019), arXiv:1908.06254 [astro-ph.CO] .
  * Kunz (2009) M. Kunz, Phys. Rev. D 80, 123001 (2009), arXiv:astro-ph/0702615 .
  * Carneiro and Borges (2014) S. Carneiro and H. Borges, JCAP 06, 010 (2014), arXiv:1402.2316 [astro-ph.CO] .
  * Escamilla-Rivera and Capozziello (2019) C. Escamilla-Rivera and S. Capozziello, Int. J. Mod. Phys. D 28, 1950154 (2019), arXiv:1905.04602 [gr-qc] .
  * Camarena and Marra (2020) D. Camarena and V. Marra, Phys. Rev. Res. 2, 013028 (2020), arXiv:1906.11814 [astro-ph.CO] .
  * Zhang  _et al._ (2014) C. Zhang, H. Zhang, S. Yuan, T.-J. Zhang, and Y.-C. Sun, Res. Astron. Astrophys. 14, 1221 (2014), arXiv:1207.4541 [astro-ph.CO] .
  * Simon  _et al._ (2005) J. Simon, L. Verde, and R. Jimenez, Phys. Rev. D 71, 123001 (2005), arXiv:astro-ph/0412269 .
  * Moresco  _et al._ (2016) M. Moresco, L. Pozzetti, A. Cimatti, R. Jimenez, C. Maraston, L. Verde, D. Thomas, A. Citro, R. Tojeiro, and D. Wilkinson, JCAP 05, 014 (2016), arXiv:1601.01701 [astro-ph.CO] .
  * Stern  _et al._ (2010) D. Stern, R. Jimenez, L. Verde, M. Kamionkowski, and S. Stanford, JCAP 02, 008 (2010), arXiv:0907.3149 [astro-ph.CO] .
  * Moresco (2015) M. Moresco, Mon. Not. Roy. Astron. Soc. 450, L16 (2015), arXiv:1503.01116 [astro-ph.CO] .
  * Ratsimbazafy  _et al._ (2017) A. Ratsimbazafy, S. Loubser, S. Crawford, C. Cress, B. Bassett, R. Nichol, and P. Väisänen, Mon. Not. Roy. Astron. Soc. 467, 3239 (2017), arXiv:1702.00418 [astro-ph.CO] .
  * Marra and Sapone (2018) V. Marra and D. Sapone, Phys. Rev. D 97, 083510 (2018), arXiv:1712.09676 [astro-ph.CO] .
  * Verde  _et al._ (2017) L. Verde, J. L. Bernal, A. F. Heavens, and R. Jimenez, Mon. Not. Roy. Astron. Soc. 467, 731 (2017), arXiv:1607.05297 [astro-ph.CO] .
  * Carter  _et al._ (2020) P. Carter, F. Beutler, W. J. Percival, J. DeRose, R. H. Wechsler, and C. Zhao, Mon. Not. Roy. Astron. Soc. 494, 2076 (2020), arXiv:1906.03035 [astro-ph.CO] .
  * Bernal  _et al._ (2020) J. L. Bernal, T. L. Smith, K. K. Boddy, and M. Kamionkowski, (2020), arXiv:2004.07263 [astro-ph.CO] .
  * Alam  _et al._ (2017) S. Alam _et al._ (BOSS), Mon. Not. Roy. Astron. Soc. 470, 2617 (2017), arXiv:1607.03155 [astro-ph.CO] .
  * Alam  _et al._ (2020) S. Alam _et al._ (eBOSS), (2020), arXiv:2007.08991 [astro-ph.CO] .
  * Pettini and Cooke (2012) M. Pettini and R. Cooke, Mon. Not. Roy. Astron. Soc. 425, 2477 (2012), arXiv:1205.3785 [astro-ph.CO] .
  * Aparicio Resco  _et al._ (2020) M. Aparicio Resco _et al._ , Mon. Not. Roy. Astron. Soc. 493, 3616 (2020), arXiv:1910.02694 [astro-ph.CO] .
  * Bengaly  _et al._ (2020a) C. A. Bengaly, C. Clarkson, and R. Maartens, JCAP 05, 053 (2020a), arXiv:1908.04619 [astro-ph.CO] .
  * Bengaly (2019) C. A. Bengaly,  (2019), 10.1093/mnrasl/slaa040, arXiv:1912.05528 [astro-ph.CO] .
  * Bengaly  _et al._ (2020b) C. A. Bengaly, C. Clarkson, M. Kunz, and R. Maartens, (2020b), arXiv:2007.04879 [astro-ph.CO] .