Effective Field Theory of Dark Energy: a Review

TL;DR

This review discusses the effective field theory approach to dark energy, a formalism that models cosmic acceleration, classifies modified gravity theories, and aids in interpreting observational data from galaxy surveys.

Contribution

It provides a comprehensive overview of the effective field theory of dark energy, highlighting its development, applications, and role in analyzing cosmological data.

Findings

Systematic classification of modified gravity theories

Predictions for background and perturbation evolution

Framework for interpreting future galaxy survey data

Abstract

The discovery of cosmic acceleration has triggered a consistent body of theoretical work aimed at modeling its phenomenology and understanding its fundamental physical nature. In recent years, a powerful formalism that accomplishes both these goals has been developed, the so-called effective field theory of dark energy. It can capture the behavior of a wide class of modified gravity theories and classify them according to the imprints they leave on the smooth background expansion history of the Universe and on the evolution of linear perturbations. The effective field theory of dark energy is based on a Lagrangian description of cosmological perturbations which depends on a number of functions of time, some of which are non-minimal couplings representing genuine deviations from General Relativity. Such a formalism is thus particularly convenient to fit and interpret the wealth of new data that will be provided by future galaxy surveys. Despite its recent appearance, this formalism has already allowed a systematic investigation of what lies beyond the General Relativity landscape and provided a conspicuous amount of theoretical predictions and observational results. In this review, we report on these achievements.

Tables

Table 1. Table 2: Table of acronyms and their definitions used throughout the review.

Acronym	Definition
BAO	Baryon Acoustic Oscillations
BICEP	Background Imaging of Cosmic Extragalactic Polarization
BOSS	Baryon Oscillation Spectroscopic Survey
CDM (DM)	Cold Dark Matter (Dark Matter)
C.L.	Confidence Level
CMB	Cosmic Microwave Background
CMB-S4	Stage-4 CMB experiment
COrE	Cosmic Origins Explorer mission
CPL	Chevallier-Polarski-Linder
CS	Cosmic shear
DE	Dark Energy
DESI	Dark Energy Spectroscopic Instrument
DHOST	Degenerate Higher Order Scalar-Tensor Theories
DoF	Degree of Freedom
EFT	Effective Field Theory
EB	Einstein-Boltzmann
FLRW	Friedmann-Lemaître-Robertson-Walker
GBD	Generalized Brans-Dicke theories
GC	Galaxy Clustering
GLPV	Gleyzes-Langlois-Piazza-Vernizzi
GR	General Relativity
GWs	Gravitational Waves
ICs	Initial conditions
ISW	Integrated Sachs-Wolfe
JLA	Joint Light-curve Analysis
KiDS	Kilo-Degree Survey
$\Lambda$CDM	$\Lambda$ Cold Dark Matter
LIGO	Laser Interferometer Gravitational-Wave Observatory
LISA	Laser Interferometer Space Antenna
LSS	Large-scale structure
LSST	Large Synoptic Survey Telescope
MCMC	Markov Chain Monte-Carlo
MG	Modified Gravity
QS	Quasi Static
RSD	Redshift-space distortions
SKA	Square Kilometer Array
SNIa	Supernovae Ia
$w$CDM	$w$ Cold Dark Matter
WEP	Weak Equivalence Principle
WL	Weak Lensing

Table 2. Table 3: Table of symbols and their definitions used throughout the review.

Symbol	Definition
$t$, $z$, $a\left(t\right)$, $k$	cosmic time, redshift, scale factor, wavenumber
${\partial }_{\mu },{\nabla }_{\mu }(=;)$	derivative and covariant derivative
$𝔣\left(t\right)$, c(t), $\Lambda \left(t\right)$	(background) EFT functions
$M_2^4,{\overline{M}}_2^2,{\overline{M}}_3^2,{\overline{m}}_1^3,m_2^2,{\mu }_1^2$	EFT functions
$M_{\mathrm{pl}}$	Planck mass
$M\left(t\right)$	effective Planck mass
${\alpha }_M\left(t\right)$	running Planck mass
${\alpha }_B\left(t\right)$	braiding function
${\alpha }_B^{GLPV}\left(t\right)$	beyond GLPV function
${\alpha }_T\left(t\right)$	tensor speed excess
${\alpha }_K\left(t\right),{\alpha }_{K_2}\left(t\right)$	kineticity, extended kineticity
${\alpha }_H\left(t\right)$	beyond Horndeski function
$\Omega \left(t\right),{\gamma }_i\left(t\right)$	EFTCAMB basis, $i=1,\mathrm{..6}$
$H\left(t\right)\left(H_0\right)$	Hubble function (today)
$B\left(t\right)\left(B_0\right)$	scalaron Compton wavelength (today)
${\delta }_m$	linear matter perturbation
$\varphi (t,x_i)$, $X={\partial }_{\mu }\varphi {\partial }^{\mu }\varphi$	scalar field, its kinetic term
$g_{\mu \nu }\left(x_{\mu }\right)$	metric tensor
$\delta A(t,x_i)$	linear perturbation of A
$R$, $R_{\mu \nu }$	Ricci scalar and tensor
${\Delta }_m(t,x_i)$	Comoving density contrast
${\rho }_m\left(t\right),{\rho }_{\mathrm{DE}}\left(t\right)$	matter and DE densities
${\Omega }_m\left(t\right)\left({\Omega }_{m,0}\right)$	matter density parameter (today)
${\Omega }_{\mathrm{DE}}\left(t\right)\left({\Omega }_{\mathrm{DE},0}\right)$	DE density parameter (today)
$w_{\mathrm{DE}}\left(t\right)$	DE equation of state
$w_0$, $w_a$	constant $w_{\mathrm{DE}}$, first derivative of $w_{\mathrm{DE}}$ today
$c_s{\left(t\right)}^2,c_t{\left(t\right)}^2$	scalar and tensor speeds of propagation
$\alpha \left(t\right)$	kinetic term in Horndeski
${\lambda }_C\sim M_C^{-1}$	Compton length scale/Mass
${\sigma }_8\left({\sigma }_{8,0}\right)$	matter power spectrum amplitude at 8 h-1Mpc (today)
$f\left(t\right),\gamma \left(t\right)$	growth rate and growth index
$f{\sigma }_8\left(t\right)$	growth function
$G_N$	Newton constant
$\Phi (t,k),\Psi (t,k)$	Newtonian and curvature potentials
$\mu (t,k)$	effective gravitational coupling
$\Sigma (t,k)$	light deflection parameter
$\eta (t,k)$	gravitational slip parameter
${\mu }_{\mathrm{sc}}\left(t\right)/{\Sigma }_{\mathrm{sc}}\left(t\right)$	$\mu /\Sigma$ at super-Compton scale
${\mu }_{\infty }\left(t\right)/{\Sigma }_{\infty }\left(t\right)$	$\mu /\Sigma$ at sub-Compton scale
${\mu }_{ff}\left(t\right)$	fifth-force contribution to ${\mu }_{\infty }$
$\pi (t,k)$	scalar field perturbation
$K$, $K_{\mu \nu }$	extrinsic curvature trace and tensor
$ℛ,{ℛ}_{ij}$	three dimensional Ricci scalar and tensor

Table 3. Table 4: Different basis used in literature to identify the EFT functions. We used the α 𝛼 \alpha -basis [ 128 ] as reference. We note that in this review we have used M = M ∗ 𝑀 subscript 𝑀 M=M_{*} . – means the basis has not been extended to the GLPV case yet.

Ref.	${𝑴}_{\ast }^{\mathrm{𝟐}}$	${𝑴}_{\ast }^{\mathrm{𝟐}}𝑯{𝜶}_{\text{M}}$	${𝑴}_{\ast }^{\mathrm{𝟐}}{𝑯}^{\mathrm{𝟐}}{𝜶}_{\text{K}}$	${𝑴}_{\ast }^{\mathrm{𝟐}}𝑯{𝜶}_{\text{B}}$	${𝑴}_{\ast }^{\mathrm{𝟐}}{𝜶}_{\text{T}}$	${𝑴}_{\ast }^{\mathrm{𝟐}}{𝜶}_{\text{H}}$
[343]	$w_1$	${\dot{w}}_1$	$\frac{2}{3}{w}_3+6H{w}_2-6H^2{w}_1$	$-w_2+2H{w}_1$	$w_4-w_1$	–
[45, 47]	$m_0^2\Omega +{\overline{M}}_2^2$	$m_0^2\dot{\Omega }+{\dot{\overline{M}}}_2^2$	$2c+4M_2^4$	$-{\overline{M}}_1^3-m_0^2\dot{\Omega }$	$-{\overline{M}}_2^2$	$2\hat{M}-{\overline{M}}_2^2$
[344]	${𝒢}_T$	${\dot{𝒢}}_T$	$2\Sigma +12H\Theta -6H^2{𝒢}_T$	$-2\Theta +2H{𝒢}_T$	${ℱ}_T-{𝒢}_T$	–
EFT basis	$M_{\mathrm{pl}}^2𝔣+{\overline{M}}_2^2$	$M_{\mathrm{pl}}^2\dot{𝔣}+{\left({\overline{M}}_2^2\right)}^{\cdot }$	$2c+4M_2^4$	$-\left(M_{\mathrm{pl}}^2\dot{𝔣}+{\overline{m}}_1^3\right)$	$-{\overline{M}}_2^2$	$2{\mu }_1^2-{\overline{M}}_2^2$
[46, 80]	$M_{\ast }^{2}f+2m_4^2$	$M_{\ast }^2\dot{f}+2{\left(m_4^2\right)}^{\cdot }$	$2c+4M_2^4$	$m_3^3-M_{\ast }^2\dot{f}$	$-2m_4^2$	$2\left({\tilde{m}}_4^2-m_4^2\right)$
EFTCAMB basis	$m_0^2\left(1+\Omega +{\gamma }_3\right)$	$m_0^2\left(\dot{\Omega }+{\dot{\gamma }}_3\right)$	$2c+4H_0^2{m}_0^2{\gamma }_1$	$-m_0^2\left(H_0{\gamma }_2+\dot{\Omega }\right)$	$-m_0^2{\gamma }_3$	$m_0^2\left(2{\gamma }_5-{\gamma }_3\right)$
[85, 80]	$M^2\left(1+{ϵ}_4\right)$	${\left(M^2\left(1+{ϵ}_4\right)\right)}^{\cdot }$	$2M^2\left(𝒞+2{\mu }_2^2\right)$	$-M^2\left(\mu -{\mu }_3\right)$	$-M^2{ϵ}_4$	$M^2\left({\widetilde{ϵ}}_4-{ϵ}_4\right)$

Table 4. Table 5: Summary of the cosmological constraints on the MG parameters discussed in Sections 4.1 , 4.2 and 4.3 for the α 𝛼 \alpha -basis and linear-de form eq. ( 193 ). We use the notation adopted in the original papers. If not specified otherwise the background assumed is Λ Λ \Lambda CDM. A (*) denotes ranges we estimated (by eye) from the marginalized contour plots in the respective papers.

Linear-de form eq. (193)
Basis	EFT functions	Constraints	Dataset	Ref.
$\alpha$-basis	${\alpha }_M$	${\alpha }_M>-1.6$ (*)	CMB + $H_0$ prior 95.4% C.L.	[75]
	${\alpha }_M$, ${\alpha }_B$ ${\alpha }_K$, $c_t^2=1$	${\hat{\alpha }}_M={0.25}_{-0.29}^{+0.19}$ ${\hat{\alpha }}_B={0.20}_{-0.33}^{+0.20}$ ${\hat{\alpha }}_K=0$ (fixed)	KiDS+GAMA 95% C.L.	[214]
		$c_M={0.20}_{-0.82}^{+1.15}$ $c_B={0.63}_{-0.62}^{+0.83}$ $c_K=0.1$ (fixed)	CMB+BAO +RSD+mPk 95%C.L.	[220]
	${\alpha }_M$, ${\alpha }_B$, ${\alpha }_T$, ${\alpha }_K$	$c_K=10$(fixed)	CMB+BAO +RSD+PK 95% C.L.	[205]
	$({\alpha }_M,{\alpha }_B)\times S\left(\frac{k}{k_V}\right)$, ${\alpha }_K\times S\left(\frac{k}{k_V}\right)$, $c_t^2=1$, see eq. (72)	${\hat{a}}_K=0.056$ (fixed) $\sigma \left({\hat{a}}_M\right)=0.065$ $\sigma \left({\hat{a}}_B\right)=0.049$	CMB+GC+CS (forecasts)	[227]
		${\hat{a}}_K=0.01$ (fixed) ${\hat{a}}_M=126\%$ ${\hat{a}}_B=41\%$	CS: 3DWL linear (forecasts)	[230]
		${\hat{a}}_K=0.01$ (fixed) ${\hat{a}}_M=158\%$ ${\hat{a}}_B=54\%$	CS: tomography linear (forecasts)	[230]
	$({\alpha }_M,{\alpha }_B)\times S\left(\frac{k}{k_V}\right)$, $({\alpha }_K,{\alpha }_T)\times S\left(\frac{k}{k_V}\right)$, see eq. (72)	$\sigma \left(c_M\right)=0.056$ $\sigma \left(c_B\right)=0.123$ $\sigma \left(c_K\right)=3.1$ $\sigma \left(c_T\right)=0.146$	S4+LSST +SKA1-IM+DESI (forecasts)	[146]
	${\alpha }_M,{\alpha }_B,$ ${\beta }_{\gamma }^2$ (see eq.(73))	$\sigma \left({\alpha }_{M,0}\right)=0.0146$ $\sigma \left({\alpha }_{B,0}\right)=0.0030$ $\sigma \left({\beta }_{\gamma }^2\right)=0.00135$	GC+WL ISW-Galaxy (forecasts)	[228]
	${\alpha }_M$, ${\alpha }_B$, ${\alpha }_K$, ${\alpha }_H$, $c_t^2=1$	(*) (*) ${\alpha }_K$ (fixed)	CMB+BAO+RSD 95% C.L.	[145]

Table 5. Table 6: The same as in Table 5 for both the α 𝛼 \alpha and EFTCAMB basis and constant, early/late transition form eq. ( 201 ) and scaling- a 𝑎 a form eq. ( 196 ).

Basis	EFT functions	Constraints	Data sets	Ref.
Constant form
$\alpha$-basis	$M^2$, ${\alpha }_B$, ${\alpha }_T$	$\sigma \left({\tilde{M}}_0\right)=0.006$ $\sigma \left({\alpha }_0^B\right)=0.02$ $\sigma \left({\alpha }_0^T\right)=0.001$	CMB-S4+DESI (forecasts)	[70]
EFTCAMB basis	$\Omega$, ${\gamma }_2$, ${\gamma }_3$	$\sigma \left({\Omega }_0\right)=0.01$ $\sigma \left({\gamma }_2^{\left(0\right)}\right)=0.05$ $\sigma \left({\gamma }_3^{\left(0\right)}\right)=0.003$	CMB-S4+DESI (forecasts)	[70]
early/late transition form eq. (201)
$\alpha$-basis	$M^2$, ${\alpha }_B$, ${\alpha }_T$	$\sigma \left({\tilde{M}}_{early}\right)=0.05$ $\sigma \left({\tilde{M}}_{late}\right)=0.007$ $\sigma \left({\alpha }_{early}^B\right)=0.04$ $\sigma \left({\alpha }_{late}^B\right)=0.08$ $\sigma \left({\alpha }_{early}^T\right)=0.02$ $\sigma \left({\alpha }_{late}^T\right)=0.002$	CMB-S4+DESI (forecasts)	[70]
EFTCAMB basis	$\Omega$	$\sigma \left({\Omega }_{early}\right)=0.03$ $\sigma \left({\Omega }_{late}\right)=0.02$	CMB-S4+DESI (forecasts)	[70]
Scaling-$a$ form eq. (196)
$\alpha$-basis	${\alpha }_M=-{\alpha }_B$	${\alpha }_{M0}=-{0.015}_{-0.017}^{+0.019}$ $\beta ={0.66}_{-0.21}^{+0.44}$	Planck18+WL +BAO/RSD 68%C.L.	[1]
	${\alpha }_M$, ${\alpha }_B$, ${\alpha }_K$, $c_t^2=1$	$c_M={0.27}_{-0.26}^{+0.54}$ $c_B={0.48}_{-0.46}^{+0.83}$ $c_K=0.1$ (fixed)	CMB+BAO +RSD+mPk 95%C.L.	[220]
EFTCAMB basis	$\Omega$		Planck13+WP +BAO+Lensing 95%C.L.	[73]
	$\Omega$, CPL	${\Omega }_0=-{0.07}_{-0.18}^{+0.17}$ $s_0>0.435$	CMB+BAO +SNIa+WL 95% C.L.	[135]
		$2\sigma \left({\Omega }_0\right)=110\%$ $2\sigma \left(s_0\right)=68\%$	GC+WL+CMB (forecasts)	[135]
	$\Omega$, ${\gamma }_1$, ${\gamma }_2$ $c_t^2=1$, CPL	${\Omega }_0={0.03}_{-0.25}^{+0.31}$ $s_0>0.215$ ${\gamma }_1^0>0.217$ ${\gamma }_2^0=-{0.9}_{-2.0}^{+1.3}$ $s_2>0.330$	CMB+BAO +SNIa+WL 95% C.L.	[135]
		$2\sigma \left({\Omega }_0\right)=128\%$ $2\sigma \left(s_0\right)=96\%$ $2\sigma \left({\gamma }_2^0\right)=240\%$ $2\sigma \left(s_2^0\right)=136\%$ ${\gamma }_1^0=5$, $s_1=1.4$ (fixed)	GC+WL+CMB (forecasts)	[135]

Table 6. Table 7: The same as in Table 5 for the e-fold form eq. ( 198 ), de-1 form eq. ( 194 ), de-density form eq. ( 197 ).

Basis	EFT functions	Constraints	Data sets	Ref.
e-fold form eq. (198)
$\alpha$-basis	${\alpha }_B=-2{\alpha }_M$	(*) $a_t=0.5$,$\tau =1.5$ (fixed)	RSD (DESI) (forecast) 68% C.L.	[192]
	${\alpha }_B=-2{\alpha }_M$ CPL		CMB+BAO +RSD+SNIa 95% C.L.	[207]
de-1 form eq. (194)
[85]	$\mu$, ${\mu }_3$, ${ϵ}_4$	$p_1=-{0.28}_{-0.20}^{+0.17}$ $p_3=0.04\pm 0.17$ $p_4=-{0.030}_{-0.035}^{+0.068}$	CMB (Planck,WMAP) 68% C.L.	[112]
		$p_1={0.10}_{-0.37}^{+0.58}$ $p_3={0.13}_{-0.40}^{+0.28}$ $p_4=-{0.18}_{-0.13}^{+0.28}$ $p_3^1={0.41}_{-0.91}^{+0.39}$ $p_4^1={0.03}_{-0.11}^{+0.18}$	CMB (Planck,WMAP) 68% C.L.	[112]
	${\mu }_1$, ${\mu }_2^2$, ${\mu }_3$, ${ϵ}_4$	$p_{10}=-{0.000}_{-0.002}^{+0.002}$ $p_{11}=-{0.127}_{-0.096}^{+0.095}$ $p_{20}={1.697}_{-2.157}^{+2.933}$ $p_{21}=-{0.926}_{-5.990}^{+5.852}$ $p_{30}={1.022}_{-0.806}^{+0.930}$ $p_{31}=-{1.447}_{-1.812}^{+1.510}$	$f{\sigma }_8+f+{\sigma }_8+{\dot{G}}_N$ 95% C.L.	[174]
de-density form eq. (197)
EFTCAMB basis	$\Omega$, CPL	${\Omega }_0=-{0.018}_{-0.019}^{+0.032}$	CMB+BAO +SNIa+WL 95% C.L.	[135]
		$2\sigma \left({\Omega }_0\right)=22\%$	GC+WL+CMB (forecasts)	[135]
	$\Omega$, ${\gamma }_1$, ${\gamma }_2$ $c_t^2=1$, CPL	${\Omega }_0={0.047}_{-0.051}^{+0.068}$ ${\gamma }_1^0>0.295$ ${\gamma }_2^0=-{0.23}_{-0.32}^{+0.26}$	CMB+BAO +SNIa+WL 95% C.L.	[135]
		$2\sigma \left({\Omega }_0\right)=48\%$ $2\sigma \left({\gamma }_2^0\right)=40\%$ ${\gamma }_1^0=4.4$ (fixed)	GC+WL+CMB (forecasts)	[135]

Equations

$$d{s}^2=-d{t}^2+a(t{)}^2\left(\frac{d{r}^2}{1-\bar{\kappa }{r}^2}+r^{2}d{\bar{\Omega }}^2\right),$$

$$\varphi (t,\vec{x})=\bar{\varphi }(t)+\delta \varphi (t,\vec{x}),$$

$$n_{\mu }=-\frac{{\partial }_{\mu }\varphi }{\sqrt{-({\partial }_{\mu }\varphi {)}^2}}\to -\frac{{\delta }_{\mu }^0}{\sqrt{-g^{00}}},$$

$$K_{\mu \nu }=h_{\mu }^{\sigma }{\nabla }_{\sigma }{n}_{\nu },K={\nabla }^{\nu }{n}_{\nu },$$

$$S=\frac{1}{2}\int d^{4}x$$

$$+M_2^4(t)(\delta g^{00}{)}^2-{\bar{m}}_1^3(t)\delta g^{00}\delta K-{\bar{M}}_2^2(t)\delta K^2$$

$$-{\bar{M}}_3^2(t)\delta K_{\mu }^{\nu }\delta K_{\nu }^{\mu }+{\mu }_1^2(t)\delta g^{00}\delta R+m_2^2(t)h^{\mu \nu }{\partial }_{\mu }{g}^{00}{\partial }_{\nu }{g}^{00}$$

$$+\dots ]+S_m[g_{\mu \nu },{\chi }_m],$$

$$c$$

$$\Lambda$$

$${\dot{\rho }}_{\mathrm{m}}+3H({\rho }_{\mathrm{m}}+p_{\mathrm{m}})=0,$$

$$H^2+\frac{\bar{\kappa }}{a^2}$$

$$\dot{H}-\frac{\bar{\kappa }}{a^2}$$

$${\dot{\rho }}_{\mathrm{DE}}+3H({\rho }_{\mathrm{DE}}+p_{\mathrm{DE}})=3M_{\mathrm{pl}}^2\dot{\mathfrak{f}}\left(H^2+\frac{\bar{\kappa }}{a^2}\right).$$

$$\displaystyle c\$$

$$\Lambda$$

$${\rho }_{\mathrm{DE}}$$

$$p_{\mathrm{DE}}$$

$$H^2+\frac{\bar{\kappa }}{a^2}$$

$$\dot{H}-\frac{\bar{\kappa }}{a^2}$$

$$t\to \tilde{t}=t+\pi (x^{\mu }),x^i\to {\tilde{x}}^i=x^i,$$

$$\mathfrak{g}(t)\to \mathfrak{g}(t+\pi (t,x^i))=\mathfrak{g}(t)+\dot{\mathfrak{g}}(t)\pi (t,x^i)+\frac{1}{2}\ddot{\mathfrak{g}}(t)\pi (t,x^i{)}^2+...,$$

$$g^{00}$$

$$\delta K_{ij}$$

$$\delta K$$

$${\mathcal{R}}_{ij}$$

$$\mathcal{R}$$

$$A\ddot{\pi }+B\dot{\pi }+(C+k^{2}D)\pi +E=0,$$

$${\mathcal{L}}_{\mathcal{Q}}\sim -\frac{1}{2}(\partial \varphi {)}^2-V(\varphi )\underset{\text{mbox}unitarygauge}{\to }-\frac{1}{2}{\dot{\bar{\varphi }}}^2{g}^{00}-V(\bar{\varphi }),$$

$$c(t)=\frac{1}{2}{\dot{\bar{\varphi }}}^2,\Lambda (t)=V(\bar{\varphi }),$$

$$d{s}^2=-N^{2}d{t}^2+h_{ij}(d{x}^i+N^{i}dt)(d{x}^j+N^{j}dt),$$

$$L=L(N,\mathcal{R},\mathcal{S},K,\mathcal{Z},\mathcal{U},{\mathcal{Z}}_1,{\mathcal{Z}}_2,{\alpha }_1,{\alpha }_2,{\alpha }_3,{\alpha }_4,{\alpha }_5;t),$$

$$\mathcal{S}=K_{\mu \nu }{K}^{\mu \nu },\mathcal{Z}={\mathcal{R}}_{\mu \nu }{\mathcal{R}}^{\mu \nu },\mathcal{U}={\mathcal{R}}_{\mu \nu }{K}^{\mu \nu },$$

$${\mathcal{Z}}_1={\nabla }_i\mathcal{R}{\nabla }^i\mathcal{R},{\mathcal{Z}}_2={\nabla }_i{\mathcal{R}}_{jk}{\nabla }^i{\mathcal{R}}^{jk},{\alpha }_1=a^i{a}_i,$$

$${\alpha }_2=a^i\Delta a_i,{\alpha }_3=\mathcal{R}{\nabla }_i{a}^i,{\alpha }_4=a_i{\Delta }^2{a}^i,{\alpha }_5=\Delta \mathcal{R}{\nabla }_i{a}^i,$$

$$S_{ADM}$$

$$\mathcal{A}$$

Full text

††journal: Physics Reports

1. 1 Introduction

2. 2 Effective Field Theory of Dark Energy

3. 2.1 Action and its formulation

4. 2.2 Modified background Friedmann equations

5. 2.3 The Stückelberg trick

6. 2.4 Mapping: a general recipe

7. 2.5 Stability conditions

8. 2.6 $\alpha$-basis: a phenomenological parameterization

9. 2.7 Couplings with matter fields

10. 2.8 Einstein-Boltzmann codes

11. 3 Novel predictions

12. 3.1 Gravitational couplings

13. 3.2 Phenomenology from $\mu$, $\Sigma$ and $\eta$

14. 3.3 Impact of stability conditions

15. 3.4 Linear growth rate of matter density fluctuations

16. 3.5 Imprints on cosmological power spectra

17. 4 Cosmological constraints

18. 4.1 Running Planck mass

19. 4.2 Pure Horndeski and GLPV models

20. 4.3 Forecasts with next generation surveys

21. 4.4 Gravitational waves

22. 4.5 Constraints on $\mu$, $\Sigma$ and $\eta$

23. 4.6 Specific modified gravity models

24. 4.7 Neutrinos and modified gravity

25. 5 Astrophysical implications

26. 5.1 Massive astrophysical bodies

27. 5.2 GW170817 and GRB170817A

28. 6 Conclusion and outlook

29. A Acronyms and symbols

30. B Alternative basis and pure EFT parameterizations

31. C Constraints summary

1 Introduction

The Universe is a physical system whose dynamics is strongly affected by gravity and that is characterized, at very large scales, by invariance under spatial translations and rotations. Other bona fide symmetries of our fundamental Lagrangians, time translations and Lorentz boosts, are spontaneously broken. One curious occurrence for the Universe is that its expansion rate has accelerated twice, once at very primordial times and high energies (inflation) and then through a second much later period that is still ongoing (“Dark Energy”, DE). When trying to address these phenomena, it looks natural to speculate about the nature of gravity itself. Indeed, the necessity of invoking accelerating expansion twice might even cast doubts on the very basic geometrical description inherited from classical General Relativity (GR), when applied to the Universe as a whole.

The standard cosmological model, $\Lambda$CDM, is based on GR and assumes the Universe is made of a DE component in the form of a cosmological constant ($\Lambda$), cold dark matter (CDM) particles and ordinary matter. Although $\Lambda$CDM gives an astonishing description of the Universe [1, 2], the model shows some shortcomings: the so-called cosmological constant problems. These arise from the mismatch between the naturally expected value of $\Lambda$ from quantum corrections and the observed one, and to the late time coincidence problem (see [3, 4, 5] for reviews). Furthermore, some mild observational tensions among different datasets emerge in this model, for instance, on the value of the Hubble constant $H_{0}\,(=100\,h$ km s*-1* Mpc${}^{-1})$ and amplitude of the matter power spectrum at present time and scale of 8 h*-1*Mpc, denoted by $\sigma_{8,0}$. At the level of the background evolution of the Universe, comparing for example Cosmic Microwave Background (CMB) radiation data by Planck [6, 1] and local measurements of $H_{0}$ based on the cosmic distance ladder leads to a tension of $4.4\sigma$ [7, 8, 9]. On the other hand, Baryon Acoustic Oscillations (BAO) measurements from the Baryon Oscillation Spectroscopic Survey (BOSS) [10] and the Sloan Digital Sky Survey (SDSS) [11] show a $2.5\sigma$ discrepancy in $H_{0}$ with Planck [12]. Regarding matter perturbations, a discordance of about $2.3\sigma$ on $\sigma_{8,0}$ can be spotted between the Kilo-Degree Survey (KiDS) [13] and Planck data [14, 13, 15, 16, 17]. This picture summarizes the motivations at the basis of speculations on the validity of the $\Lambda$CDM model and the search for new physics beyond the standard model.

When pondering about the theory of gravity, one immediately faces a general deep lesson: no other theory than GR is compatible with the basic requisites of a single massless spin two field (the graviton) and recovering Lorentz invariance. As corollary, Lovelock’s theorem [18, 19] implies any infrared departure from GR must bring in new degrees of freedom (DoFs). Usually one refers to these proposals as modified gravity (MG) theories [5, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33]. The simplest modification at the basis of both DE and Inflationary models is to introduce an extra scalar DoF to GR thereby accessing the realm of scalar-tensor theories [34, 35, 36, 26, 24, 37, 5, 28, 38, 32]. It is possible to think of this scalar field as the Goldstone field of broken time translations characteristic of an expanding Universe [39]. This specific symmetry breaking pattern is well described within the unitary gauge which is equivalent to a choice of the time coordinate. Spacial diffeomorphisms are therefore left unbroken. By setting the time coordinate proportional to the value of the scalar field, its fluctuation disappears from the dynamics and remains encoded in the metric DoFs. There comes an important advantage from this procedure and the unitary gauge: it allows to write down the most general action for cosmological perturbations with relative ease, providing an effective framework that does not rely on any specific model, yet being a genuine description of departures from GR: a so-called effective field theory (EFT). Regarding expanding universe, this formalism was first applied to inflation [40, 41, 42], then to quintessence models without conformal couplings [43], and finally to DE [44, 45]. The latter case enables one to treat any DE/MG models with one additional scalar DoF in a model-independent approach described through a variety of geometrical operators compatible with the symmetry imposed and accompanied by time dependent functions, namely the EFT functions. By requiring that the linear perturbation equations contain derivatives up to second order, theories of the “Horndeski” type [34, 36] can be recovered quite straightforwardly in the unitary gauge [46, 47]. Even more, working in the unitary gauge has made very natural to explore an entire new class of scalar-tensor theories beyond Horndeski, the so-called Gleyzes-Langlois-Piazza-Vernizzi (GLPV) theories [37]. The EFT approach can also encode the cosmology of Lorentz violating theories [44, 48] such as Hořava gravity [49]. The landscape of theories that, despite the presence of higher derivatives in the equations of motions, avoid Ostrogradsky instabilities [50] has been exhausted by the Degenerate Higher Order Scalar-Tensor Theories (DHOST) [38, 51, 52, 53, 54, 55]. In parallel, the original EFT formalism has been also extended to include those theories with second-order derivative equations of motion with a vector or tensor additional field [56, 57], such as Generalised Proca [58], generalized Einstein-Aether [59] and massive bi-gravity [60]. Finally, while the EFT formalism was originally developed to describe general patterns of DE/MG models at linear cosmological scales, there has also been great progress in extending the framework to include non-linear perturbative effects [61, 62, 63, 64, 65, 66].

A key goal of future surveys such as Euclid [67], Dark Energy Spectroscopic Instrument (DESI) [68], Square Kilometre Array (SKA) [69], Stage-4 CMB experiment (CMB-S4) [70], Large Synoptic Survey Telescope (LSST) [71] is to test gravity on cosmological scales with an exquisite precision to shed light on the phenomenon of late time cosmic acceleration. Cosmological probes such as Galaxy Clustering (GC), CMB, Weak Lensing (WL), Redshift-Space Distortions (RSD), Supernovae Ia (SNIa) and BAO give access to precious physical information about the expansion history and growth of large-scale structures (LSS). The final aim is to confirm the cosmological standard model with even more accuracy, or, eventually, to single out signatures of DE/MG. The systematic exploration of DE/MG models against cosmological data is facilitated by the EFT approach. It offers indeed the possibility to relate the phenomenology of large classes of models directly to cosmological observables. This unified framework allows to explore the realm of gravitational theories in a coherent fashion identifying clear testable patterns. The synergy between theory and observations is further exploited thanks to Einstein-Boltzmann (EB) codes constructed on top of the EFT framework. The latter are EFTCAMB [72, 73], hi_class [74], COOP [75] and EoS_class [76] and they have been used to obtain cosmological constrains using a wide selection of data sets.

In this review, we aim to provide a comprehensive overview of the many developments occurred in the field of MG when treated in the EFT of dark energy framework. The detailed structure of the review is as follows. Section 2 is dedicated to its theoretical aspects which include: the construction of the EFT action, a description of background and linear perturbation equations, the mapping procedure to write any specific scalar-tensor theory in the EFT language and the theoretical stability requirements to guarantee the viability of a gravity model. We also discuss the so-called $\alpha$-basis which reshuffles the EFT functions according to physical properties of the Universe and the extension of the EFT framework to direct couplings between gravity and matter fields. Finally, we review the EB codes which implement the EFT framework. Section 3 describes the links between theory and observations where we highlight the novel predictions brought by the EFT formulation. We discuss the construction of a road map to select stable DE/MG models depending on crucial observations and their implications. We review the studies of large samples of DE/MG models within the EFT framework which identified clear patterns in observables such as the growth of LSS and power spectra. We also summarize the results about the impact of stability conditions on the parameter space. Section 4 is devoted to a discussion of cosmological constraints using present day data and forecasts from future surveys. We attempt to cover all the models analyzed in literature within the EFT framework, these include direct parameterizations of the EFT functions as well as the mapping of specific MG theories in the framework. Section 5 contains an overview about the implications astrophysical constraints have on the parameter space identified by the EFT functions. In particular, we discuss constraints from massive astrophysical bodies such as dwarf, neutron stars, pulsars and galaxy clusters and the recent bounds derived from the joint detection of the GW170817 and GRB170817A events. Section 6 is devoted to final remarks and discussion on the prospect of future inputs into this field.

We provide a guide to acronyms and symbols used in this manuscript respectively in Table 2 and Table 3 (Appendix A). We show the relations among the different EFT basis adopted in literature and a detailed list of the relevant EFT parameterizations used in the analyses in Appendix B. In Appendix C, we give tables summarizing the observational constraints discussed in Section 4.

2 Effective Field Theory of Dark Energy

The EFT is a model-independent approach encompassing all single-field DE/MG models. It describes both the evolution of the cosmological background and linear perturbations. In the following, we review the building blocks of this approach and the mapping procedure which allows to write any specific scalar-tensor theories in the EFT language. We also discuss a phenomenological EFT basis, the so-called $\alpha$-basis, and its extension to include the direct coupling between DE and matter fields. Furthermore, we present the stability requirements a gravity model has to satisfy to be considered theoretically viable. These are particularly needful when dealing with model-independent parameterizations. We conclude by reviewing the existing EB codes which implement the EFT framework.

2.1 Action and its formulation

The EFT action is constructed in the unitary gauge and written in terms of operators compatible with the residual symmetries of unbroken spatial diffeomorphisms. The operators are organized in powers of the number of perturbations and spatial derivatives and these symmetries allow for time-dependent functions multiplying each operator. We review these features:

The background: in the EFT framework, the perturbations are assumed to evolve on a homogeneous and isotropic background, thus a Friedmann-Lemaître-Robertson-Walker (FLRW) line element of the following form is considered:

[TABLE]

where $a(t)$ is the scale factor and $t$ is cosmic time, $\bar{\kappa}$ is the spatial curvature constant and $d\bar{\Omega}^{2}=d^{2}\theta+sin^{2}\theta d^{2}\varphi$.

• 2.

The unitary gauge: the unitary gauge is also known as the velocity orthogonal gauge  [77]. It corresponds to the choice of the basis in which the perturbation of the extra DoF, responsible for the spontaneous symmetry breaking, vanishes. In details, let us consider a scalar field, $\phi(t,\vec{x})$, and its decomposition in a perturbed FLRW Universe as follows

[TABLE]

where $\vec{x}$ are the spatial coordinates, $\bar{\phi}(t)$ is the homogeneous background value of the scalar field and $\delta\phi$ its perturbation. In order to apply the unitary gauge, one has to choose the time coordinate such that $\delta\phi=0$, thus $t$ becomes a function of $\phi$, $i.e.$ $t=t(\phi)$ 111This statement is valid at all perturbative orders as long as $\phi$ is a monotonic function of time.. According to this choice $\phi$ defines a preferred time slicing ($\phi=const$) and constant time hypersurfaces coincide with constant scalar field hypersurfaces. The gradient of the scalar field is thus assumed to be time-like. In a fluid description, the latter implies that the velocity of the scalar field is orthogonal to the constant time hypersurfaces, hence the name “velocity orthogonal gauge”. The assumption of the unitary gauge implies the action will not show any explicit dependence on the scalar field, as such, the latter is said to be “eaten” by the metric, making the manifest number of DoFs minimal. In conclusion, using this gauge implies two major consequences. Firstly, the EFT action will not be constructed in terms of $\phi$ and its perturbations, but will be written only in terms of the metric and geometrical quantities. Secondly, it breaks the full diffeomorphism invariance while leaving unbroken the subgroup of time-dependent spatial diffeomorphisms.

• 3.

The operators: according to the unitary gauge, the operators in the EFT action are constructed with all the invariants under the residual symmetries of unbroken spatial diffeomorphisms. The full EFT action is thus derived considering the metric $g_{\mu\nu}$ and the unit vector $n_{\mu}$, perpendicular to the time slicing, defined as 222Note that in Ref. [45] $n_{\mu}$ is defined with the opposite sign. The different definition can change the sign of some operators in the EFT action, see for instance [78].

[TABLE]

where $g^{00}$ is the time-time component of the inverse metric and $\partial_{\mu}$ is the four dimensional derivative. The covariant derivatives of $n_{\mu}$ must be considered as well, or equivalently, its projection orthogonal to the constant time hypersurfaces, $i.e.$ the extrinsic curvature tensor and its trace

[TABLE]

with the induced metric defined as $h_{\mu\nu}=g_{\mu\nu}+n_{\mu}n_{\nu}$ and $n^{\sigma}\nabla_{\sigma}n_{\nu}\propto h_{\nu}^{\ \;\;\mu}\partial_{\mu}g^{00}$ where $\nabla_{\sigma}$ is the four dimensional covariant derivative. Additional operators to include are also the Ricci scalar $R$ and any curvature invariants and contractions of tensors with $n_{\mu}$, $g_{\mu\nu}$ and derivatives $\nabla_{\mu}$. Additionally, the symmetry of the action is still satisfied if each operator is accompanied by a time dependent function. These free functions of time which scale the perturbations on the FLRW background are called EFT functions. Unlike standard covariant approaches, the EFT action is constructed using a perturbative approach thus each operator is expanded. For example, the operator $A$ will be written as $A(t,x_{i})=\bar{A}(t)+\delta A(t,x_{i})$, where $\bar{A}$ is the background value and $\delta A$ its perturbation. Note this perturbation scheme enables one to write down an action with operators expanded at any order. In a cosmological context, linear perturbations are generally the most focused upon and therefore only operators up to quadratic order are considered. Already at this order, a large number of operators and combinations can be constructed, especially if one considers an arbitrary number of spatial derivatives acting on the perturbed quantities. However, one has to keep in mind that additional spatial derivatives increase the scaling dimension of an operator and the higher the derivatives an operator contains, the less it becomes relevant on large linear scales.

• 4.

Couplings with matter fields: the first version of the EFT framework for dark energy assumes the validity of the weak equivalence principle (WEP) and hence the existence of a metric $g_{\mu\nu}$ universally coupled to the matter fields $\chi_{m}$, making the Jordan Frame the natural choice [44, 45]. This assumption can be relaxed considering a frame where the gravitational interaction between the additional scalar DoF and the matter fields is explicit [79, 80, 81, 82]. In this Section, we follow the historical path and we review the extension of the EFT framework to include direct couplings with matter fields in Section 2.7.

• 5.

Regime of validity and limitations: the regime of applicability of the EFT framework spans from largest cosmological scales down to the ultraviolet (UV) cut-off, $i.e.$ the energy scale at which a classical description of gravity breaks down. In the cosmological framework, one can assume this cut-off $M_{\rm cut-off}$ must be larger than the Hubble parameter at present time $H_{0}$ in order to describe the background cosmology and observable perturbation modes [83]. Extensions of the linear EFT framework to mildly non-linear scales are possible by including non-linear operators and appropriate additional EFT functions to the action (5) [62, 63, 64, 65, 66]. The EFT approach to linear and non-linear scales is thus regulated by: the strong coupling energy scale of the dark energy fluctuations and the screening scale. Physically, the strong coupling scale corresponds to the scale at which the non-linear interactions of the EFT action exit perturbative unitarity [64] thus setting a cut-off above which the EFT looses predictability and where UV completion becomes important. Such scale is characteristic of the specific model under consideration [84, 64]. While the screening scale is associated to the non-linear scale of the extra scalar field. If the screening is weak or fails to operate (see Section 5), the linear EFT approach is still valid, on the other hand if the screening is strong non-linear terms dominate. Other shortcomings of the EFT formulation exist. For example, it does not describe strong gravity nor higher-dimensional regimes. In parallel, the initial limitation of the EFT framework to include only DE/MG models based on a single scalar DoF was overcome recently by a proposal including additional vector and tensor fields [57]. We focus only on the original proposal encompassing DE/MG models with an extra scalar DoF in this review.

Several versions of the EFT action exist in literature [44, 45, 46, 85, 86, 87] which differ in notation, in number of operators and order of perturbations. In this Section, we consider the EFT action for DE/MG up to second order in perturbations, including the most relevant operators encompassing well known DE/MG models. We use the notation as presented for the first time in Ref. [44] and we discuss the other formulations in Appendix B. This EFT action reads [44]

[TABLE]

where $M_{\mathrm{pl}}^{2}$ is the Planck mass, $g^{00}=-1+\delta g^{00}$, $g$ being the determinant of the metric, $\delta R$ and $\delta R_{\mu\nu}$ are the perturbations of the Ricci scalar and tensor, respectively. $S_{m}$ is the matter action for all matter fields $\chi_{m}$ and the EFT functions are $\mathfrak{f}$, $\Lambda$, $c$, $M_{i}$, $m_{i}$, $\bar{M}_{i}$, $\bar{m}_{i}$ and $\mu_{i}$. The action is organized in a specific way: the first line contains the operators contributing to both the background evolution of the Universe and the linear perturbation equations. The corresponding EFT functions are therefore denoted as background EFT functions. The operators in the second and third lines enter only in the linear perturbations equations. The ellipsis stand for additional second order operators which can be included, $e.g.$ $(\delta R)^{2}$ [46] or higher order terms in both derivatives and perturbations. For example Refs. [88, 48, 78] include the operators with higher spatial derivatives to describe theories like high-energy Hořava gravity [89] in the EFT description. Higher order operators, such as $(\delta g^{00})^{3}$ and $(\delta g^{00})^{2}\delta K$, are relevant for an extension of the action to non-linear scales [62, 64, 65]. Additional three operators need to be included [54] in the above action to describe the recently found DHOST [55, 90].

[FIGURE:]

As anticipated, the power of this model-independent construction of a gravitational action lies in the direct connection it has with most of the well known cosmological theories proposed in the last decades. It captures for example the physics at linear cosmological scales of theories such as $f(R)$-gravity [97, 98], Horndeski [34, 36], GLPV theories [37, 80], low-energy Hořava gravity [49] and more, as illustrated in Table 1. The four major sub-sets of theories encoded in the EFT action (5) are:

Scalar-tensor theories à la Brans-Dicke, hereafter dubbed Generalized Brans-Dicke theories (GBD), which are described only by the background EFT functions $\{\mathfrak{f},\Lambda,c\}$; 2. 2.

Horndeski theories for which the relation $\bar{M}_{2}^{2}=-\bar{M}_{3}^{2}=2\mu_{1}^{2}$ (and $m_{2}^{2}=0$) must be imposed; 3. 3.

GLPV theories which require the condition $\bar{M}_{2}^{2}=-\bar{M}_{3}^{2}$ (and $m_{2}^{2}=0$); 4. 4.

Lorentz violating theories, such as Hořava gravity, for which $m^{2}_{2}\neq 0$.

The second and third conditions allow to eliminate from the EFT action spatial derivatives higher than second order.

The advantage of the EFT formalism is twofolds: on one hand it allows to perform model-independent explorations of DE/MG models without assuming any particular theory. This is done by selecting all the EFT functions in action (5) or a sub-set of them usually according to the above itemized list. This procedure is known as pure EFT approach. The investigation of the above four sub-sets of theories in a pure EFT fashion allowed to identify quite general features of the gravity force as further discussed in Section 3. On the other hand, the EFT can be used to investigate specific DE or MG models, $e.g.$ $f(R)$. In this case, the EFT functions assume characteristic forms [44, 45, 47, 46, 104, 48, 78]. We illustrate the mapping recipe to encode specific DE/MG theories in the EFT formulation in Section 2.4. One refers to this procedure as mapping approach. The latter is particularly useful to test a specific theory against cosmological data as explained in Section 2.8.

2.2 Modified background Friedmann equations

In this Section, we proceed further by illustrating the model-independent formulation of the background Friedmann equations for DE/MG models. In order to obtain such equations, one must vary the action (5) with respect to the metric. This procedure yields the modified Friedmann equations

[TABLE]

where the dots correspond to derivatives with respect to cosmic time $t$, $H(t)$ is the Hubble function defined as $H\equiv\frac{1}{a}\frac{da}{dt}$ and $\rho_{\mathrm{m}}(t)$ and $p_{\mathrm{m}}(t)$ are respectively the background energy density and pressure of matter. To complete the set of background equations, one supplies the above system with the matter continuity equations as follows

[TABLE]

where the perfect fluid approximation is assumed. The modified Friedmann equations can be rewritten following the fluid description:

[TABLE]

where $\{\rho_{\rm DE},p_{\rm DE}\}$ are the density and pressure of the dark fluid. Now, differentiating eq. (9) with respect to time and combining it with eqs. (8)-(10), one obtains a “non-conservation” equation for the dark fluid, which reads

[TABLE]

The non-conservation of the dark fluid is thus regulated by the time derivative of the EFT function $\mathfrak{f}$. Using eqs. (9)-(10) in eqs. (6)-(7) allows the simplification

[TABLE]

One could equivalently define an effective density and pressure for the dark fluid $\rho_{\rm DE}^{\mathrm{eff}}$ and $p_{\rm DE}^{\mathrm{eff}}$ which have the forms

[TABLE]

They allow to recover a more standard form of Friedmann equations

[TABLE]

where the dependency on the non-minimal coupling $\mathfrak{f}$ has been hidden in the dark component. From the first line above, one recovers the usual constraint relation: $\Omega_{\rm m}+\Omega_{\rm DE}+\Omega_{\bar{\kappa}}=1$ where the density parameters $\Omega_{i}$ are defined as $\Omega_{\rm m}=\rho_{\mathrm{m}}/3M_{\mathrm{pl}}^{2}H^{2}$, $\Omega_{\rm DE}=\rho^{\mathrm{eff}}_{\rm DE}/3M_{\mathrm{pl}}^{2}H^{2}$ and $\Omega_{\bar{\kappa}}=-\bar{\kappa}/a^{2}H^{2}$. In the following we indicate with $\Omega_{i,0}$ their present day values.

Note that the EFT approach provides a general model-independent description to study the evolution of cosmological perturbations and, in this spirit, one can assume any background expansion. The modified Friedman equations contain three unknown EFT functions $\{\mathfrak{f},\Lambda,c\}$ and the unknown expansion history, $H$. Thus, one has to fix two out of the four free functions, then the remaining can be obtained from the Friedman equations. Usually, the procedure adopted is to assume a chosen form for $\mathfrak{f}$ and an equation of state for the dark fluid, $i.e.$ $w_{\rm DE}$ ($p_{\rm DE}=w_{\rm DE}\rho_{\rm DE}$) in order to fix $H$. Then both $c$ and $\Lambda$ are fixed through eqs. (12)-(13). In this respect, the most common choice is to fix $H(t)$ to that of $\Lambda$CDM. An attempt to find an explicit form for $\mathfrak{f}$ was made in Ref. [105], where the study of the dynamical system associated to the background equations allowed to find appropriate ans$\ddot{\mbox{a}}$tz for $\mathfrak{f}$. One might claim that fixing $H$ to follow a specific expansion history might induce strong hypotheses on the behavior of the underlying gravity theory. Using too specific parametric forms could indeed bias the generality of predictions since one is already selecting a branch of models at the level of the background. To avoid such drawbacks, one can alternatively parameterize $\Lambda$ and $c$ and solve eqs. (6)-(7) to derive $H$ and $\mathfrak{f}$ [106, 107]. Both the approaches are used in the phenomenological analysis presented in Section 3 and cosmological constraints in Section 4.

2.3 The Stückelberg trick

The unitary gauge used to write the action (5) is useful from a theoretical point of view to identify the main operators which introduce modifications at large scales. Furthermore, as presented before, it allows to select relevant sub-classes of models for cosmological purposes straightforwardly, yet it is not convenient to study the evolution of the extra scalar DoF and that of the metric perturbations separately. The reason is because the extra DoF is hidden inside the metric thus an explicit evolution equation for the scalar field cannot be obtained. It is possible to make such field appear explicitly in the action by restoring the full diffeomorphism invariance upon application of the Stückelberg trick 333The Stückelberg trick is a common procedure to study theories with broken gauge symmetries. The introduction of new fields is used to reveal a symmetry of the gauge-fixed theory.[44, 45]. To do so, one has to force back the broken gauge transformation on the field in the Lagrangian by imposing the following time coordinate transformations

[TABLE]

where $\pi$ is the perturbation of the extra DoF. Time translation invariance is thereby restored. The above transformations induces time dependent functions in the action to transform as

[TABLE]

while scalars do not transform. The transformations of the quantities of interest in action (5) are [44, 45]

[TABLE]

where $\mathcal{R}$ is the three dimensional Ricci scalar and $\mathcal{R}_{ij}$ is the corresponding tensor.

Once the above transformations have been inserted in the action (5), one obtains the following action [44, 45]:

[TABLE]

where $\nabla^{i}$ is the spatial covariant derivative.

It is now possible to obtain the dynamical equation for the extra DoF explicitly by varying the action with respect to $\pi$. This yields

[TABLE]

where the spatial part has been Fourier transformed and $k$ is the wavenumber. $A,B,C,D,E$ are functions of time and scale $k$. $E$ also includes couplings with the metric linear perturbations. Their expressions can be found in [86]. Let us note that in order to obtain observable predictions this equation needs to be coupled with the other linear perturbative equations found by varying the above action with respect to each metric component [45, 46].

2.4 Mapping: a general recipe

The relevance of the EFT approach relies in its capability to encompass DE/MG models with a single additional scalar DoF. Along the model-independent explorations of DE/MG (pure EFT approach), it is possible to map specific theories with one additional scalar DoF in the EFT language (mapping approach) [44, 45, 47, 46, 104, 48, 78]. Hereafter we discuss the latter.

In order to map a theory in the EFT language one can follow two paths. The first consists in starting from the covariant action of a specific theory and then impose the unitary gauge. One has then to identify each term in the action with the corresponding one in EFT action. Let us present a simple example. We use the covariant quintessence Lagrangian and we impose the unitary gauge as follows

[TABLE]

where $\delta\phi=0$ is considered according to the definition of unitary gauge. Then, it is straightforward to identify the following correspondence:

[TABLE]

from the EFT action (5). Deriving the mapping relations from a covariant Lagrangian requires to apply these steps for every theory and it can become cumbersome in more involved cases. This lack of generality is solved by the second option [88, 104, 78], $i.e.$ working out the mapping for a general Lagrangian written in unitary gauge with all the relevant operators. As a result each EFT function can be written in terms of this general Lagrangian. Once this general recipe is derived, any model can be translated in the EFT formalism. In the following we illustrate this general approach in details.

For this purpose one introduces the Arnowitt-Deser-Misner (ADM) formalism [108], for which the line element can be written as:

[TABLE]

where $N(t,x^{i})$ is the lapse function, $N^{i}(t,x^{i})$ the shift and $h_{ij}(t,x^{i})$ is the three dimensional spatial metric. Using the ADM formalism, a general Lagrangian describing scalar-tensor theories can be written as function of the following operators [88]:

[TABLE]

where in details:

[TABLE]

with $a_{\nu}=n^{\mu}\nabla_{\mu}n_{\nu}$ being the acceleration of the normal vector and $\Delta=\nabla_{k}\nabla^{k}$. The operators considered in the Lagrangian (30) allow to describe gravity theories containing up to sixth order spatial derivatives. According to the EFT action (5), in the following we consider only some of the operators introduced above, $i.e.$ $L=L(N,{\mathcal{R}},\mathcal{S},K,\alpha_{1},t)$. However, an extended EFT action and the general mapping including the whole set of operators listed in eq. (30) can be found in Refs. [88, 78].

Expanding the Lagrangian (30) up to quadratic order in perturbations of these operators yields the following action [46]:

[TABLE]

where $\bar{L}$ is the background expression of the Lagrangian, $L_{\mathcal{S}}\equiv\partial L/\partial\mathcal{S}$ and equivalently for the others and

[TABLE]

The EFT action (5) has to be written in ADM form as well in order to be compared with eq. (32). After some manipulations (see [46, 88, 78] for details) which include the use of the Gauss-Codazzi relation [109] and the transformation

[TABLE]

where $\delta N$ is the perturbation of the lapse function, from which one can deduce $(\delta g^{00})^{2}=4(\delta N)^{2}$ at second order, one obtains:

[TABLE]

At this point, it is very easy to deduce the following mapping by identification:

[TABLE]

Let us show a practical example by considering the $f(R)$-gravity theory [97, 98]. The mapping of the latter into the EFT language was first derived in Refs. [44, 104], but here we use the above recipe to find the mapping relations. The $f(R)$-gravity action is:

[TABLE]

where $f(R)$ is a general function of the four dimensional Ricci scalar. We expand this action around the background value of the Ricci scalar, $\bar{R}$:

[TABLE]

where $f_{R}\equiv\frac{df}{dR}$. Now we use the Gauss-Codazzi relation [109] to write the above action in the ADM formalism:

[TABLE]

Finally, using eqs. (2.4), one finds the following mapping relations

[TABLE]

and the other EFT functions are zero.

The general recipe of eqs. (2.4) is very handful when implementing a specific model in an EB code using the mapping procedure as illustrated in Section 2.8.

2.5 Stability conditions

We reviewed the construction of a very general framework enclosing many DE/MG models in the previous Sections. Because of the wide generality of the EFT approach, it is of crucial importance to ensure that the theory of gravity under consideration is free from pathological instabilities, such as ghosts, gradient and tachyonic instabilities [110]. When testing gravity models with cosmological data using statistical tools [111, 72, 73, 74], these viability criteria can reduce the parameter space to explore [85, 112, 113] or even dominate over the constraining power of data [73, 48, 113].

Such instabilities are related to the evolution of the extra scalar DoF and when matter fields are involved, they can also contribute to the stability conditions and, as such, alter the viability space of the theory [114, 115, 116, 117, 88, 80, 118, 119]. Hence, the latter need to be consistently considered when analyzing the stability of the whole system. Fundamental then becomes the choice of the matter action. Recently, it has been shown [120, 118] that among the models describing the matter action [121, 122, 114, 115, 116, 123, 117, 88, 80, 79, 82, 119], the more appropriate choice is the Sorkin-Schutz action [121, 122]. This action describes general matter fluids and the canonical field characterizing the matter DoFs is the matter density perturbation, $\delta_{m}$. The latter allows to avoid the problem of a divergent action when including pressure-less matter fluids, such as baryon and dark matter.

One can then construct an action $S$ made by the EFT action (5) ($S_{EFT}$) and the Sorkin-Schutz one ($S_{m}$), $S=S_{EFT}+S_{m}$. The action $S$ then includes one DoF for the gravity sector, namely $\zeta$, defined as the scalar perturbation of $h_{ij}$ 444The ADM metric perturbations for the scalar and tensor components reads: $ds^{2}=-(1+2\delta N)dt^{2}+2\partial_{i}\psi dtdx^{i}+a^{2}\left[(1+2\zeta)\delta_{ij}+h_{ij}^{T}\right]dx^{i}dx^{j}$, where $\delta N(t,\vec{x})$ is the perturbation of the lapse function, $\partial_{i}\psi(t,\vec{x})$ is the scalar perturbation of the shift function, $\zeta(t,\vec{x})$ of the three dimensional metric and $h_{ij}^{T}(t,\vec{x})$ are the perturbed metric components which contribute to tensor modes., as many DoFs, $\delta_{m,i}$, as matter fluids considered and the tensor modes $h_{ij}^{T}$. As a result one obtain an action for scalar modes ($S^{s}$) and one for tensor modes ($S^{T}$). In Fourier space they have the following compact forms:

[TABLE]

where $\vec{\chi}^{t}=(\zeta,\delta_{i})$ is a dimensionless vector and $\textbf{A},\textbf{G},\textbf{B},\textbf{M}$ are matrices whose coefficients are combinations of EFT functions and some matrices also manifest a $k$ dependence 555A similar analysis for the no-ghost and no-gradient conditions can be performed by starting from the action (5) and restoring the broken symmetry by means of the Stückelberg trick (see Section 2.3). In this case the propagating DoF associated to the gravity sector is $\pi$. The stability conditions have been derived with this approach for Horndeski-like models [45, 85].. We refer the reader to [118] for further details and the complete expressions of the matrices. Finally, $M^{2}(t)$ and $c_{t}^{2}$ are respectively the effective Planck mass and the speed of propagations of tensor modes. In section 2.6 we provide their expressions in terms of EFT functions.

Let us now discuss the three main sources of instability related to the above actions:

The Ghost instability corresponds to having modes with negative kinetic energy. In this case the high energy vacuum is unstable to the spontaneous production of particles [124, 125]. Such a pathology is regulated by demanding for a positive kinetic term if only one field is involved, or a positive kinetic matrix if more fields define the system. In the scalar modes action such condition corresponds to requiring A to be positive definite, i.e. all eigenvalues must be strictly positive. The condition is imposed only in the high energy regime because the ghost instability generated in the infrared regime corresponds to the physical phenomenon of the Jeans/tachyonic instability [126], which can be controlled demanding for specific conditions, as we will discuss in the following. Let us note that in the case of the EFT framework, one has to consider that such approach is valid up to a certain cut-off scale, namely $\Lambda_{\rm cut-off}$. Then, when performing high-$k$ expansions one has to assume the following relations between the scales involved, $H\ll k/a\ll\Lambda_{\rm cut-off}$. The no-ghost condition when applied to the tensor modes action (42) reads $M^{2}>0$.

• 2.

The Gradient or Laplacian instability occurs when the DoFs propagate with negative speeds, $i.e.$ $c_{s,i}^{2}<0$ and $c_{t}^{2}<0$. This signals the presence of exponentially growing modes. The regime in which the gradient instability manifests itself is in the high-$k$ regime. The speeds of propagation of each DoF can be identified by computing the field equations associated to the actions (41) and (42) and considering their high-$k$ expansions. To avoid the gradient instabilities, which would be catastrophic for the system, one has to require that the speeds of propagation are positive, $i.e.$ $c_{s,i}^{2}>0$ and $c_{t}^{2}>0$.

• 3.

The Tachyonic and Jeans instabilities are the less severe instabilities and they appear when the DoF has a negative mass squared. In particular, they arises when the perturbations are not computed about the true vacuum of the theory [5]. In order to account for this pathology, one can look at the boundedness of the Hamiltonian at low momenta. These conditions are less explored with respect to the no-ghost/no-gradient conditions and a full and general derivation in the context of the EFT framework is done in Ref. [118]. Starting from the above action (41), it is possible to obtain the associated Hamiltonian, namely $\mathscr{H}(\Phi_{i},\dot{\Phi}_{i})$ of canonical fields $\Phi_{i}$, which in the case of one fluid assumes the form

[TABLE]

where $\mu_{1}$ and $\mu_{2}$ are the mass eigenvalues. The Hamiltonian is unbounded from below if the mass eigenvalues are negative, $i.e.$ $\mu_{i}(t,0)<0$. Requiring $\mu_{i}>0$ would result in a too stringent condition. A less severe request, if the $\mu_{i}$ are negative, is to demand they satisfy the condition $|\mu_{i}(t,0)|\lesssim H^{2}$. In this case the time scale of evolution of the instability is larger than the Hubble time so that it will not affect the stability of the system. Such condition will allow to have $\mu_{i}$ negative at some times. This behavior is known as Jeans instability and it is necessary in order for structures to form.

The pathologies listed above concern both the scalar propagating modes as well as the tensor ones. Let us note that in the case of the scalar action (41) a unique condition for each of the above points, encompassing all viable theories included in the EFT framework, cannot be computed. This is because of the large number of operators involved and the different $k$ dependence they carry, in particular the operators accompanying $\bar{M}^{2}_{3},\bar{M}^{2}_{2}$ and $m^{2}_{2}$. Because the stability requirements discussed above demand for high-$k$ or low-$k$ limits, a certain number of sub-cases need to be considered to account for different powers of $k$ appearing in the action and equations [118]. The three most relevant cases are [118]: 1) the general case for which $m_{2}^{2}\neq 0$ and $\bar{M}_{2}^{3}+\bar{M}^{2}_{2}\neq 0$; 2) the GLPV case which includes Horndeski as sub-case; 3) the Hořava gravity like case. Regardless of the sub-cases considered, the matter fields do not modify the ghost conditions. The same is not true for the case of gradient and tachyon instabilities where the combination of matter and gravity fields is non trivial. For example, in the case of GLPV, it has been found that the speed of propagation of the gravity DoF and that of the radiation fluid do not decouple even at high-$k$ [37, 88, 82, 118]. Most of the investigations we review in the next Sections are limited to the case of Horndeski theories. For this case, the stability requirements concerning the absence of ghosts and positive speeds of propagation become very simple because the $k$ dependence in the matrices disappears and the matter fields do not affect any of them. On the contrary the mass eigenvalues are strongly modified by the presence of matter fields. Regardless of the sub-case considered, if direct couplings between the extra scalar field and matter are included, the no-ghost and no-gradient conditions are modified [79, 81].

Finally, it has been noticed that the conditions for the absence of ghost and gradient instabilities do not dependent on the gauge choice, on the contrary the expressions for the mass eigenvalues do [127]. The Hamiltonian (43) indeed is written in terms of fields which may not have a clear physical interpretation. Thus in order to look for the proper tachyonic conditions, it would be more appropriate to consider the following gauge invariant quantity describing the linear density perturbation of the DE field [127]:

[TABLE]

where $\psi$ and $\gamma$ are the scalar perturbations respectively of the shift function and of the metric tensor of the three dimensional spatial slices. Standard bars stand for background values. This definition is very general and applicable both in the presence of matter fields and in the late time Universe. The action (41) can be written in terms of $\delta_{\phi}$. The resulting stability requirements in the de Sitter limit, which arise by imposing the absence of ghost and gradient instabilities, do not change when considering the $\zeta$ field nor $\delta_{\phi}$. On the contrary, the mass term for the two fields is distinctively different [127]. Further analysis in this direction is still necessary in order to extend the results of [127] by adopting the procedure illustrated in [118]. This would allow to definitively identify the mass of the physical perturbation field $\delta_{\phi}$ in presence of matter fields.

The relevance of imposing physical motivated stability requirements is not limited to demanding a viable gravity theory but they also have power to constrain models parameters. We review and discuss some considerable results in Section 3.3.

2.6 $\alpha$-basis: a phenomenological parameterization

An alternative parameterization of the EFT action, dubbed the $\alpha$-basis, was developed in [128] in order to describe specific physical properties of the Horndeski theory. In that case any departure from GR is described by four time dependent phenomenological functions, namely $\alpha_{M}(t)$, $\alpha_{B}(t)$, $\alpha_{K}(t)$, $\alpha_{T}(t)$. The original $\alpha$-basis was generalized later to include GLPV models by adding an additional function, $\alpha_{H}(t)$ [80, 104], and finally it was further developed to include higher spatial derivatives operators accounting for Lorentz violation, $e.g.$ $\alpha_{K_{2}}(t)$, $\alpha_{B}^{GLPV}(t)$ [78]. The $\alpha$-basis has the benefit of relating the evolution of the coupling functions to clear physical effects, hence it is a more phenomenological approach.

The quadratic action in the $\alpha$-basis encompassing Horndeski, GLPV and low-energy Hořava gravity can be written in ADM formalism and Fourier space as follows [78]

[TABLE]

where the geometrical quantities with tildes are the Fourier transforms, $\delta_{2}$ refers to taking the expansion at second order in perturbations, $h$ is the determinant of the spatial metric $h_{ij}$. The identification with the EFT basis reads

[TABLE]

where $M^{2}(t)=M_{\mathrm{pl}}^{2}\mathfrak{f}-\bar{M}^{2}_{3}$ is the effective Planck mass and $c_{t}$ is the speed of propagations of tensor modes or equivalently of gravitational waves (GWs). One can complement the above functions with the running of the effective Planck mass,

[TABLE]

which characterizes the evolution rate of the effective Planck mass.

Let us now discuss the physical interpretation of the above basis [128, 80, 104, 78, 129]:

$\alpha_{M}$ is the running Planck mass. As specified above, this function parametrizes the time evolution of the effective Planck mass. A running Planck mass modifies the growth of structures, introduces anisotropic stress and modifies the friction term in the GW equation (74).

• 2.

$\{\alpha_{B},\alpha_{B}^{GLPV}\}$: $\alpha_{B}$ is the braiding function [128] 666The definition of $\alpha_{B}$ presented here follows the one in Ref. [128]. We note that it differs by a minus sign and a factor 2 from the one defined in Ref. [104, 78].. It describes the mixing between the metric and the DE field. $\alpha_{B}$ is different from zero for all the theories showing non-minimal coupling to gravity and/or possessing the $\delta N\delta K$ operator in the action, $i.e.$ $f(R),L_{3}^{H},L_{4}^{H},L_{5}^{H},L_{4}^{GLPV},L_{5}^{GLPV}$, where $L_{i}^{H}$ and $L_{i}^{GLPV}$ are respectively the Lagrangians of Horndeski and GLPV theories. This operator does not appear when one considers quintessence and K-essence models ($L_{2}^{H}$). The additional function $\alpha^{GLPV}_{B}$ extends the braiding effect to scalar-tensor theories beyond GLPV. Both braiding functions take place in the kinetic term and the speed of propagation of the scalar mode, hence impacting the clustering properties of DE.

• 3.

$\{\alpha_{K},\alpha_{K2}\}$: $\alpha_{K}$ is called kineticity and is purely a kinetic function and $\alpha_{K2}$ is its extension to Lorentz violating theories. They both enter the definition of the kinetic term. They affect the speed of propagation of the DE field hence the condition for the absence of a scalar ghost. In particular, large positive values of these functions suppress the sound speed of scalar perturbations. The $\alpha_{K}$ function is characteristic of theories belonging to GLPV models for which $\alpha_{K2}$ vanishes. The opposite holds for low-energy Hořava gravity. $\alpha_{K}$ is the only coupling present in quintessence or perfect-fluid DE models.

• 4.

$\alpha_{T}$ is the tensor speed excess and describes the deviation of the speed of propagation of GWs from the speed of light. This function is present in Horndeski, GLPV and low-energy Hořava gravity. It affects the evolution of the scalar gravitational potentials leading to anisotropic stress (see Section 3.1).

• 5.

$\alpha_{H}$ characterizes the departure from Horndeski theories. It contributes to the speed of propagation of the scalar DoF and couples the gravitational field to the velocity of matter [37]. This function is present in GLPV and low-energy Hořava gravity models. In particular, in the latter case this function can be further extended with additional functions associated with higher order spatial derivatives terms in high-energy Hořava gravity [78] for example.

The above basis was carefully built to consider the different phenomenological aspects of the DE fluid. However, let us notice that the desired correspondence between the $\alpha$-functions and physical effects becomes weaker when going beyond the Horndeski class.

2.7 Couplings with matter fields

The EFT formalism described through the action (5) assumes the matter fields, $\chi_{m}$, are minimally coupled to gravity through a unique metric $L_{m}(g_{\mu\nu},\chi_{m})$. This frame dubbed Jordan frame is the standard frame where the interpretation of cosmological measurements are performed. The reason for using the Jordan frame relies on the fact that stringent constraints exist on the couplings between the extra DoF and the standard matter species, $i.e.$ baryons and photons [130, 131]. As a consequence such couplings are chosen to be minimal, $i.e.$ the matter fields are not coupled to the scalar curvature. However, in the case of dark matter and neutrinos, observational constraints are less severe and one has more freedom. One can consider therefore a frame where the gravitational interaction between the additional scalar DoF and the matter fields is explicit. The so-called Einstein frame. In this frame, the metric $\hat{g}_{\mu\nu}$ is related to the Jordan frame metric $g_{\mu\nu}$ by a conformal/disformal transformation as follows 777Let us note that the transformation (48) can be more general, for instance $C_{i}$ can also depend on $X$ or even it can include an extended disformal term given by a rank-two symmetric tensor [132]. :

[TABLE]

where $C_{i}$ and $D_{i}$ are respectively the conformal and disformal couplings for each matter species $i$ and $X=\partial_{\mu}\phi\partial^{\mu}\phi$. The condition $C_{i}>0$ is necessary to preserve the signature of the Jordan frame metric. Let us note that the transformation (48) preserves the structure of the Lagrangian in GLPV theories [80], while in the case of Horndeski theories the structure is preserved if and only if $D_{i}(\phi)$ [133, 132].

Couplings with matter fields have been investigated in the EFT framework for Horndeski theory [79] and later generalized to GLPV [80, 81, 82]. We treat hereafter the general case of GLPV [82]. Using the unitary gauge and the ADM formalism, the disformal transformation in eq. (48) reads

[TABLE]

Then, one can add to the $\alpha$-basis in action (45) for each matter species three coupling functions defining the conformal and disformal interactions [82]

[TABLE]

Note that for the case of Horndeski theories $\alpha_{X,i}=0$. One has a total of $3N_{i}+5$ time-dependent functions if $N_{i}$ matter species are present where the +5 are the standard coupling functions characterizing GLPV theories. However, the arbitrariness in the choice of the gravitational metric used to define the gravitational and matter sectors makes three of these functions redundant. The structure of the action is preserved under transformations of the reference metric thus the number of physically relevant functions of time reduces to $3N_{i}+2$ [82].

The inclusion of disformal couplings has been shown to impact the no-ghost and no-gradient stability requirements [79, 81]. The disformal couplings indeed contribute to the kinetic term of the scalar modes, thus the condition for the absence of ghost is modified and in turn the speed of propagation as well. Furthermore, the corresponding conditions are frame independent [79, 81].

2.8 Einstein-Boltzmann codes

In order to perform explorations of cosmological observables and constrain cosmological and model parameters, one option is to modify existing EB codes based on $\Lambda$CDM. The latter allow to numerically evaluate the linear evolution of relevant perturbed quantities in DE/MG scenarios. As a result, general purpose codes have been developed using the EFT formalism in order to encompass a wide range of DE/MG models and to allow for model-independent exploration of their properties (pure EFT approach) as well as for a full implementation of specific models once the mapping has been worked out (mapping approach). In the pure EFT approach, the user has to specify the expansion history, $i.e.$ choose the functional form of two out of the four background functions ($\mathfrak{f}$, $w_{\rm DE}(a)$, $c$ and $\Lambda$) and the functional forms of the other EFT functions. In the mapping approach, the user has to implement a background solver to find the expansion history, $H$, for the desired model and provide the mapping relations for the EFT functions as discussed in Section 2.4. A variant in the background solver is the so-called designer approach, $i.e.$ the expansion history is fixed ($w_{\rm DE}(a)$) and the dynamical equation for the gravity field is solved. The EFT functions are then specified using the mapping recipe. Such EB codes allow to evolve the full set of linear perturbative scalar and tensor equations without relying on any Quasi Static (QS) approximation. The QS approximation consists in neglecting the time derivatives of linear perturbations and it is usually applied within the sound horizon of the DE mode, $k/aH>c_{s}$ [134, 135]. However, the latter is a necessary condition for the QS approximation to hold but not sufficient to exploit the full dynamics of the extra DoFs [136].

A collection of EB solvers based on the EFT framework has been compared and cross-calibrated recently [137, 76]. Among them, there are: the Effective Field Theory for CAMB (EFTCAMB) [72, 73], Horndeski in CLASS (hi_class) [74], Cosmology Object Oriented Package (COOP) [75] and Equation of State for CLASS (EoS_class) [76]. All are publicly available 888 EFTCAMB webpage: www.eftcamb.org; hi_class: webpage: www.hiclass-code.net; COOP webpage: www.cita.utoronto.ca/~zqhuang; EoS_class webpage: https://github.com/fpace.. Here we briefly describe their main features:

EFTCAMB is a patch for the public EB solver CAMB [138] which implements the EFT approach using the formulation presented in [45]. The set of EFT functions used in EFTCAMB is in Appendix B (we refer to it as EFTCAMB basis hereafter). The code has built-in models which include specific theories such as designer $f(R)$-gravity [73], $f(R)$-Hu Sawicki model [139], minimally coupled quintessence [86], low-energy Hořava gravity [48], covariant Galileon [140], K-mouflage [141], Galileon Ghost Condensate [142], beyond Horndeski model [143] as well as several model-independent parameterizations of the DE equation of state and choices for the EFT functions. It allows to use the $\alpha$-basis as well [128]. A novelty introduced in the EFTCAMB patch is the built-in stability module which checks for the viability of the underlying theory of gravity by imposing the full set of physical conditions discussed in Section 2.5, $i.e.$ no-ghost and no-tachyonic conditions and a positive speed of propagation. The resulting viable parameter space is supplied as prior when using the Markov Chain Monte-Carlo (MCMC) code named EFTCosmoMC [73].

• 2.

hi_class is a modified version of the public EB solver CLASS [144] which implements Horndeski theory by using the formalism of the $\alpha$-basis, $i.e.$ $\alpha_{M},\alpha_{K},\alpha_{T},\alpha_{B}$ [128]. It has been recently extended to include also $\alpha_{H}$ which is characteristic of GLPV models [145]. The code includes the EFTCAMB basis as well. The code is comprised of several built-in specific theories such as quintessence, with different choices for the potential, Brans-Dicke theory [146], covariant Galileon [147] and also allows for a model-independent exploration of the $\alpha$-basis with built-in models. A stability check ensures the viability of the theory by implementing the no-ghost and positive speed conditions. hi_class is interfaced with MontePython [148, 149] to compute cosmological constraints.

• 3.

COOP is an EB code which solves the linear cosmological perturbations for GLPV theories and sub-classes by implementing the formalism of the $\alpha$-basis. It also includes a stability module verifying the absence of ghost and gradient instabilities. The likelihoods to perform cosmological constraints are embedded within the code itself.

• 4.

EoS_class is a modified version of the public EB solver CLASS [144]. It implements the $\alpha$-basis description of the Horndeski theory [128]. Unlike the previous EB codes, EoS_class is based on the equation of state approach [150]. This formalism encloses the modifications to GR in an effective fluid described by a non-trivial stress-energy tensor and the coefficients multiplying the fluid perturbations are written in terms of $\alpha$-functions. The code also includes as specific theory the designer $f(R)$-model (CLASS_EOS_FR) [151].

Authors in Ref. [137] compared the shapes of the CMB angular power spectrum and of the dark matter power spectrum predicted by EFTCAMB, hi_class and COOP. The results agree at a level of $0.1\%$ for the matter power spectrum at all scales and for the TT, EE, TE spectra of CMB for angular scales $\ell>100$. Deviations of up to $0.5\%$ arise for $\ell<100$ due to known lack of convergence issues already detected when comparing results of CAMB and CLASS in $\Lambda$CDM. The analysis of COOP, instead, shows that COOP achieves the required precision only for $k<1h{\rm Mpc}^{-1}$. EoS_class was cross-checked with hi_class for the $\alpha$-basis [76] and with EFTCAMB for the designer $f(R)$-model [137] showing a sub-percent agreement. These results strengthened the confidence on these numerical codes making them efficient for precision constraints on cosmological and gravitational parameters.

The initial conditions (ICs) for N-body codes can be fixed by using the linear EB codes or using the Zel’dovich’ approximation [152]. Due to the complexity of modeling the physical phenomena approximations are sometimes employed to set ICs, which for example do not account for dynamical perturbations of DE that instead might be significant. In this regard, ICs for N-body simulations were discussed within the EFT framework [153]. The model considered for the analysis is the $f(R)$-gravity, modeled using the designer approach [154], $i.e.$ the expansion history is fixed and the dynamical equation is then solved for $f(R(a))$. In this specific case the expansion history is chosen to closely mimic the $\Lambda$CDM and $w$CDM. As result, ICs set at early time leave imprints up to 5% at Mpc scales. This notably implies that one must go beyond the $\Lambda$CDM ICs for a proper implementation of $N$-body simulations in MG theories. A public code for the generation of ICs exists: FalconIC [153] 999FalconIC webpage: http://falconb.org. It can be linked to any version of both EB codes CAMB and CLASS, including EFTCAMB.

Many of the phenomenological investigations and cosmological constraints we present in the next Sections are obtained using the above linear EB packages.

3 Novel predictions

Ongoing and future cosmological surveys offer an unprecedented insight into gravity on cosmological scales. The EFT framework revealed to be a powerful theoretical tool to systematically identify clear patterns and predictions of MG and DE proposals. In this Section, we review the novel predictions obtained with this approach, namely by using the systematic enforcement of stability conditions and the straightforward computation of cosmological predictions and observables. Let us note that until now the widely investigated class of models within the EFT framework is Horndeski gravity, which encloses a large class of well known scalar-tensor theories. Thus, most of the results presented in this Section apply to Horndeski models even though in some cases they can be straightforwardly extended.

3.1 Gravitational couplings

Cosmological probes can be schematically divided into two subsets: surveys observing the smooth expansion of the Universe, such as SNIa and BAO, will constrain the Hubble rate $H(t)$ or the DE equation of state $w_{\rm DE}$; GC, RSD, CMB and WL surveys scrutinize the clumpy nature of Universe at large scales leading to measurements of gravitational potentials, matter density and temperature fluctuations power spectra. The latter data can be seen as hybrid probes since they contain information both on the evolution of the background and perturbations. These cosmological probes can be used to test DE and MG proposals at cosmological scales. Regarding the evolution of perturbations in scalar-tensor theories, the two powerful and handy phenomenological functions $\mu(t,k)$ and $\Sigma(t,k)$ prove useful to interpret theoretical predictions in light of observations [155, 156, 157, 158, 159]. The former, known as the effective gravitational coupling or the effective Newton constant, characterizes the modifications of gravity on the clustering of matter. The latter, called the light deflection parameter, describes the modifications of gravity on null geodesics, $i.e.$ how light travels on cosmological distances. To define these phenomenological functions, let us consider linear scalar perturbations for which the line element in Newtonian gauge reads

[TABLE]

where $\{\Psi(t,x^{i}),\Phi(t,x^{i})\}$ are the gravitational potentials. The $\mu(t,k)$ and $\Sigma(t,k)$ functions are defined in Fourier space as

[TABLE]

where $G_{N}$ is the Newton gravitational constant and the comoving density contrast is defined as $\Delta_{\rm m}=\delta_{\rm m}+3Hv/k$, where $\delta_{\rm m}=(\rho_{\rm m}-\bar{\rho}_{\rm m})/\bar{\rho}_{\rm m}$ is the density contrast and $v$ is the irrotational component of the peculiar velocity. GC and RSD data, being statistics of the matter perturbation stochastic field, are direct probes of $\mu$. $\Sigma$ measuring the deviation in the Weyl potential ($\Phi+\Psi$) can be probed with measurements sensitive to the lensing of light. A third quantity, although its connection to observations is less obvious, called the gravitational slip parameter, is often considered

[TABLE]

The three phenomenological functions are thus linked by the relation

[TABLE]

Note that the GR limit is recovered when the three phenomenological functions are equal to unity.

Finding analytical forms of any of the above functions is generally not possible but requires numerical solving. However, a direct connection between them and a specific theory of gravity can be derived if the QS approximation is assumed. This approximation amounts to neglecting terms involving time derivatives in the Einstein equations for linear perturbations. It has been proved to be a valid assumption within the sound horizon of the DE mode [134, 135]. The advantage offered by the EFT approach is twofold. On the one hand, it grants the possibility to explore modification of gravity beyond the simplifying assumption that is the QS approximation. In this case $\mu$, $\Sigma$ and $\eta$ can be computed numerically using EB codes which evolve the full linear perturbative equations of the EFT formulation. On the other hand, it allows to obtain explicit and algebraic expressions in the QS limit for each phenomenological function, enabling a direct and neat connection between observables and the relevant couplings of the gravitational interaction, $i.e.$ the EFT functions. In this case the analytical expressions are obtained for large classes of models ($e.g.$ Horndeski, GLPV, Hořava gravity) instead of being derived for any specific theory such as $f(R)$-gravity. As we will see further on, most of the discussions on the phenomenology of $\mu$, $\Sigma$ and $\eta$ in literature deals with Horndeski models and sub-classes. We thus focus on these cases here while we refer the reader to [79, 54, 64, 82, 145] for the expression of the phenomenological functions in scenarios beyond this landscape. For the case of Horndeski models and sub-classes, the analytical expressions of the phenomenological functions in the QS approximation can be written as follows

[TABLE]

where $f_{i}$ and $M_{\mathrm{C}}^{2}$ depend on the EFT functions (see [160] for details). $M_{\mathrm{C}}$ sets a crucial transition scale: the one below which the scalar field mediates a fifth force [5]. A fifth force is characteristic of MG theories with extra DoFs and the transition scale depends on the dynamical mechanism which screens the strength of the scalar fifth force in local environments [5]. Astrophysical scales are typical examples of screened environments. In other words, the transition scale relates to the Compton wavelength of the scalar field, $\lambda_{\mathrm{C}}$, as $M_{\mathrm{C}}\propto\lambda_{\mathrm{C}}^{-1}$. Theories of the chameleon type [161, 162] display a small Compton wavelength $\lambda_{\mathrm{C}}\lesssim 1\;\rm{Mpc}$, whereas models exhibiting self acceleration, and in a more general sense, models where the extra DoF sources cosmic acceleration, bare a very large Compton wavelength $\lambda_{\mathrm{C}}\propto H^{-1}$. Let us look into both regimes in more details.

On super-Compton scales, $i.e.$ $k/a\ll M_{\mathrm{C}}$, which we identify with subscript “${\rm sc}$”, the eqs. (55)-(57) reduce to [160]

[TABLE]

The above are thus representative of the modification of gravity which remain in a screened environment. For scalar-tensor theories such environment amounts to a medium where the scalar field is decoupled from the matter fields. The only mediators left to transmit long-range interactions are the tensors modes.

On the other hand, in the sub-Compton regime, $i.e.$ $k/a\gg M_{\mathrm{C}}$, which we denote with subscript “$\infty$”, the eqs. (55)-(57) yield [160]

[TABLE]

where

[TABLE]

with $\alpha$ being the no-ghost condition and $c_{s}^{2}$ the speed of propagation of the scalar mode. These represent then the full modifications of gravity on large cosmological scales.

The phenomenological functions written as in eq. (3.1) allows a deeper understanding of the origins of modifications of gravity. For example the effective gravitational coupling $\mu$ can be understood as

[TABLE]

where $\mu_{\mathrm{ff}}=\beta_{\xi}^{2}/(1+\alpha_{T})$ characterizes directly the strength of the fifth-force mediated by the scalar field [163]. This contribution must be positive for a viable model since the gravitational interaction induced by a healthy spin-0 field is always attractive. This is indeed the case thanks to the stability conditions: $\{c_{s}^{2},\alpha\}>0$ for scalar modes and $1+\alpha_{T}>0$ for tensor modes.

We must emphasize that care must be taken regarding the normalization of the effective Planck mass. Defining a screened version of the phenomenological functions amounts implicitly to deciding the normalization scheme adopted for the effective Planck mass. Experiments measuring the value of the Newton constant are performed in a screened environment and thus probe the gravitational coupling which remains in the EFT action once the scalar field is decoupled from the gravitational potentials: $G_{\rm sc}=(1+\alpha_{T})/M^{2}$ [163]. Therefore, at present time, $t_{0}$, its value must coincide with the Newton constant, $i.e.$ $G_{\rm sc}(t_{0})\equiv G_{\mathrm{N}}$, which implies $\mu_{\mathrm{sc}}$ must be normalized as [163]

[TABLE]

Upon identification with eq. (58) one deduces $M^{2}_{\mathrm{pl}}=M(t_{0})^{2}/\left(1+\alpha_{T}(t_{0})\right)$. Importantly, this also implies that $\mu_{\mathrm{sc}}(t_{0})=1$. Note that writing a simple and analytical definition of $\mu_{\mathrm{sc}}$ might not be straightforward when the QS approximation is not applied. In this respect, this normalization is not always adopted in literature. In the rest of the review, if not stated otherwise, we assume the convention to set the effective Planck mass today to the Planck mass is the one adopted101010The previous being a matter of definition, one might be led to believe a simple rescaling of results in one normalization into the other to be sufficient to establish a fair comparison. This however has some caveats because of the underlying stability conditions. Let us consider an illustrative example. In the context of Horndeski theories, stability conditions and observations generically push towards $M^{2}(t)\gtrsim M^{2}_{\mathrm{pl}}$. Assuming this strict prior and considering the normalization $M^{2}(t_{0})=M^{2}_{\mathrm{pl}}$, an MCMC analysis would therefore produce models satisfying $dM^{2}(t_{0})/dt\gtrsim 0$ only. On the contrary, if the normalization of $M^{2}$ is left free, provided $M^{2}(t_{0})\neq M^{2}_{\mathrm{pl}}$, both outcomes $dM^{2}(t_{0})/dt\gtrsim 0$ and $dM^{2}(t_{0})/dt\lesssim 0$ would appear and be stable. As a result, the posterior distribution on $M^{2}$, if not all parameters of the analysis, will be different than in the former case. Both the size and the shape thus depend on the normalization and a simple rescaling would not suffice..

3.2 Phenomenology from $\mu$, $\Sigma$ and $\eta$

The phenomenological functions $\mu,\Sigma,\eta$ introduced in the previous Section can be used for a systematic interpretation of observations with the aim of constraining or ruling out classes of DE/MG models. Analyzing their expressions in terms of EFT functions already gives an insight about the phenomenology of models. Furthermore, numerical investigations of these phenomenological functions give also precious information since the planes identified by $\mu$ - $\eta$ and $\mu$ - $\Sigma$ have been proven to be sound benchmarks to highlight peculiar features of DE/MG models. In this regards the use of a Monte-Carlo approach to systematically generate large samples of pure EFT models ($\sim 10^{4}$) under stability requirements enables to draw conclusions of wide applicability.

Figure 1 [160] is a practical example of how model-independent measurements can help in the diagnostic of Horndeski models. Self-accelerating models [164, 165] have a very small scalar mass, comparable to the Hubble rate, then a detection of $k$-dependence by LSS surveys in either of the phenomenological functions would rule these models out. The condition $\mu_{\infty}>\mu_{\rm sc}$ follows from the presence of attractive fifth force in Horndeski models, hence a detection of the opposite, $\mu_{\infty}<\mu_{\mathrm{sc}}$, would rule out the whole class of models. In the super-Compton limit, $\alpha_{T}$ regulates the deviation in the slip parameter. Thus an observation of $\eta_{\rm sc}=1$ or $\mu_{\mathrm{sc}}=\Sigma_{\rm sc}$ would imply $\alpha_{T}=0$. The opposite also holds, a measurement of $\alpha_{T}=0$ requires $\eta_{\rm sc}=1$ otherwise Horndeski theories would be ruled out. GWs are weakly constrained on cosmological scales with CMB but recently the Laser Interferometer Gravitational-Wave Observatory (LIGO)/Virgo and the INTernational Gamma-ray Astrophysics Laboratory (INTEGRAL)/Fermi collaborations constrained the deviation in the speed of GWs to be smaller than $10^{-15}$ [166, 167]. This has severe implications for MG and in particular on the form of $G_{4}$ and $G_{5}$ in the Horndeski Lagrangian [168, 169, 170, 171, 172] (see Section 5.2 for a more detailed discussion). Assuming therefore that $\alpha_{T}$ is negligible, a measurement of $\mu_{\infty}\neq\Sigma_{\infty}$ would give a signature of $\alpha_{M}\neq 0$. On the contrary it would constrain the form of $G_{4}$ and $G_{5}$ in the Horndeski Lagrangian [160].

The QS relations of $\mu$ and $\Sigma$ discussed in Section 3.1 allowed to deduce another important feature characterizing the Horndeski models, the so-called $\mu$-$\Sigma$ conjecture [160]

[TABLE]

which states that measurements of $\mu-1$ and $\Sigma-1$ of opposite signs at any redshift, $z$, and scale would strongly disfavor Horndeski models. A value of $\mu_{\mathrm{sc}}<1$ is predicted by self-accelerating models due to an increasing effective Planck mass. Then one does not expect to observe $\Sigma_{\rm sc}>1$. The latter is extremely unlikely since it would require a large positive $\alpha_{T}$ to change the trend of the effective Planck mass only in $\Sigma_{\rm sc}$. In the sub-Compton regime by comparing $\mu_{\infty}$ and $\Sigma_{\infty}$ from eqs. (3.1), it is straightforward to deduce that the conjecture breaks down if the following inequalities hold [136]

[TABLE]

The sign agreement between $\mu_{\infty}-1$ and $\Sigma_{\infty}-1$ thus depends on the competition between the EFT functions $\{M^{2},\,\alpha_{T}\}$ and the quantities $\{\beta_{\xi},\,\beta_{B}\}$. This physically translates into the competition between screened and fifth force contributions. The failure of the conjecture on sub-Compton scales would then require a significant fine-tuning. The sign agreement in $\mu_{\infty}-1$ and $\Sigma_{\infty}-1$ is more obvious because they have the same pre-factor and the additional contributions from $\alpha_{T}$ and the fifth force are of the same order [160]. The existence of a correlation between $\mu$ and $\Sigma$ was also found using a Monte-Carlo approach to generate a large sample of viable Horndeski model with a fixed $w$CDM background 111111$w$CDM is the cosmological model where $w_{\rm DE}=w_{0}$ and $w_{0}$ is a free constant parameter. and the QS approximation on sub-Compton scales to model $\mu_{\infty}$ and $\Sigma_{\infty}$ [173]. This diagnostic suggests that Horndeski theories would be strongly disfavored if [173]: $i)$ $\mu_{\infty}$ and $\Sigma_{\infty}$ are observed to have opposite sign for $z\gtrsim 1.5$; $ii)$ $\mu_{\infty}<1$ at $z=0$. A complementary analysis allowed to investigate the validity of the conjecture not only at different redshift but also at different scales as shown in Figure 2 [136]. The results showed that the conjecture holds very well within the QS approximation, but the exact behaviors of $\mu$ and $\Sigma$ obtained with EFTCAMB violate the conjecture at $k=0.001$ h/Mpc, where the full dynamics of the Horndeski models allows $\Sigma>1$ and $\mu<1$ [136]. The validity of the conjecture at all times is in tension with the previous work where the sign agreement is lost for $z<1.5$ [173]. The reasons stem from how viable are the models considered: using only stability conditions (ghost and gradient) leads to the violation of the conjecture at small redshifts [173], while requiring additional priors based on observational constraints allows for the conjecture to be restored [160]. In particular, Solar System bounds on the variation of the Newton’s constant translate into a stringent prior on $\mathfrak{f}$ at present time [174]. Importantly, the breaking of the conjecture at $k=0.001$ h/Mpc was shown to occur despite such scale being well within the sound horizon of the scalar field, highlighting the fact that time derivatives of the scalar field and the metric can no longer be neglected. Thus the condition of validity of the QS approximation within the sound horizon of the dark field has to be seen as a necessary condition but not sufficient [136].

A way to go beyond the QS approximation at large scales while keeping analytical expressions of the phenomenological functions is to adopt a semi-dynamical treatment [175], which amounts to consider time derivatives of the metric potentials and velocity fields on large scales [175]. This extension of the QS approximation was designed by evolving the perturbations at a given pivot scale and extrapolating the relations between the perturbations to other scales. Doing so a more precise treatment of large scales is encapsulated while the connection between theory and observables remains thus straightforward thanks to analytical expressions of $\mu$ and $\eta$. They read [175]:

[TABLE]

where $\mu_{\pm i}$ and $\eta_{+i}$ are expressed in terms of EFT functions and corrective terms arising from the semi-dynamical treatment. Corrections near the Hubble scale, $k_{H}=k/aH$, up to order $k_{H}^{6}$ in $\mu$ and up to $k_{H}^{4}$ in the gravitational slip parameter appear. In the small scales limit ($k\rightarrow\infty$), the semi-dynamical expressions for $\eta_{\infty}$ and $\mu_{\infty}$ for Horndeski models match naturally the ones in the QS approximation. However in the same limit, the velocity field and time derivative of the spatial metric potential in GLPV theories should not be neglected. The semi-dynamical expressions differ from the QS ones by factors proportional to $\alpha_{H}$. Thus for these theories a semi-dynamical treatment should be preferred over the QS one even at small scales. The validity of this result is shown in Figure 3 [175].

With the perspective of large amount data from future missions, one might also tailor efficient procedures to reconstruct the time scaling of the phenomenological functions. As such, one can capture the evolution of $\mu$, $\Sigma$ and $w_{\rm DE}$ in joint prior covariance matrices [107]. These matrices can hence be considered as priors which encapsulate the definite trends on $w_{\rm DE}$, $\mu$ and $\Sigma$ of Horndeski theories and sub-classes that are shown in Figure 4. The Monte-Carlo procedure highlighted that for instance these models can display $\mu(z=0)<1$ as soon as the normalization of the effective Planck mass today is left free. In parallel, this analysis shows quantitatively how the redshift evolution of $\mu$ and $\Sigma$ is strongly correlated in scalar-tensor theories, while their correlation with $w_{\rm DE}$ is progressively washed out when going beyond GBD models. The correlations to use in future model reconstructions can be effectively mapped into the CPZ parameterization [176] where for $\mu$ and $\Sigma$ they are shown to scale as $|a-a^{\prime}|$ and for $w_{\rm DE}$ as $|\mbox{ln}a-\mbox{ln}a^{\prime}|$, where $a$ and $a^{\prime}$ are two bins in the scale factor [107].

3.3 Impact of stability conditions

In the EFT formulation, the conditions which guarantee the stability of a theory are expressed in terms of the EFT functions. This allows to gauge straightforwardly whether a gravitational model is physically viable. We have discussed in Section 2.5 the required physical stability conditions a model has to satisfy, namely, the avoidance of ghosts, gradient and tachyonic instabilities. Here, we review the implications of such conditions on the parameter space of models.

From this stability perspective, the “naturalness” of the $\Lambda$CDM scenario among Horndeski theories can be addressed [85]. In this case the Horndeski class of models is parametrized within the EFT framework using the pure EFT approach. We refer to them as pure Horndeski models hereafter. As displayed in Figure 5, the ghosts and gradient stability conditions shrink the allowed space of EFT parameters up to the point where $\Lambda$CDM corresponds to a corner of the stable region. This happens quite naturally since $\Lambda$CDM corresponds to a non-propagating extra DoF which follows from a vanishing kinetic term and speed of propagation. In other words, the realm of Horndeski theories with an acceleration of the background expansion compatible with current observational constraints has $\Lambda$CDM at its border. Interestingly, this analysis also revealed the parameter space location of models which violate the null energy condition in a stable manner, $i.e.$ “super-accelerating” (or “phantom”) models with $w_{\rm DE}<-1$. The volume of the parameter space occupied by the stable pure Horndeski models gradually shrinks as the background equation of state decreases (see Figure 5). In particular, the $\Lambda$CDM model does not fall within the set of phantom models, highlighting the fact that at least one non-minimal coupling must be invoked for models to violate the null energy condition in a stable way. In parallel, there is one EFT function which does not play a role in the expression of $\mu$ nor $\eta$ when the QS approximation is considered: $M_{2}^{4}$ or equivalently $\alpha_{K}$. This coupling acts solely on the ghost stability condition and regulates hence the speed of sound of the scalar perturbations. Figure 5 shows also how large values of this coupling open up the space of stable models and thereby helps to oppose the shrinking effect due to considering $w_{\rm DE}<-1$. This feature allowed to identify a large set of viable null energy condition violating theories.

The impacts of ghost and gradient conditions on GBD class of models were investigated using the pure EFT approach with an effective Planck mass defined by the so-called “linear EFT model” [73]:

[TABLE]

where $\Omega_{0}^{\mathrm{EFT}}$ is a constant. As a consequence, all MG effects are captured by the parameter $\Omega_{0}^{\mathrm{EFT}}$ and the constant DE equation of state parameter $w_{\rm DE}$ defining a $w$CDM expansion history. This analysis was the first to highlight how the stability conditions act as a strong prior on the constraints of cosmological parameters. In the top left panel of Figure 6 one can notice for instance that a stable linear EFT model implies $\Omega_{0}^{\mathrm{EFT}}>0$. Once the effective equation of state parameter becomes smaller than $-1$, the likelihood contours produced by a statistical analysis would therefore be significantly tightened. This result independently confirms the shrinking theory space effect in the super accelerating regime discussed previously [85].

Another specific model with a running Planck mass is $f(R)$-theory. In this case, one can study the impact of stability conditions on the parameters of the model using a designer approach [154]. The resulting family of models are parametrized by the present value of the Compton wavelength of the scalaron, $f_{R}=df/dR$, [177], defined as [178]

[TABLE]

As shown in the top right panel of Figure 6, the stability conditions which include the no-ghost, no-gradient and a positive mass ($m^{2}=f_{RR}>0$) induce a stringent cut of the parameters space [73]. One can note that as opposed to the linear EFT model, the designer $f(R)$ model does not allow stable regions beyond $w_{\rm DE}<-1$ as expected theoretically, because this amounts to a minimally coupled scalar field for which violating the null energy condition can only be at the expense of it being a ghost. These viable parameter spaces have been used to put stringent constraints on the parameters of both linear EFT model and $f(R)$-theory showing in some cases to have a higher constraining power than that of data [73]. We will further discuss this point in Section 4.

The selection of the viable parameter space induced by the stability conditions has non trivial implications for EFT predictions. We already discussed in Section 3.2 the $\mu-\Sigma$ conjecture [160] and its implications [173, 136]. Other general behaviors can be identified using a Monte-Carlo approach to generate a vast sample of healthy pure Horndeski models and explore the phenomenology of $\mu$, $\eta$, $\Sigma$ and $f\sigma_{8}$ across $z$. $f\sigma_{8}$ is the so-called growth function, $i.e.$ the product of $\sigma_{8}$ and the linear growth rate of matter density perturbations $f$ (see eq. (70)). For instance, definite features can be found for models exhibiting sub-luminal propagation of scalar and tensor perturbations and DE contributes only at low redshifts [163]. As shown in Figure 7, the effective gravitational coupling displays a characteristic $S-shape$ pattern, a regular alternating succession of epochs where gravity is stronger and weaker than predicted by GR [163]. This behavior significantly affects the growth function $f\sigma_{8}$. The amplitude of this observable is generically suppressed, compared to the value expected in $\Lambda$CDM models, at intermediate redshifts ($0.5\lesssim z\lesssim 1$), the opposite being true at all other cosmic epochs, when a phenomenon of super growth is highlighted. The gravitational slip parameter $\eta$ is, instead, predicted to be bounded from above. Notably, for GBD theories, the gravitational slip parameter is at most unity and a larger deviation from GR is found only at intermediate redshifts for the Horndeski case.

Historically only the no-ghost and no-gradient conditions are employed when exploring the stability of a gravity theory in literature. As discussed in Section 2.5, such requirements are high-$k$ dependent statements thus they cannot guarantee stability on the whole range of cosmic scales. In EB codes, to prevent exponentially growing modes at low-$k$, ad-hoc mathematical conditions are implemented to complete the set of physical conditions. In the EFT formulation, these additional requirements are worked out at the level of the perturbation equation for the scalar field $\pi$ [72, 86]. The stability module of EFTCAMB includes the tachyon conditions as a function of the mass eigenvalues $\mu_{i}$, $i.e.$ the mass condition [179]. The evaluation of the impacts of switching on this additional condition reveals the mass condition to be very efficient for substituting the mathematical condition in practically all cases of cosmological interest [179]. We display the results obtained thanks to the Monte-Carlo generation of healthy pure Horndeski and GBD models in the $\mu-\Sigma$ plane at $z\approx 0.1$ in Figure 8 121212Note that in this investigation the effective Planck mass today was not normalized to $M^{2}_{\mathrm{pl}}$.. One can notably observe that as opposed to the general case of pure Horndeski, the action of the mass condition becomes very peculiar in the pure GBD case and cuts out models from the $\Sigma\,,\,\mu>1$ region at $z\approx 0.1$. The tail of models that instead are present when only the no-ghost, gradient and mathematical stability were active are thus singled out. This cut translates directly in a lower bound of the effective Planck mass which has to be positive at all redshifts. From these investigations, the mass requirements together with the no-ghost and no-gradient conditions prove a reliable and complete set of robust and physically motivated theoretical priors that guarantees the soundness of the theories on all linear scales.

A common thread linking a lot of the analyses presented throughout this review is the processing of the background expansion rate as a fixed prior when studying the effects of MG theories in the perturbed sector. As anticipated in Section 2.2, a complementary approach can be investigated. Using the free background approach, where $H(t)$ is computed from the time evolving $w_{\rm DE}$ expressed in terms of the background EFT functions $\Lambda$ and $c$, reveals that universal behaviors could also be found on quantities characterizing the background expansion [106]. The Monte-Carlo procedure employed shows that the stability conditions favor the emergence of a tracking behavior in the scaling of $w_{\rm DE}$ of Horndeski models with the value of this function becoming close to zero deep in matter domination while approaching the value $w_{\rm DE}=-1$ at present time (see Figure 6). This renders the possibility to map this evolution onto the Chevallier-Polarski-Linder (CPL) parameterization [180, 181], $i.e.$ $w_{\rm DE}=w_{0}+w_{a}(1-a)$ where $w_{0}$ and $w_{a}$ are constants. The same exploration technique for quintessence models reveals a tracker behavior with however $w_{\rm DE}<-1$ excluded by stability, as expected theoretically. The definite features of the background evolution are translated into a theoretical prior covariance matrix for $w_{\rm DE}$, correlating its values at different redshifts to simplify the reconstructions of $w_{\rm DE}$ from future data. This has been recently used also for the reconstruction of $w_{\rm DE}$ through Gaussian process in quintessence and Horndeski theories [182].

Stability conditions are usually used on top of a chosen parameterization for the EFT functions. However it is possible to construct a stable EFT basis which evades ghost and gradient instabilities or accommodates further theoretical priors such as a luminal or sub-luminal scalar propagation speed [183, 184, 185]. Such basis instead of using the usual EFT functions (or $\alpha$-basis) employs directly the stability functions $\alpha,M^{2},c_{s}^{2}$ with a boundary condition $\alpha_{B0}=const.$. An advantage of this basis is that it avoids $\Lambda$CDM to be confined to a narrow corner of the stability space. This corner would be difficult to sample and could lead to false evidence against the $\Lambda$CDM model. The enforcement of theoretical stability conditions when performing cosmological constraints is definitively an advantage as they ensure the resulting observational bounds lay in a physical parameter space. This is now common practice and most of the available codes have built-in modules which automatically check for the stability of the considered model. A further benefit in using the stable EFT basis would be secured when running MCMC algorithms: Markov chain would no longer randomly explore unstable regions of the parameter space making the whole process much faster. We note however that the stable EFT basis does not include the tachyonic condition discussed in Section 2.5, thus the stability of the underlying theory of gravity might not be guaranteed at all scales.

3.4 Linear growth rate of matter density fluctuations

The EFT formulation contributed to understand how linear structures are likely to grow in scalar-tensor theories of gravity and what theoretical patterns a viable theory should have to better accommodate cosmological observations. RSD, GC and WL independent measurements all consistently detected a lower growth rate of matter density perturbations than predicted by the $\Lambda$CDM model. This is at variance with results extrapolated from CMB studies [14, 13, 15, 16, 186, 187, 17]. Scalar-tensor theories are not expected to produce lower growth relative to $\Lambda$CDM because of the fifth force they exhibit. The EFT has contributed to showing that this is not necessarily the case [85, 188, 163].

The effective gravitational coupling modifying the Poisson equation (52) alters the growth of matter perturbations $\delta_{m}$ as follows 131313Let us note that in GLPV models, an additional modification proportional to $\alpha_{H}$ is present in the friction term in eq. (69). This is due to the fact that terms involving matter velocity cannot be neglected on sub-horizon scales, see for instance [145].

[TABLE]

The linear growth rate $f$ is then defined by [189]

[TABLE]

In the case of Horndeski models, using the sub-Compton limit one can deduce from eq. (61) that in order to have lower growth the condition $\mu_{\infty}<1$ has to be satisfied. One necessary condition is to require $\mu_{\mathrm{sc}}<1$ which originates only from the tensor part, in particular it depends on the ratio between $c_{t}^{2}/M^{2}$, $i.e.$ the tensor speed and the effective Planck mass ($M^{2}$) [188]. However, the $\mu_{\mathrm{sc}}<1$ condition is not sufficient to guarantee a weaker gravity because of the interaction between the scalar field and matter which always enhances $\mu_{\infty}$, $i.e.$ the fifth force. Numerically the growth of structure in Horndeski models was investigated using a large samples of Monte-Carlo simulated healthy Horndeski models and sub-classes with sub-luminal propagation of scalar and tensor perturbations [163]. The results are displayed in Figure 7 where it is clear that the majority of models satisfy $M^{2}>M^{2}_{\mathrm{pl}}$ and $c_{t}^{2}<1$, hence producing lower growth relative to $\Lambda$CDM at low redshifts. At higher redshifts instead it is possible to note a period of super-growth with a linear growth rate larger than that predicted in the standard model. GLPV theories are even more flexible in describing the time dependence of the linear growth rate because $\alpha_{H}$ provides extra freedom on top of $M^{2}$ and $c_{t}$ to modulate the stronger/weaker gravity pattern [188]. In particular, for GLPV it is possible to obtain a suppression in the matter power spectrum with respect to $\Lambda$CDM [145, 82].

In the standard model of cosmology, the linear growth rate can be efficiently parametrized as $f=\Omega_{\rm m}(a)^{\gamma(a)}$ [189] where the parameter $\gamma(a)$ is the so-called growth index [190]. MG models can predict a slight deviation from $\Lambda$CDM parameterization yet detectable. For example, [191] showed that the growth rate of a large class of MG models can be accurately described by

[TABLE]

where $\gamma_{0}$ and $\gamma_{1}$ are parameters which depend, in a predictable way, on the adopted gravitational theory. The EFT basis can be directly mapped onto the $\gamma_{0}-\gamma_{1}$ parameterization [85]. Thus once a specific parameterization of the EFT functions is chosen and a background is fixed one can obtain constraints on $\gamma_{i}$ parameters. This phenomenological approach allows to gain further insights on the mechanisms that tend to suppress the growth of structures on large scales. Assuming a $\Lambda$CDM background evolution and no-ghost and no-gradient conditions, the space of theories leading to weaker gravity in the $\gamma_{0}-\gamma_{1}$ plane is found to be much larger than that of stronger gravity [85] (see Figure 9). In other words, statistically, Horndeski theories are more likely to produce models which exhibit weaker gravity, in the local Universe, as compared to $\Lambda$CDM. In particular, once the cosmic expansion rate is fixed to mimic $\Lambda$CDM, no viable theory can show a value for $\gamma_{0}$ which is larger than that of standard scenario [85]. Figure 9 displays also forecast for a Euclid like survey in this approach.

The tendency of viable Horndeski theories to produce lower growth relative to $\Lambda$CDM is not necessarily limited to late-time DE (LDE) models but it is also present in scenarios where DE contributes throughout matter domination [173]. The former are models for which the GR limit at early times is recovered, the latter are distinguished between early DE (EDE) models, $i.e.$ models with the effective Planck mass contributing to the total energy moment tensor even at early time, and early MG (EMG) models, $i.e.$ models where on top of the previous the EFT function do not vanish at early times. These scenarios are constructed using different asymptotic behaviors of the EFT functions. The $f\sigma_{8}$ - $\Sigma_{\infty}$ plane proves to be instrumental in discriminating between the aforementioned behaviors embedded within Horndeski theories and their capability to generate a lower growth. The diagnostic reveals that LDE models are strongly disfavored if $f\sigma_{8}<(f\sigma_{8})_{\Lambda CDM}$ at $z\gtrsim 1.5$; EDE models if $f\sigma_{8}>(f\sigma_{8})_{\Lambda CDM}$ at $z\gtrsim 1.5$ or, simultaneously, $f\sigma_{8}<(f\sigma_{8})_{\Lambda CDM}$ and $\Sigma_{\infty}>1$ at $z\gtrsim 1.5$. The only possibility to display a large enhancement of growth is to use the EMG model. The redshift evolution of these models in the $f\sigma_{8}$ - $\Sigma_{\infty}$ plane is displayed in Figure 10.

Within the class of Horndeski theories, a specific phenomenological model known as the No Slip Gravity model [192] produces lower growth than $\Lambda$CDM at low redshifts as shown in Figure 11 (left panel). The model is characterized by a gravitational slip parameter equal to unity, an unmodified speed of propagation for tensor modes and by the relation $\alpha_{B}=-2\,\alpha_{M}$. As consequence the effective gravitational coupling $\mu$ is lower than GR in the past inducing a lower growth of structure.

On the other hand, suppose future observational data do not allow: $i)$ lower growth of structure relative to the standard model; $ii)$ a gravitational slip different from unity; $iii)$ changes in the damping ($\alpha_{M}\neq 0$) and speed ($\alpha_{T}\neq 0$) of GWs. This does not induce all alternative scenarios to GR to be ruled out. Indeed, modifications of gravity due to a non-vanishing braiding function are still possible. An example of such a model is the No Run Gravity model characterized only by the $\alpha_{K}$ and $\alpha_{B}$ functions [193]. Consequently, the deviation in the effective gravitational coupling depends only on the braiding term. In particular, stability conditions allow only positive values of $\alpha_{B}$ and the model therefore always predicts an enhancement in the growth of structures relative to $\Lambda$CDM as depicted in Figure 11 (right panel) [193].

3.5 Imprints on cosmological power spectra

Non-standard gravitational scenarios are best inspected by analyzing the shape of large scale observables, such as the matter density power spectrum, the CMB angular power spectrum, the lensing spectrum, lensing B-mode contribution and the bi-spectrum. They can for instance change the lensing potential ($\Phi+\Psi$) when additional perturbative terms are included [194], change the growth of structure modifying the Poisson equation and affect the shape of the temperature-temperature CMB power spectrum at low multipole-$\ell$ through the integrated Sachs-Wolfe (ISW) effect sourced by $\dot{\Phi}+\dot{\Psi}$ [195, 196], shift its high-$\ell$ peaks due to a modified expansion history [197] and change the ratio between odd and even peaks in models where DE is coupled to dark matter [198].

The EFT approach allows thus to investigate the phenomenology of these effects at linear scales for large class of models [199, 61, 112, 200, 74, 82, 63, 201, 202, 135, 203, 145, 204, 76]. In this case the connection with a specific operator in the action (5) or a physical effect enclosed in the $\alpha$-basis can be easily identified by switching on/off single EFT functions per time. In the following we report on some of them. We note that the results in literature might depend on the chosen parameterization. Values of the effective Planck mass larger than $M_{\mathrm{pl}}^{2}$ are shown to suppress the lensing power spectrum with respect to the standard scenario and to modulate the shape of the TT power spectrum at low multipoles $l\lesssim 30$ giving rise to an enhanced ISW tail [112]. The deviations in the speed of propagation of GWs ($\alpha_{T}$) are responsible for changing the location of the inflationary peak of the BB spectrum [199]. The kineticity coupling, $\alpha_{K}$, has been found to modulate the low-$\ell$ CMB TT power spectrum due to the late-time ISW effect, however this effect is unmeasurable being dominated by the cosmic variance [135]. A non-vanishing positive braiding term, $\alpha_{B}$, leads to both an enhanced lensing auto-correlation function and matter power spectrum and, depending on its magnitude, it can generate either a suppressed ISW tail or an enhanced one [76]. The GLPV function, $\alpha_{H}$, is responsible for the damping of the matter power spectrum. Increasing positive values of $\alpha_{H}$ enhances the CMB TT power spectrum at low-$\ell$ [145, 204, 82] as it can be noticed in Figure 12 (left panel). On the same line the lensing potential decreases as a function of $\alpha_{H}$ [145, 82]. For the particular case of the No Slip Gravity model, the running of the Planck mass, $\alpha_{M}$, is shown to mildly affect B-mode reionization and recombination bumps on low multipoles $l\lesssim 10$, while at $l\gtrsim 100$ its impacts on the lensing B-mode polarization are more significant (of order four times the maximum value of $\alpha_{M}$) [201], see Figure 12 (right panel).

Given the ever-increasing precision and scales probed by surveys, integrated effects on ultra-large scales such as ISW, Doppler, Shapiro time-delay, etc, might soon enter the observational window. It is thus interesting to explore their effective potentiality in probing departures from standard gravity. Such effects have been scrutinized within the EFT formulation for Horndeski and GLPV models [200, 204]. For instance, the amplitude of the ISW signal when forecasted for Horndeski theories may deviate by up to $\mathcal{O}(1000\%)$ from standard GR expectations, while the same models only induce $\mathcal{O}(10\%)$ modifications of the amplitude of local observable such as $f\sigma_{8}$ [200]. In parallel, taking into consideration the lensing convergence or the Shapiro time delay leads to additional contributions highlighting different sensitivities to the EFT functions as displayed for the galaxy number count angular power spectrum in Figure 13 (top panel) [200]. The ISW contribution remains the most sensitive probe to MG effects, yet including other integrated effects can therefore improve the total sensitivity of cosmological spectra to measure any gravity departures from GR. Note that the contribution of each ultra-large scale relativistic effects changes across redshift. For instance, the Doppler effect is shown in Figure 13 (bottom panel) to be significant at all epochs in GLPV theories while the ISW, the time-delay and velocity potential effects become significant for $z\gtrsim 3$ [204].

4 Cosmological constraints

In this Section, we review the cosmological constraints on DE/MG models obtained using the EFT approach. The key characteristic of these analyses is the use of MCMC procedures or Fisher forecasts to constrain gravity at large scales combined often with the theoretical priors exposed in the previous Section. We distinguish between models described by a pure approach and specific gravity models implemented in EB codes using the mapping approach (see definitions in Section 2.8). We include in Appendix B the list of the pure EFT parameterizations used in this Section and tables summarizing their observational constraints in Appendix C.

4.1 Running Planck mass

Modifications of gravity induced by a running Planck mass have several distinct effects on observables. In the scalar sector, a Planck mass varying across time has been shown to impact, for instance, the background evolution, the growth of matter perturbations, lensing potential and the ISW tail of the CMB power spectrum. In the tensor sector, it alters the friction term in the GW equation (eq. (74)) affecting the amplitude of the primordial polarization peak in B-modes. We reviewed these features in Section 3. In the EFT formalism, the running of the Planck mass is encoded by $\alpha_{M}$ in the $\alpha$-basis and a linear combination of the derivatives of $\mathfrak{f}(t)$ and $\bar{M}^{2}_{3}(t)$ in the EFT basis (see Section 2.6 for details).

Generally, it is not possible to consider only $\alpha_{M}$ and set all other EFT functions to zero. This choice would give a vanishing kinetic term (see eq. (60)) leading to ghost instabilities and strong coupling problems. One way to proceed is to tweak the function $\alpha_{K}$ and fix it to a positive value to prevent the appearance of any pathological behavior. Note that in such a configuration any noticeable effects on cosmological observables are only due to $\alpha_{M}$ because $\alpha_{K}$ has no detectable impact [205, 206, 135]. We comment further about this aspect in Section 4.2.

One choice for the form of $\alpha_{M}$ often adopted is to model it as proportional to the DE density parameter, $i.e$ $\alpha_{M}=\alpha_{M,0}\Omega_{\rm DE}(a)$, where $\alpha_{M,0}$ is a constant (eq. (193)). Using this model with a $\Lambda$CDM background, the combination of Planck CMB data and the $H_{0}$ prior from Riess et al. [7] is found to favor a positive $\alpha_{M}$ at $\sim 2\sigma$, however, such preference is weakened when geometrical probes are included [75].

The effects of a running Planck mass have been further investigated when a specific link with the braiding function $\alpha_{B}$ is retained. In particular, the relation $\alpha_{M}=-\alpha_{B}$ is typical of conformally coupled models such as $f(R)$-gravity, Brans-Dicke and chameleon theories. Investigations of this class of models has been performed by the Planck collaboration in 2015 [208]. There, $\alpha_{M}$ is modeled as a power law in the scale factor with amplitude $\alpha_{M,0}$ and scaling $\beta$ with a $\Lambda$CDM background. The combination of datasets involving CMB, WL, BAO, RSD measurements leads to stringent constraints on the amplitude of the running Planck mass, $\alpha_{M,0}<0.097\,\,(95\%\mbox{C.L.})$, while the statistical power to constrain its scaling is much weaker, $0.92^{+0.53}_{-0.25}\,\,(95\%\mbox{C.L.})$. For this model, WL data prefer higher values of the expansion rate than CMB ones and because of that WL data lead to weaker constraints [208]. Using the latest Planck data (2018), the results favor $\alpha_{M,0}$ to be negative when only the CMB data are employed, while the inclusion of WL, BAO, RSD measurements give the bounds $\alpha_{M,0}=-0.015^{+0.019}_{-0.017}$ at $68\%\mbox{C.L.}$ (with CMB lensing) [1].

Following on the constraints where a relation between $\alpha_{M}$ and $\alpha_{B}$ is imposed, No Slip Gravity models are representatives with $\alpha_{B}=-2\alpha_{M}$ [192]. For this class of models, stability conditions require $\alpha_{M}\geq 0$. Using the e-fold form eq. (198) to parameterize alternatively $\alpha_{M}$ or $M^{2}$ and a $\Lambda$CDM background, Planck 2015 data favor values of the running Planck mass much lower than what found with LSS data [192, 201]. This is justified by $\alpha_{M}>0$ which induces a higher CMB lensing power and the loss of power on large scale in the matter power spectrum. A time dependent DE equation of state allows to explore larger portions of the parameter space. In this case, for example, the running Planck mass can even assume negative values [207]. CMB, BAO, SNIa and RSD datasets constrain indeed $\alpha_{M}$ to be negative as shown in Figure 14. Remarkably a negative running of the Planck mass at $z\gtrsim 1$ does not prevent the model to produce lower growth relative to the standard model nor a positive ISW-galaxy cross-correlation [207].

The EFT basis (action (5)) offers a natural environment for constraining the running of the Planck mass as well. For example using the linear EFT model eq. (67), the cosmological constraints shown in Figure 15 yield the bound $\Omega_{0}^{\mathrm{EFT}}<0.061$ at 95% C.L. for a $\Lambda$CDM background and $\Omega_{0}^{\mathrm{EFT}}<0.058$ for a $w$CDM background [73], in agreement with what is found later by the Planck collaboration (2015) [208]. The linear EFT model has also been explored to gauge whether this proposal could reconcile the tension between Planck CMB estimation of the lensing amplitude parameter $A_{L}$ and the result from the lensing reconstruction [209]. Simulations of CMB anisotropy and CMB lensing spectra, assuming Planck 2015’s best-fit values and Planck blue book on beam and noise specifications, show that models with an effective Newton constant stronger than $G_{\mathrm{N}}$ can have a modulating effect similar to that of $A_{L}$ [210]. Nevertheless, this induces higher values of $\sigma_{8,0}$ making the tension with WL surveys more severe [13, 15, 14, 16, 186, 187]. The $A_{L}$ tension remains an open issue as confirmed by the latest Planck 2018 results [1].

4.2 Pure Horndeski and GLPV models

The pure Horndeski and GLPV models are characterized in the EFT language by a selected set of free functions of time beyond $H(t)$. For pure Horndeski models they are $\{\mathfrak{f},\bar{m}^{3}_{1},M^{4}_{2},\bar{M}^{2}_{2}\}$ in the EFT basis and $\{\alpha_{M},\alpha_{K},\alpha_{B},\alpha_{T}\}$ in the $\alpha$-basis. An additional function is required to select pure GLPV models, i.e $\mu_{1}^{2}$ or $\alpha_{H}$ depending on the basis. In this description, instead of choosing the functional form of the $G_{i}$ or $F_{i}$ functions respectively in the Horndeski Lagrangians [211] and GLPV ones [37], one can fix the form of each EFT function and describe the linear part of the theory. These sets of functions allow to explore and constrain the primary features of Horndeski and GLPV theories in a pure EFT fashion. It is noteworthy to say that choosing an appropriate model for the EFT functions is a challenging task. Sometimes the chosen forms of the EFT functions result in simplified behaviors of relevant physical quantities, such as $\mu(t,k)$ or $\Sigma(t,k)$, if compared to the complex behavior they have when the mapping approach is considered [212, 213]. The consequence is that one might underestimate the modification of the gravity force or even miss its signatures. In this regard one has to be careful in generalizing the results obtained in the EFT framework recalling that they might be strictly related to the chosen parameterization. On the other hand, it is worth to note that relevant results have been obtained and common trends have been identified using this approach with different parameterizations. We review them in the following.

The first observational constraints of pure Horndeski models have been derived in the $\alpha$-basis, where each $\alpha$-function is parametrized using the linear-de form (eq. (193)). The constraints on the $\alpha_{i}$ parameters using cosmological datasets, such as CMB, BAO, RSD and the power spectrum of galaxies from the WiggleZ survey, are shown in Figure 16 [205]. Notably, the data favor models with the effective Planck mass $M^{2}$ larger than the Planck mass and stability conditions induces a hard prior on the initial value of the effective Planck mass, $i.e.$ $M^{2}_{\rm ini}\geq M^{2}_{\mathrm{pl}}$. The constraints display a neat preference for a positive value for the braiding, negative for the running Planck mass and a sub-luminal propagation of tensor modes. This parameterization was further investigated in light of KiDS+GAMA data considering luminal propagation of tensor modes [214]. While the constraints on cosmological parameters are compatible with $\Lambda$CDM, a preference for positive values of $\alpha_{M}$ and $\alpha_{B}$ is found. Furthermore, the clustering quantity $S_{8}=\sigma_{8}\sqrt{\Omega_{\mathrm{m}}/0.3}$ is highlighted to be better in agreement with the Planck estimate when considering Horndeski theories than $\Lambda$CDM. Just as considering stability or sub-luminal priors improves observational constraints significantly, so does the consideration of “positivity bounds” [215]. The latter arise from requiring basic principles such as a unitary, causal, local UV completion [216, 217, 218, 219]. Such positivity bounds when applied to Horndeski theory imply indeed some constraints on the $G_{i}$ functions which can be rewritten in terms of EFT functions [215]. The inclusion of these additional bounds as theoretical priors when performing parameter estimation analysis, led to a reduction of over 60% of the allowed parameter space for the sub-classes Horndeski models parametrized with the linear-de form (eq. (193)). The additional constraint of sub-luminal propagation of GWs in the same models narrowed down further the viable parameter space to less than a 1%.

A parameterization of the EFT basis also in terms of linear-de form (eq. (193)) has been considered with a $\Lambda$CDM background. For this model [112], the constraining power of viability priors with CMB data are put in perspective. Setting the kineticity to zero so as to be in the most restrictive setup in terms of stability, the posterior distributions of the EFT parameters are understood to be mostly driven by the stability priors (see Figure 2 in Ref. [112]). On top of this, asking for sub-luminal propagation of scalar and tensor perturbations reduces drastically the marginalized contours of the parameters (see Figure 17 top panels). Despite the different parameterization with respect to the model previously discussed [205, 214], the tendency of the data to favor a lower effective Planck mass and a lower propagation speed for tensor modes is recovered. Furthermore, with the aim to be more general, one can be led to increase by one the number of free parameters characterizing each EFT function. This is done expanding in powers of $\Omega_{\rm m}-\Omega_{{\rm m},0}$ for instance and retaining up to two free parameter per EFT function following the de-1 parameterization (eq. (194)). This has the consequence of loosening the posterior constraints as depicted in Figure 17 (middle and bottom panels). From this figure it is also clear that the sub-luminal propagation of perturbations can impact the data constraints and even limit the loosening of the constraints when more free parameters are considered.

Constraining pure Horndeski models with a large set of cosmological probes tends to show that CMB data is only next to RSD data in terms of constraining power [206, 220]. Within CMB data itself, the constraints on the EFT functions are mostly driven by the ISW effect while the inclusion of RSD has the important effect of breaking degeneracies between $\alpha_{M}$ and $\alpha_{B}$ (when $\alpha_{T}=0$) [220]. This result is illustrated in Figure 18 [220] where the constraints are obtained using the scaling-a parameterization (eq.(196)) of the $\alpha$-functions. Assuming luminal propagation of tensor modes, RSD data sets tight bounds on the evolution of the running of the Planck mass disfavoring large positive values. Interestingly, the integrated components of the CMB probe tend to provide more sensitivity to the scaling of the EFT functions while combining with RSD data ameliorates the stringent bounds on their present value. When $\alpha_{T}\neq 0$, data prefer sub-luminal values for the speed of GWs and a larger viability region for the running Planck mass is possible. The intricate link between $\alpha_{M}$ and $\alpha_{B}$ can be further constrained when considering the radiative stability of Horndeski theories [221].

Focusing further on RSD constraints on MG, the recent release of measurements of $f$ and $\sigma_{8}$ separated thanks to galaxy-galaxy lensing by VIPERS [222] and SDSS [223] give constraints on pure Horndeski models as competitive as the full set of $f\sigma_{8}$ data [174]. The models investigated display $c_{t}^{2}=1$ and are modeled using the de-1 form (eq. (194)) with a background expansion set to $\Lambda$CDM. The combination of this new set of data with the full $f\sigma_{8}$ recollection yields a gain on the precision of the constraints on the EFT parameters of at least 20% with respect to the $f\sigma_{8}$ set alone. In parallel, a stringent bound from Solar System tests on the evolution of the Newton constant, $|\dot{G}_{\rm N}/G_{\mathrm{N}}|<0.002\;H_{0}$, is shown to translate into a prior on the unscreenable contribution to gravity, $\mu_{\mathrm{sc}}$, or equivalently the running Planck mass, which in turns increases the constraining power of $f$ and $\sigma_{8}$. In parallel, this prior stands as drastic cut on the weaker gravity the effective Planck mass can induce. Nevertheless, pure Horndeski theories are still able to produce lower growth relative to $\Lambda$CDM with $\sigma_{8,0}$ found to be lower than the Planck $\Lambda$CDM cosmology at more than 2$\sigma$.

Regardless of the parameterization chosen, the kineticity $\alpha_{K}$ or equivalently $\gamma_{1}$ (see Appendix B) has been shown to have a weak impact on observables. This coupling is thus mostly unconstrained even though it bears a fundamental role for the viability of the model considered [205, 206, 135, 174]. The kineticity does not affect the constraints on the other parameters even when considering a large combination of data sets as displayed in Figure 16 [205]. The difficulty in constraining the kineticity coupling has been further investigated [206, 135]. It has been shown in particular that it cannot be constrained by data directly because its contribution to observables is below the cosmic variance [135]. In principle, a way to overcome the cosmic variance limitation could be to use sophisticated multi-tracer techniques [224]. Despite this feature, the kineticity coupling has a significant role in defining the stable parameter space as it enters in the definition of the no-ghost condition for the scalar sector and thus regulates the speed of propagation of scalar modes. As consequence, it changes the viable parameter space used in the MCMC explorations and it cannot be simply discarded [206, 135]. Furthermore, bounds on $\alpha_{K}$ or $\gamma_{1}$ can be obtained using on top of stability conditions the prior on $G_{\mathrm{N}}$, $f$, $\sigma_{8}$ and $f\sigma_{8}$ data sets [174]. These significantly restrict the space of viable models which in turn put indirect constraints on $\alpha_{K}$ or $\gamma_{1}$.

Moving beyond the Horndeski landscape, constraints on pure GLPV models with $c_{t}^{2}=1$ are derived using the linear-de form (eq. (193)) with both $\Lambda$CDM and CPL backgrounds [145]. A combination of CMB, BAO and RSD datasets shows that the GLPV parameter $\alpha_{H}$ is degenerate with the braiding parameter $\alpha_{B}$ and the running Planck mass $\alpha_{M}$ but not with either the cosmological parameters nor the CPL ones. The marginalized distributions of $\hat{\alpha}_{H}$ (the constant parameter of $\alpha_{H}$) excludes the GR limit ($\hat{\alpha}_{H}=0$) in most of the cases analyzed, the data constrain $\hat{\alpha}_{H}$ to be of $\mathcal{O}(1)$ and favoring generally positive values. However, no statistically significant preference of this pure GLPV model over $\Lambda$CDM is found.

The previous constraints are de facto parameterization dependent. Interestingly, it is possible to minimize such dependence by exploiting data-driven reconstruction techniques. In these approaches indeed the number of parameters increases largely, but the constraining power of data can be explored more faithfully. Such an approach is used to constrain Horndeski models and sub-classes (GBD, Quintessence, K-essence and Kinetic Gravity Braiding) by making use of the EFT formulation [225]. The EFT functions are modeled on a fixed time grid and interpolated thanks to piece-wise fifth order spline. The information from the specific theories are encoded thanks to a Gaussian smoothing kernel defined by the CPZ correlation [226]. Doing so, the EFT functions identifying the models are reconstructed across cosmic times using cosmological data and could thus be used to derive specific model properties. This data-driven approach notably highlights that the constraining power of present cosmological probes such as CMB, WL, BAO, SNIa and local measurements of $H_{0}$, is not only high for Horndeski but is enhanced as the complexity of the theory considered increases. In other words, the more EFT functions are considered, the more freedom is allowed, the more the data can express their constraining efficiency. In parallel, the analysis also shows that Horndeski theories might alleviate the tension in $H_{0}$ at low redshifts while the sub-classes do not.

4.3 Forecasts with next generation surveys

Next generation surveys will probe the Universe at extended and complementary redshifts and scales delivering highly accurate data and offering an unprecedented insight into gravity on cosmological scales. Exploiting their potentiality and investigating the improvement in constraining cosmological parameters is one of the main goals of current investigations. Since the EFT framework does not rely on a specific model but allows to make general prediction on large classes of models, it provides a powerful benchmark for forecasting the cosmological signals to which the future missions will give access to. In the following, we review the cosmological forecasts analyses performed using the pure EFT approach [228, 146, 229, 70, 227, 230, 135]. Also in this case, the class of models that has been largely explored belongs to the pure Horndeski models.

The future Stage IV photometric redshift surveys LSST, the radio galaxy survey SKA and CMB-S4 experiments have been used to obtain forecasts for the Horndeski models in the $\alpha$-basis parametrized as follows [146]

[TABLE]

where $\Omega_{\rm DE}(t)$ is chosen to follow $w$CDM and $\alpha_{i}^{0}$ are constant. This modeling includes a phenomenological assumption for screening on small scales. The scale of the Vainshtein screening mechanism $k_{V}$, above which GR is recovered, has been set to be 0.1 h/Mpc according to numerical simulations [231]. Including screening effects helps, for instance, not to overestimate the surveys capacity to test gravity on such scales. The results of the analysis display an improvement in the constraints by a factor of 5 with respect to previous results based on present day surveys [205]. In particular, in Figure 19 (top panel) one can observe the effects of BAO and RSD data by an independent DESI-like experiment on top of CMB-S4, LSST and SKA measurements on the parameter constraints. Tighter constraints are obtained by setting $c_{t}=1$ [227] as shown in Figure 19 (bottom panel). Such constraints exclude variations of the effective Newtonian constant larger than 10% over the age of the Universe [227]. However, the contours obtained from the MCMC analysis are larger than those obtained from a Fisher analysis. This points to the danger of using Fisher forecasts for non-Gaussian likelihoods in forecast investigations. Finally, the model is investigated using two cosmic shear methods [230]: the tomographic method, where the correlations between the lensing signal in different redshifts bins allow to keep track of the redshift information, and the 3D approach, where all the redshift information is carried throughout the analysis. For an Euclide-like experiment, the 3D analysis is shown to constrain the model better than the tomographic approach by about 20% thanks to the increased redshift information (see Figure 20). The size of the constraints on $\Omega_{\rm m}$, the sum of neutrino mass and EFT functions is further improved by including non-linear corrections in the power spectrum as displayed in Figure 20 (bottom panel).

The synergy between future CMB-S4 and DESI has been explored on both the EFTCAMB and $\alpha$-basis for the Horndeski class of models [70]. Two parameterizations of the EFT functions are considered: constant and a time varying behavior with a smooth transition between early and late time values (see eq. (201)). When modeling the $\alpha$-functions as constants, CMB-S4 surveys are found to provide the tightest bounds on the effective mass Planck and tensor speed excess respectively by a factor 1.5 with respect to DESI and 2.5 with respect to Planck, while DESI and Planck show the same sensitivity. On the contrary, the Planck measurements are slightly stronger than CMB-S4 on the constant $\alpha_{B}$. When considering the time varying parameterization, CMB-S4 measurements are sensitive to both early and late time values while the former are, as expected, more efficiently constrained by the CMB surveys than LSS ones. In the EFTCAMB basis and for the constant parameterization, Ref. [70] finds that the sensitivity of CMB probes are unmatched when constraining the parameter $\Omega_{\mathrm{0}}^{\mathrm{EFT}}$. Additionally, CMB-S4 measurements better constrain $\gamma_{0}^{(2)\mathrm{EFT}}$ and $\gamma_{0}^{(3)\mathrm{EFT}}$. For the time varying parameterization, late time values of the $\gamma_{i}$ parameters do not change as compared to the constant case, while the forecast bounds on $\Omega_{\mathrm{0}}^{\mathrm{EFT}}$ slightly degrades. At early time, the constraints are mostly constrained by physical viability requirements.

Different time behaviors for the EFT functions $\{\Omega,\gamma_{1},\gamma_{2}\}$ of the EFTCAMB basis have been explored further by considering also the constraints on GWs which suggest $\gamma_{3}=0$ [135]. The remaining three EFT functions are modeled firstly with the scaling-a parameterization and then with the de-density form (respectively eq. (196) and eq. (197)) both on a CPL background. A constraint analysis is performed using present day observational data (Planck+BOSS DR12+$H_{0}$+JLA+KiDS) and forecasts from combinations of GC and WL for a prototype of next generation galaxy surveys with specifications like DESI and SKA2. These future surveys will be able to increase the precision on constraints of model parameters by one order of magnitude in both the parameterizations [135].

We conclude the section discussing the constraints on modifications of gravity when a direct coupling between DE and DM is present. In Ref. [228], the linear-de form (eq. (193)) of the $\alpha$-functions is considered on top of a $w$CDM expansion history. The dark matter coupling function is parametrized as follows:

[TABLE]

where $\beta_{\gamma}$ is a constant and $c_{s}(t)$ and $\alpha(t)$ correspond to the DE speed of propagation and kinetic coupling. The forecasts analysis is performed using the specifications of an Euclid-like survey for which three probes are considered: the galaxy power spectrum in redshift space, tomographic weak-lensing shear power spectrum and the correlation spectrum between the ISW effect and the galaxy distribution. This is done selecting three fiducial models: $\Lambda$CDM, a braiding and an interacting model. We display the results for the $\Lambda$CDM fiducial in Figure 21 where one can observe how the parameters are degenerate, yet some degeneracies can be effectively broken by combining all three observational probes. $1\sigma$ constraints on the MG parameters are of order $\sim 10^{-2}-10^{-3}$ (68%C.L.) for the first two fiducial models and one order better in the interacting fiducial model. The error on the dark matter coupling parameter is $\sim 10^{-4}$ in all cases. One can also appreciate from Figure 21 how an Euclid-like survey gives an order of magnitude tighter constraints than CMB-S4 experiments. Finally, the non-minimal coupling of DE enhances the effects of modification of gravity and reduces the statistical errors accordingly [79].

4.4 Gravitational waves

Alternatives to GR can lead to two specific modifications in the propagation of GWs: the first is the friction term affected by the running Planck mass, $i.e.$ $\alpha_{M}$, the second is the speed of propagation, $c_{t}^{2}=1+\alpha_{T}$. The former alters the amplitude of the GWs and the latter its phase. In detail, the evolution equation of the GWs on cosmological scales reads in Fourier space 141414A mass term can be included in the equation of GWs [232] but this does not hold for the classes of theories considered in this review.

[TABLE]

where $h_{ij}^{T}$ is the tensor component of the spatial metric.

$\alpha_{M}$ and $\alpha_{T}$ can lead to observable effects on CMB, on both the temperature and the polarization spectra at any time. In particular, information on primordial GWs is encoded in the B-modes at larger scales [233, 234, 235, 236]. The modification in the speed of propagation at early time can affect the position of the inflationary and of the reionization peaks in the B-modes [199, 237]. Furthermore, along with the expansion rate $H$, the running Planck mass can damp GWs: for the class of Horndeski models, if $\alpha_{M}<0$ the GW amplitude is smaller than that predicted by the standard scenario, the opposite holds if $\alpha_{M}>0$ [238]. A modified amplitude then introduces a difference in the GWs and electromagnetic luminosity distance, which can be tested by GWs experiments and standard sirens [239, 240, 241]. The next generation Laser Interferometer Space Antenna (LISA) [242] would constrain a constant $\alpha_{M}$ at redshift 1.5 with a precision varying between 0.03 and 0.13 independently of the underlying cosmological model [243]. A recent analysis combining LISA with CMB+BAO+SNIa data shows that $\alpha_{M}$ can be measured to an accuracy reaching 1.1% [244]. Joint measurements with standard sirens will not only directly constrain the running Planck mass but also give constrains on $H_{0}$ which might also help in resolving the $H_{0}$ tension [7, 8, 186, 187, 12]. We refer the reader to [232] for a review on GWs astronomy.

The tensor-to-scalar ratio $r$ is degenerate with $\alpha_{M}$ as they both affect the amplitude of the primordial peak [245]. Considering the Background Imaging of Cosmic Extragalactic Polarization 2 (BICEP2) data [246] in order to have $r$ close to zero, $\alpha_{M}$ must assume negative values, in particular it goes towards $-2$. However, negative values of the running Planck mass enhance the CMB BB spectra at large $\ell$, where the BICEP2 data would favor instead smaller values. This feature needs to be considered in light of other probes on scalar perturbations, indeed as extensively discussed in previous Sections $\alpha_{M}$ enters in the growth of structure, lensing and ISW effect. For example, considering the Brans Dicke theory, the background and scalar perturbations give $0\leq\alpha_{M}\leq 0.01$ and the degeneracy with $r$ is removed [245]. General cases do not allow for such tight bounds hence the B-modes can provide useful constraints on early-time MG. Planck and BICEP2 datasets constrain $c_{t}^{2}=1.30\pm 0.79$ and $c_{t}^{2}<2.85$ at 95%C.L. by assuming a power law primordial tensor power spectrum and $c_{t}^{2}<2.33$ at 95% C.L. if the running of the spectral index is allowed [237]. Forecasts for the next generation CMB satellites Cosmic Origins Explorer mission (COrE) [247] and Polarized Radiation Imaging and Spectroscopy Mission (PRISM) [248], will be able to constrain the speed of GWs at percent level [237].

As discussed previously, the modification of the speed of propagation of GWs generates a gravitational slip $\eta$ different from unity. The relationship between an anomalous speed of GW and $\eta$ is investigated in three cases [249]: Bi-metric, Einstein-Aether gravity and Horndeski theories. Horndeski theories are studied in their EFT formulation and are shown to a have enough freedom to hide dynamically the modification in the gravitational slip whenever $\alpha_{T}\neq 0$, but at the cost of making the perturbations evolve towards a divergent kinetic term. On the contrary, the other theories do not offer this possibility. This result can be used as theoretical argument for the interpretation of future observations if they favor $\eta=1$.

Finally, the anomalous speed of GWs can be used to break the degeneracy between MG and DE behaviors within Horndeski theories, $i.e.$ singling out self-accelerating models [250]. The linear shielded Horndeski model [251] retains this property but recovers a $\Lambda$CDM background expansion history and exhibits $\mu(t,k)=\eta(t,k)=1$ at the linear perturbations level in the QS regime. Within this scenario, the region allowing for self-acceleration is very narrow. In particular, for these models to have cosmic acceleration attributed to genuine MG, the present tensor speed has to be $<95\%$ than that of light and $\alpha_{M}$ $5\%$ less efficient than in $\Lambda$CDM [250]. Note that when imposing $c_{t}^{2}=1$, these models produce a $3\sigma$ poorer fit to cosmological observables as compared to $\Lambda$CDM [252], challenging the concept of self-acceleration within Horndeski theories. This model leads to the conclusion that surveys probing LSS and the background expansion alone are not sufficient to ultimately discriminate between $\Lambda$CDM and such modifications of gravity. Nevertheless it highlights how the addition of current and future GWs detections can break crucial degeneracies induced by MG and offer complementary means to test gravity.

4.5 Constraints on $\mu$, $\Sigma$ and $\eta$

Figure 22: Top panels: Figure 7 in Ref. [112]. 68% and 95% C.L CMB constraints on $\mu_{MG}(z=0)$, $\Sigma_{MG}(z=0)$ and $\gamma_{MG}(z=0)$ (in this review $\mu,\Sigma,\eta$ respectively). The EFT functions follow the linear-de form (eq. (193)) on $\Lambda$CDM background. Different combinations of stability conditions are imposed: stable (no ghost and no gradient) (green), plus the sub-luminal propagation for scalar mode ($c_{s}<1$) (red) and plus sub-luminal propagation of tensor mode ($c_{t}<1$) (blue). Bottom panels: Figure 5 in Ref. [174]. RSD constraints are compared with (red) and without (blue) the prior on the variation of the Newton constant for pure Horndeski models with $c_{t}=1$ parametrized using the de-1 form (eq. (194)) on a $\Lambda$CDM background.

A systematic investigation of possible DE/MG extensions to GR can be performed using the $\mu,\Sigma,\eta$ phenomenological approach [155, 156, 157, 158] reviewed in Section 3.1. Signatures of deviations from the cosmological standard model can then be captured by constraining these functions directly. Despite this approach being model-independent, the way in which constraints on these functions are obtained might instead be model-dependent: one has to choose somehow the underlying gravity model in order to obtain cosmological constrains on either $\mu,\Sigma$ or $\eta$. One can indeed make use of the EFT formulation to compute the constraints. In this case, given a parameterization of the EFT functions one can then translate the constraints on the EFT parameters into the phenomenological functions using the relations (52)-(53) in Section 3.1. Alternatively, one can parameterize directly the phenomenological functions. In the latter case, the connection with a specific model is lost, however depending on the chosen parameterization, $\mu,\Sigma,\eta$ can mimic the predictions of certain class of theories by appropriately choosing their time and scale dependencies. A public EB code which allows to do so is Modified Growth with CAMB (MGCAMB)151515 MGCAMB webpage: www.aliojjati.github.io/MGCAMB/ [111, 253]. Let us note that in both cases the constraints are highly sensitive to the particular parameterizations and priors assumed. A third approach instead considers $\mu,\Sigma,\eta$ as free functions (e.g. using a principal component analysis [254]) and as such the resulting constraints are completely model-independent. In this Section we focus on the first two approaches.

Firstly let us consider the approach where the constraints are derived within the EFT formulation. In principle, starting from each of the results discussed in Sections 4.2-4.3 one can deduce the corresponding constraints on $\mu,\Sigma$ or $\eta$. Major advantages in using the EFT approach to constrain these phenomenological functions come from the possibility to: impose the appropriate stability conditions to guarantee the viability of the chosen model (see Section 3.3); use observational and experimental priors such as those on the effective Planck mass, speed of GWs and Hubble parameter. Their impact on the $\mu-\Sigma$ or $\mu-\eta$ planes is largely demonstrated in Section 3 and when such requirements are used as priors in MCMC analysis they show a strong constraining power [73, 48, 112, 174]. We show an example in Figure 4.5 (top panel) [112] where the $c^{2}_{s}\leq 1$ and $c^{2}_{t}\leq 1$ priors push notably the constraints for pure Horndeski models in the $\eta(z=0)<1$ and $\Sigma(z=0)<1$ quadrant. This is in agreement with the sign conjecture on $\mu$ and $\Sigma$ [160] discussed in Section 3.2. The conjecture is also confirmed for pure Horndeski models with $c^{2}_{t}=1$ using the Solar-System prior on the variation of the present day value of the effective Planck mass combined with RSD data and a high redshift CMB prior [174]. The Solar-System prior allows notably to obtain a stringent bound on the present value of the gravitational slip parameter $\eta=1.0000\pm{9.3\times 10^{-4}}$ (95% C.L.). Having the effective Planck mass stringently bounded at low redshifts by the Solar-System prior and at high redshifts by the CMB one, the redshift evolution of $\eta$ is constrained by RSD data to be close to unity across matter domination. This implies $\mu$ and $\Sigma$ to be almost equal across redshifts as shown in Figure 4.5 (bottom panels). The fifth force contribution $\mu_{\mathrm{ff}}$ is favored at more than $2\sigma$ at present time from the bound $\mu(z=0)=\mu_{\mathrm{ff}}(z=0)=1.321^{+0.370}_{-0.284}$ (95% C.L.).

In terms of next generation galaxy surveys with specifications like DESI and SKA2, forecasts on $\mu$ and $\Sigma$ obtained from pure Horndeski models with $c_{t}^{2}=1$ parametrized as in eqs. (196)-(197) on a CPL background are found to reach the 1% level [135]. In Figure 23 [135], we show the forecasted 2$\sigma$ errors on $\Sigma(z)$. For the scaling-a models (M1), the errors are of $\mathcal{O}(10^{-3})$, and decrease towards $\mathcal{O}(10^{-4})$ for $z>1.5$. This is because $\Sigma(z)$ asymptotically tends to unity, independently of the cosmological parameters, which implies very small errors. For the de-density form models (M2), the errors are constant in $z$. One can notice that future data will allow to distinguish these EFT models from $\Lambda$CDM at more than 3$\sigma$, assuming the best-fit values obtained with present day data (Planck+BOSS DR12+$H_{0}$+JLA+KiDS) hold [135].

Let us now concentrate on the second approach, $i.e.$ a direct parameterization of the phenomenological functions. A parameterization commonly used in the literature is [255, 208, 256, 1]

[TABLE]

where $\Xi$ stands for either $\{\mu,\Sigma,\eta\}$. The latter has been considered by the Planck Collaboration [208, 1] and the best fit resulting from the combined analysis of several cosmological probes deviates from $\Lambda$CDM by more than $2\sigma$. The constraints obtained suggest that $\mu(z=0)<1$ and $\eta(z=0)>1$. This is explained by the preference of data for lower values of $\sigma_{8,0}$ relative to the standard model (see discussion in Ref. [112] and related references therein). Note that using a different data set combination this tension was not recovered [256]. Other examples of direct parameterizations and corresponding constraints can be found in Refs. [257, 258]. An important aspect to note is that each chosen direct parameterization is differently sensitive to the redshift evolution according to the time scaling of the model thus the tightness of the constraints obtained can indeed vary significantly [257].

The direct parameterization of the phenomenological functions has some advantages: it provides a null test of the standard cosmological model and enables to clearly identify the tendencies of data with a low amount of free parameters. However, being not based on a gravitational theory, stability and normalization requirements or additional observational priors cannot be considered or are hard to implement. These instead are easily included in the EFT approach. Comparing the results from Planck [208, 1] and the results obtained from the Monte-Carlo generation of viable pure Horndeski models discussed in Sections 3.2-3.3, it is clear that the marginalized posterior distributions obtained with the model in eq. (75) dwell in a region where Horndeski theories are mostly absent. Further imposing the normalization of the effective Planck mass to the Planck mass, $i.e.$ $\mu_{\mathrm{sc}}(z=0)=1$ prevents to recover models with $\mu(z=0)<1$ since $\mu_{\mathrm{ff}}(z)>1$ as shown in Figure 4.5. Additionally, in a given theory, $\mu,\Sigma,\eta$ are not independent as eq. (54) highlights and are hence bound to be correlated. Direct parameterizations of these phenomenological functions do not allow to retain the link they share affecting the constraints. Using instead the EFT approach allows to overcome this issue. However in both approaches the time evolutions of these phenomenological functions might result in simplified behaviors which poorly capture the ones obtained from full covariant theories, $e.g.$ scalar-tensor theories [212, 163, 213]. The risk here is to miss information in the data. A way to bypass this issue would be to extract model-independent information about the phenomenological functions from current and future data. An example to minimize the model dependence is to consider them as unknown functions of both redshift and scale, and bin them into a large number of narrow bins on a grid in the $(z,k)$ space [259, 254, 260, 176]. This approach would however miss most of the advantage described above and has its own limitations.

4.6 Specific modified gravity models

The EFT formalism has a twofold face when implemented in EB codes, indeed the background dynamics can be solved either with a pure approach as discussed in the previous Sections or by implementing a specific background solver for a chosen theory, $i.e.$ the mapping approach. In the latter case thanks to the mapping procedure discussed in Section 2.4, the EFT functions are fully specified. Let us note that the mapping approach can also involve a designer procedure to solve the background. In the following we review specific DE/MG models implemented in EB codes using the mapping approach.

Quintessence [261]: the action for quintessence reads

[TABLE]

where $V(\phi)$ is the potential of the scalar field $\phi$. It is possible to show that the freedom in choosing the functional form of the potential can be replaced by the choice of the equation of state, indeed in the EFT formalism one gets

[TABLE]

Then, one has only to specify the expansion history and use a designer approach to investigate minimally coupled quintessence models. The model is implemented in EFTCAMB. Using the CPL background, the stability conditions alone drastically limit the allowed range of variation of the CPL parameters $\{w_{0},w_{a}\}$ (removing the $w_{\rm DE}<-1$ space) and drive the observational constraints [113]. As a result, the region of the CPL plane that would be favored by WL data is eliminated, thus quintessence becomes significantly disfavored with respect to the standard cosmological model. An extension of the viability region of CPL could be performed with multifields quintessences models which might provide a better fit to data. Furthermore, the no-ghost, no gradient and mathematical161616Mathematical conditions are a set of requirements implemented in EFTCAMB worked out at the level of the dynamical equation for the perturbation of the scalar field $\pi$. They guarantee the stability of the theory in absence of the full set of physical conditions. Such condition can be replaced by the no-tachyonic condition which instead is a well physical motivated condition [179]. help reducing the tension between KiDS and Planck data sets [113].

• 2.

Designer $f(R)$-gravity [178, 262, 154, 98]: the action for $f(R)$-gravity is shown in eq. (37) with the corresponding mapping into the EFT framework. The model is implemented in EFTCAMB using the designer approach [73]. On a $w$CDM background, the free parameters are $B_{0},w_{0}$, $i.e.$ the present day values respectively of the scalaron and the equation of state for the DE. The results show that using a combination of data from Planck temperature and lensing potential spectra, WMAP low-$\ell$ polarization spectra (WP), and BAO, the constraints on the DE equation of state are $w_{0}\in(-1,-0.9997)$ at 95% C.L.. In particular, the lower bound is the result of the requirement of stability conditions which induce a strong correlation between $B_{0}$ and $w_{0}$ (see Figure 6 top panel), indeed when $Log_{10}B_{0}<-4$ follows $w_{0}\rightarrow-1$. Furthermore, the combination of the viability priors and the Planck lensing data allowed to obtain a stringent constraint also on $B_{0}$: $Log_{10}B_{0}=-3.35^{+1.79}_{-1.77}\;(95\,\%\rm C.L.)$. This result was confirmed later in Ref. [210].

• 3.

$f(R)$* Hu-Sawicki model* [263]: this model is characterized by a specific form of the $f(R)$- function which reads

[TABLE]

where $m=\rho/3M_{\mathrm{pl}}^{2}$, $n>0$ and $c_{1,2}$ are constant. The above functional form is introduced to mimic $\Lambda$CDM in the high-energy regime and to give rise to an accelerated expansion which is not driven by a true cosmological constant at low-energy. The model is implemented in EFTCAMB using the mapping approach [139]. In the latest analysis [139], three cases are considered: $n=1$, $n=4$, $n$ free, with only one extra free parameter, $i.e.$ the scalaron Compton wavelength $log(-f_{R}^{0})$, where $f_{R}^{0}\equiv df/dR(z=0)$. The constraints for the baseline datasets, D1 (Planck15, JLA, BAO) are weak in the case $n=1$ ($log(-f_{R}^{0})<-2.7$ at 95% C.L. ), in the other cases the constraints on $log(-f_{R}^{0})$ are not statistically significant. The inclusion of the WiggleZ data drives the bound on this parameter away from $\Lambda$CDM, for $e.g.$ in the $n=1$ case $log(-f_{R}^{0})=-3.4^{+1.4}_{-1.2}$, at 95% C.L., while considering D1+CFHTLenS the value of $\log_{10}(-f_{R}^{0})$ is driven back to its $\Lambda$CDM limit, $log(-f_{R}^{0})<-4.5$ at 95% C.L. (n=1). This analysis reveals for the first time a degeneracy between $\sigma_{8}$ and $f_{R}^{0}$. This was possible thanks to a full implementation of the background solver which is not forced to be $\Lambda$CDM. Interestingly, when $n=4$ and $log(-f_{R}^{0})>-2$, this degeneracy changes in direction, a feature also found in the posterior distribution of $\sigma_{8}$ and $H_{0}$. Furthermore, for $log(-f_{R}^{0})>-2$, this parameter shows a degeneracy with $H_{0}$ also, due to the fact that the effect of the background modification is not negligible within this bound.

• 4.

Hybrid-metric Palatini $f(\hat{\mathcal{R}})$ gravity [264]: the action is constructed by adding a Palatini correction $f(\hat{\mathcal{R}})$ to the usual Hilbert-Einstein term as

[TABLE]

where $\hat{\mathcal{R}}\equiv g^{\mu\nu}\hat{\mathcal{R}}_{\mu\nu}$ is the Palatini curvature and $\hat{\mathcal{R}}_{\mu\nu}$ is defined in terms of an independent connection $\hat{\Gamma}^{\alpha}_{\mu\nu}$. The mapping relations using the $\alpha$-basis read [265]

[TABLE]

where prime is the derivative with respect to $\ln\,a$. The background evolution is solved using the designer approach with $\Lambda$CDM background [265]. In this model, early modifications of gravity become significant after recombination while they decay towards the present. Background probes and Planck measurements do not show evidence for such effects. The constraints for the scalar field value are $|f_{\hat{\mathcal{R}}}(z=z_{\mathrm{on}})|\lesssim 10^{-2}$ where $z_{on}$ is the redshift at which the decaying early-time modification is introduced ($z_{\mathrm{on}}\sim(500-1000)$), and $|f_{\hat{\mathcal{R}}}(z=0)|\lesssim 10^{-8}$ (95% C.L.) at present time.

• 5.

Jordan Brans-Dicke theory (JBD) [95]: the JBD action reads

[TABLE]

where $\omega_{BD}$ is the constant Brans-Dicke parameter and $V(\phi)$ is a potential. GR is recovered for $\omega_{BD}\rightarrow\infty$. The theory is very well constrained at all scales: $\omega_{BD}>4\times 10^{4}$ from Shapiro time delay [266], $\omega_{BD}>6\times 10^{2}$ from Planck [267]. The model is investigated further in Ref. [146] with $V=0$ and a tracker solution for the background is found where $\phi=\phi_{0}a^{1/(\omega_{BD}+1)}$ and $\phi_{0}=(2\omega_{BD}+4)/(2\omega_{BD}+3)$ [268]. In the $\alpha$-basis the JBD model can be written as

[TABLE]

where the tracker solution implies $\alpha_{M}=1/(\omega_{BD}+1)$. This background evolution and the mapping relations are implemented in hi_class and the analysis reveals that the combined future Stage-IV surveys will be able to place the bound $\omega_{BD}>1.7\times 10^{4}$ [146], which is comparable to Solar-System and astrophysical tests.

• 6.

Covariant Galileon (CG) [103]: the Galileon field $\phi$ was first introduced in flat space with the Galileon symmetry $\partial_{\mu}\phi\rightarrow\partial_{\mu}\phi+b_{\mu}$ which guarantees that the resulting theory posses second order equations of motion [102]. Its generalization to a dynamical space-time leads to the breakdown of such symmetry when demanding the field preserves second order equation of motion [103]. The resulting theory is called Covariant Galileon (CG). Even though it does not satisfy the original Galileon symmetry it preserves the shift symmetry: $\phi\rightarrow\phi+c$. The CG action reads

[TABLE]

where $G_{\mu\nu}$ is the Einstein tensor, $X=\phi^{;\mu}\phi_{;\mu}$ and $\{;\}$ stands for the covariant derivative. Moreover, $c_{i}$ are constant dimensionless parameters and $M^{3}=M_{\mathrm{pl}}H_{0}^{2}$. One can fix a canonical normalization by choosing $c_{2}=-1$ [269]. Note the original version of the action also includes a $c_{1}\phi$ term, but if one wants the late time acceleration to be driven by the field kinetic energy, one should set $c_{1}=0$. The above model has been widely studied in literature because of its rich phenomenology which makes it a possible candidate to explain late time acceleration [165, 269, 147]. It was recently implemented both in EFTCAMB and hi_class thanks to the mapping procedure. Thus the EFT functions are fully determined once the background is solved. Three cases are explored [147, 140]: Cubic model ($G_{3}$), $i.e.$ $c_{3}\neq 0\,,\{c_{4},c_{5}\}=0$; Quartic model ($G_{4}$), $\{c_{3}\,,c_{4}\}\neq 0\,,c_{5}=0$; Quintic model ($G_{5}$), $\{c_{3},c_{4},c_{5}\}\neq 0$. Using data from the CMB (including lensing), BAO and ISW, $G_{3}$ is excluded as a viable candidate because it shows 7.8$\sigma$ tension with the data [147]. The $G_{4}$ and $G_{5}$ cases show a reduced viable parameter space due to the ISW data and in this region the goodness-of-fit is comparable to $\Lambda$CDM [147]. When using the combination of CMB (temperature and polarization), BAO, $H_{0}$ and WL measurements, the $G_{3}$, $G_{4}$ and $G_{5}$ models are also statistically ruled out [140]. This result is obtained using only cosmological data as such it is independent from any assumption coming from GWs [168, 169, 170]. We further discuss the implication of GW on DE/MG in Section 5.2.

• 7.

*Low-energy Hořava gravity [49, 270]: * Hořava gravity is considered a candidate for the ultra-violet completion of GR [49, 270]. The idea behind this theory is to add only higher order spatial derivatives in the action in order to modify the graviton propagator. To this purpose the theory is formulated using the 3+1 decomposition of the ADM formalism. Considering the above arguments, the action of Hořava gravity at low-energy can be written as follows [89]

[TABLE]

where $\left\{\lambda,\xi,\eta\right\}$ are dimensionless running coupling constants and $\bar{\Lambda}$ is the “bare” cosmological constant. The model parameters are constrained using both cosmological [271, 272, 273, 274, 275, 48] and astrophysical [276, 277, 278, 279] data. Furthermore, the theory is still viable [280] after the detection of GW170817 event [166, 167]. The theory is mapped into the EFT formalism and implemented in EFTCAMB [48]. Using the combination of CMB (temperature and lensing), galaxy power spectrum, local Hubble measurements, SNIa and BAO data, two cases are analyzed [48]: H3, where the three extra free parameters $\{\xi,\eta,\lambda\}$ are free to vary and H2 where the Parameterize post-Newtonian (PPN) bounds are imposed [281], which imply $\xi=\eta/2+1$. Thus only two parameters of the theory are left to vary $\{\eta,\lambda\}$ 171717Let us note that after the detection of GW170817 [166, 167] the viability region identified by the PPN has been revised and the condition $\xi=\eta/2+1$ becomes less relevant [280].. In the H3 case the constraints on the parameters are $\xi-1=-0.01^{+0.01}_{-0.02}$ and $log_{10}(\lambda-1)<-4.31$ (at 99.7%C.L) while $\eta$ is unconstrained. Furthermore, the constraint on derived parameters such as $G_{cosmo}/G_{N}$, where $G_{cosmo}$ is the cosmological gravitational constant appearing in the Friedmann equation, improves by one order of magnitude with respect to previous results [271]. Indeed, $G_{cosmo}/G_{N}-1<0.028$ at 99.7%C.L.. For the H2 case, the results are $log_{10}(\lambda-1)<-4.39$, $log_{10}(\eta)<-4.51$ and $G_{cosmo}/G_{N}-1<6.1\times 10^{-5}$ (at 99.7%C.L.). The latter is several orders of magnitude better than the Big Bang Nucleosynthesis bound [271].

• 8.

K-mouflage [282, 283]: these theories are characterized by a screening mechanism which acts through the derivative of the scalar field $\phi$. They are constructed by including a universal coupling of the scalar field to matter fields in the K-essence action. The action can be written as follows:

[TABLE]

where $\mathcal{M}^{4}$ is the energy scale of the scalar field, $g_{\mu\nu}$ obeys to the transformation $g_{\mu\nu}=A^{2}(\phi)\tilde{g}_{\mu\nu}$ between the Jordan frame and Einstein frame metric $\tilde{g}_{\mu\nu}$. The kinetic term in the Einstein frame is defined as $\tilde{\chi}=-\tilde{g}^{\mu\nu}\partial_{\mu}\phi\partial_{\nu}\phi/2\mathcal{M}^{4}$, and $K$ is a general function of the scalar field. The quantities with tilde are defined in the Einstein frame. The deviation with respect to $\Lambda$CDM is modeled through two functions:

[TABLE]

The model is implemented in EFTCAMB and two cases are investigated [141]: the first in which the background is solved in the mapping approach exploiting the full dynamics of the K-mouflage model and a second one where the scalar field is forced to reproduce a background evolution degenerate with $\Lambda$CDM (K-mimic). Using CMB, CMB lensing, SNIa and different galaxy catalogues and forecasts for future CMB probes, the $H_{0}$ bound alleviates the tension between Planck and low-redshift probes. Moreover, because the model predicts a suppression in the growth of matter perturbations relative $\Lambda$CDM, this results in a lower value for $\sigma_{8}$, easing the tension between Planck and WL measurements [14, 13, 15, 16, 186, 187]. Finally, the constraints on the parameters of the model are at 95% C.L.: $-0.04<\epsilon_{2,0}<0$ for K-mouflage and $0<\epsilon_{2,0}<0.002$ for K-mimic. They could improve by approximately one order of magnitude with future surveys such as COrE [284].

• 9.

Galileon ghost condensate (GGC) [285]: the model extends the cubic covariant Galileon [103] by including an additional higher-order field derivative $X^{2}$ following the Ghost-Condensate model [286]. The action reads

[TABLE]

where $a_{1,2,3}$ are constants. The model is implemented in EFTCAMB [142] and for convenience the following dimensionless variables are considered on a FLRW background

[TABLE]

A phenomenological analysis shows that interestingly the low-$\ell$ ISW tail is lower than in the standard $\Lambda$CDM scenario for $x_{3}\ll x_{2}$ [142]. Another relevant feature is that the equation of state of DE can be in the region $-2<w_{\rm DE}<-1$ at low redshifts [285]. The constraints on the present day values of the above parameters using a combination of CMB, BAO, RSD and SNIa data [142] yield: $x_{1}^{(0)}$ and $x_{2}^{(0)}$ of order 1, with $x_{1}^{(0)}$ negative and $x_{2}^{(0)}$ positive, $x_{3}^{(0)}$ smaller than $x_{2}^{(0)}$. Finally, according to the Deviance Information Criterion (DIC) [287] and the Bayesian evidence factor ($\log_{10}B$) [288, 289], the GGC model is found to be statistically preferred over $\Lambda$CDM even with additional model parameters. In details, the above quantities computed with respect to $\Lambda$CDM give: $\Delta{\rm DIC}=-0.6$ (negative value supports the GGC model) and $\Delta\log_{10}B=5.1$ (value $>2$ favors the GGC model) for the complete dataset [142].

• 10.

Beyond Horndeski dark energy model [285]: The model (hereafter BH) belongs to the quartic-order GLPV theories and extends the GGC model in action (87) by including a beyond Horndeski term. The action reads

[TABLE]

where $a_{1,2,3,4}$ are constants. Along with the dimensionless parameters in eq. (88) on a FLRW background, one can consider an additional parameter which reads

[TABLE]

The above parameter defines the deviation with respect to Horndeski model, the so-called beyond Horndeski parameter $\alpha_{H}$, which is given by

[TABLE]

The model is implemented in EFTCAMB and investigated using CMB, BAO, RSD and SNIa data [143]. A tight upper bound $|\alpha_{H}^{(0)}|\leq{\cal O}(10^{-6})$ on the present day value of $\alpha_{H}$ was found [143]. This is mostly due to the shift of CMB high-$\ell$ peaks induced by the early-time modification in the cosmological background and linear perturbations arising from the dominance of $\alpha_{H}$ in the DE density. Another bound on $\alpha_{H}$ comes from the GW decay to DE and is of order $10^{-10}$ [171]. However let us notice that it is still a matter of debate whether GWs bounds from LIGO can be applied to MG at large linear scales [290] (see Section 5.2 for further discussion). The complete analysis of the BH model also reveals it to suppress the low-$\ell$ ISW tail of the CMB TT power spectrum with respect to $\Lambda$CDM due to the existence of the cubic term ($x_{3}$) and Ghost condensate one ($x_{2}$). The modified background expansion also affects the high multipoles CMB power spectrum. These features allow for the BH model to fit the data better than $\Lambda$CDM according to the $\chi^{2}$ statistics, however the DIC criterion slightly favors the latter [143]. From this result follows that since the BH model with $\alpha_{H}=0$ corresponds to the GGC model which is favored over $\Lambda$CDM [142] as discussed before, there are no particular signatures for deviations from Horndeski theories in current data.

4.7 Neutrinos and modified gravity

Massive neutrinos leave an imprint on cosmological observables: they impact the matter power spectrum and the shape of CMB anisotropies significantly while neutrino masses in the sub-eV to eV range alter the expansion history at the epoch of radiation-matter equality [291, 292]. Massive neutrinos change the height of the first acoustic peak of the CMB TT power spectrum due to the early ISW effect [293]: during the transition from radiation to matter epoch, the evolution for the metric perturbation changes hence the photon geodesics as well. This happens near the epoch of photon decoupling leaving a signature in the CMB anisotropies. A further contribution of massive neutrinos to the early ISW effect comes from their transition from relativistic to non-relativistic regime if their mass is of order 1 eV. Furthermore, neutrinos have large thermal velocity so they do not fall into the potential wells at $k$ larger than the neutrino free-streaming comoving wavenumber. This weakens the gravitational potential wells and suppresses the growth of structure on small scales. Some of these effects can therefore be degenerate with those of DE and MG and generally the constraints on neutrinos mass depend on the cosmological model assumed [294, 295, 296, 297, 298, 253, 299]. The degeneracy between alternatives to GR and massive neutrinos has been revisited in the EFT framework [300, 301].

In the context of designer $f(R)$-gravity as implemented in EFTCAMB using a combinations of Planck 2013, BAO measurements and LSS data from WiggleZ, a lower degeneracy with respect to previous literature results is found [300]. This is motivated by the fact that the dynamics of the $f(R)$ model is fully exploited with EFTCAMB (the code does not use QS approximation). The bounds on the Compton wavelength parameter and the neutrino mass read: $\log_{10}B_{0}<-4.1$ for a fixed $\Sigma m_{\nu}=0.06$ eV, which set a new upper limit on $B_{0}$; $\log_{10}B_{0}<-3.8$ for varying neutrino mass and $\Sigma m_{\nu}<0.32$ at 95% C.L.. The improved bounds on these parameters are driven mostly by the WiggleZ data which are highly sensitive to changes in $B_{0}$ and thus are able to partially break the degeneracy between these two parameters. The second model considered is a non-minimally gravitational coupling model parametrized through the linear EFT model (67) with a $\Lambda$CDM background. In this model, no sizable degeneracy is found and the bound obtained on the free parameter is $\Omega_{0}^{\mathrm{EFT}}<0.05$ (95% C.L.). The constraint is slightly improved with respect to what was previously obtained [73] where no massive neutrinos were considered. Finally, the constraint on the sum of neutrino mass is $\Sigma m_{\nu}<0.26$ at 95% C.L.. In this case the CMB lensing drives the constraints on the coupling constant $\Omega_{0}^{\mathrm{EFT}}$ and the bound on $\Sigma m_{\nu}$ is slightly looser. The combined dataset without the CMB lensing results in an slightly improved value for $\Sigma m_{\nu}(<0.25)$.

The degeneracy has been explored in the $\alpha$-basis with the $\alpha$-functions parametrized as the z-transition form in eq. (200) on a $\Lambda$CDM background [301]. The neutrino mass is found to be partially degenerate with the free parameters $c_{i}$. In particular, using forecast CMB and galaxy power spectrum datasets, one observes the parameter $c_{B}$, characterizing $\alpha_{B}$, to dominate the correlation with the total neutrino mass. Furthermore, $\alpha_{B}$ can cancel the power suppression due to the massive neutrinos at a given redshift. The breakdown of such degeneracy depends on the cosmological LSS data used at different redshifts. Next generation surveys such as Euclid would limit but not fully break the degeneracy between these two parameters [146], where no apparent degeneracy between MG and the sum of neutrino masses is obtained given the forecasted precision on the $\alpha_{i}$ parameters.

5 Astrophysical implications

In this Section, we review the implications astrophysical constraints have on the parameter space identified by the EFT functions. The EFT framework discussed in this review holds in the linear regime. In order to connect the EFT functions with physical quantities describing astrophysical processes, one can assume the validity of a linear treatment only if modifications of gravity are not screened or weakly screened, i.e. in a regime where non-linearities are subdominant. In this context, we review the bounds which constrain the strength of gravity inside massive astrophysical bodies such as dwarf, neutron stars, pulsars and galaxy clusters and the constraint derived from the detection of the GWs event GW170817 and its electromagnetic counterpart GRB170817A.

5.1 Massive astrophysical bodies

A common feature of MG models is the property to screen the fifth force on small scales or high density environments, notably where Solar-System and astrophysical tests constrain gravity to be that of GR with astonishing precision [302, 303]. The GR limit is then recovered in such models thanks to screening mechanisms [5]. In scalar-tensor theories, the latter can be classified using the type of interaction as a phenomenological criterion. It can depend on the local field value as in the cases of symmetron [304], chameleon [161, 162] and dilaton [305] mechanisms; or on the first derivative $\partial\phi$, $e.g.$ K-mouflage screening mechanism [306] and finally if it acts through the second derivative $\partial^{2}\phi$, typical of the Vainshtein mechanism [307, 102, 308, 309].

The Vainshtein screening mechanism is characteristic of Horndeski and GLPV theories. A peculiarity of GLPV theories is that the screening of the extra DoF is not complete however. Outside an extended object the Vainshtein mechanism allows to fully reproduce GR equations, yet it exhibits a “partial breaking” inside astrophysical bodies [310]. The perturbed potentials of the Minkowski metric, $\phi(r)$, $\psi(r)$, inside the objects obey the equations [310, 311]

[TABLE]

where $\Upsilon_{i}$ are dimensionless constants depending on the specific theory and $\tilde{M}(r)$ is the mass inside the object. Inside a matter overdensity the gravitational interaction becomes dependent on local matter density and the gradient of the gravitational potentials are no longer equal. This is transcribed by the fact that both the $\Upsilon_{i}$ are different functions of the MG couplings. This “breaking” of the Vainshtein screening effect opens a new window to constrain bounds on GLPV models through massive objects, such as dwarf, neutron, hyperon and quark stars and galaxy clusters [312, 313, 314, 311, 315]. In the EFT formulation, the $\Upsilon_{i}$ are related to the EFT coupling functions by [314, 313]

[TABLE]

Note that the breaking occurs if and only if $\alpha_{H}\neq 0$, hence it does not apply to Horndeski models. Thus, bounds on these two parameters can be seen as constraints on the EFT functions and eventually they can be used to rule out competitors to $\Lambda$CDM.

Observations of several low mass red dwarf stars give a conservative upper bound on the parameter $\Upsilon_{1}\lesssim 0.4$ at redshift zero [312]. Stable spherically static stellar solutions also require $\Upsilon_{1}>-2/3$ [313], which is satisfied when considering the lower bound $\Upsilon_{1}>-0.51$ from the consistency of the Chandrasekhar mass with the lowest mass white dwarf [318]. Constraints on $\Upsilon_{i}$ parameters can be also obtained from extragalactic measurements of galaxy cluster profiles [314]. Using X-ray and lensing profiles of galaxy clusters from XMM Cluster Survey [319] and CFHTLenS [320], a stringent constraint on $\Upsilon_{1}$ and the first estimation of $\Upsilon_{2}$ are obtained at high redshift ($0.1<z<1.2$), $i.e.$ $\Upsilon_{1}=-0.11^{+0.93}_{-0.67}$ and $\Upsilon_{2}=-0.22^{+1.22}_{-1.19}$ at $2\sigma$ [311]. These constraints directly translate into bounds on the EFT functions. One can further restrict the allowed parameter space for the EFT functions by considering the additional bound on the speed of propagation of GWs which leads to $c_{t}^{2}=1$ ($i.e.$ $\alpha_{T}=0$) [167]. Then, only $\alpha_{H}$ and $\alpha_{B}$ enter in the definitions of $\Upsilon_{i}$. The remaining viable regions in the $\alpha_{H}$-$\alpha_{B}$ plane are showed in white in Figure 24 [316]. Let us note that the second line in eq. (92) is further modified with the inclusion of an additional constant $\Upsilon_{3}$ if theories such as DHOST are considered [321, 322, 323].

The Vainshtein mechanism has another peculiar characteristic around massive objects, the so-called piercing effect [317]. Within the screened region, although the background value of the scalar field is suppressed its gradient is not bound to vanish for theories bearing the shift symmetry, as it is the case of a sub-class of GLPV theories. Theories with anomalous speed of GWs and coupling to matter for GWs ($i.e.$ an effective Planck mass) can therefore be tested at astrophysical scales. For example, observations of the Hulse-Taylor pulsar led to a model-independent constraint on the local value of $c_{t}$ to the level of $10^{-2}$ [317] through the bound

[TABLE]

where $G_{\mathrm{GW}}$ is the GWs coupling to matter. Combining this bound with the PPN constraint from the Cassini spacecraft experiment on the screened remnant of the gravitational slip parameter $\eta_{\mathrm{sc}}-1=(2.3\pm 2.1)\times 10^{-5}$ [266] induces the constraints on the EFT functions, $\alpha_{T}=c_{t}^{2}-1$ and $\alpha_{H}$, shown in Figure 24 (right panel).

The bounds reviewed in this Section can be used to complement those obtained from cosmological scales in order to further improve our knowledge about the gravitational interaction.

5.2 GW170817 and GRB170817A

The first observed merger of a binary neutron star system on August 17, 2017 occurred through two channels: first the LIGO and Virgo collaborations detected the GWs signal from this event, known as GW170817 [166] and after $1.74\pm 0.05$ s the Fermi and INTEGRAL gamma-ray telescopes observed the gamma-ray burst GRB170817A [167]. The time delay between the two detections constrained the difference between the speed of GW and the speed of light to be $-3\times 10^{-15}\leq c_{t}-c\leq 7\times 10^{-16}$ [167].

This result has a severe impact in selecting viable MG theories compatible with such tiny bound [168, 169, 170]. Applying this constraint on the speed of propagation of GWs in the EFT formalism implies $|\alpha_{T}|<10^{-15}$. If one imposes exactly such condition at any time, broadly speaking it leads to the conclusion $c_{t}^{2}=1$. Now considering the EFT action (5), one obtains three separate conditions depending on the chosen sub-class of models:

[TABLE]

The cases of GLPV and Horndeski theories received particular attention [168, 169]. Using the above relations, the GLPV Lagrangian reduces to [168]

[TABLE]

from which one can notice that the quintic GLPV Lagrangian vanishes. Horndeski theories with $c_{t}^{2}=1$ can be obtained from the above Lagrangian considering $B_{4}$ to be solely a function of the scalar field, $i.e.$ $B_{4}(\phi)$. It is clear that also in this case the quintic Horndeski Lagrangian is completely ruled out. Note that the relations leading to action (97) do not represent a fine tuning in the theory because the choice $c_{t}^{2}=1$ is protected against large quantum corrections [324, 325]. According to action (97) some well known MG models were ruled out [170], $e.g.$ Quartic and Quintic galileon [102, 103], Fab four [326], de Sitter Horndeski [327], $G_{\mu\nu}\phi^{\mu}\phi^{\nu}$[328], $f(\phi)$-Gauss-Bonnet [329]. Considering the beyond Horndeski models, one can exclude quartic and quintic GLPV [37], quadratic DHOST (with $A_{1}\neq 0$) [38] and cubic DHOST [330]. Additionally, in a recent work the decay of GWs into DE fluctuations in presence of Lorentz breaking was investigated [171]. The dominant decay channel is the decay of GWs into two scalar fluctuations $\gamma\rightarrow\pi\pi$ and a second channel leads to $\gamma\rightarrow\gamma\pi$. In both cases the decays are driven by a coupling proportional to $\alpha_{H}$. Values of $\alpha_{H}\neq 0$ lead to a large decay rate of the GWs implying no wave would reach the detector. As consequence $\alpha_{H}$ is forced to be of order $10^{-10}$ and theories such as GLPV and some sub-classes of DHOST models are further ruled out. Many other models passed the bound on the GWs speed, such as Cubic Horndeski models discussed in Refs. [102, 103, 331, 332], the Kinetic Gravity Braiding model [99], the shift symmetric GLPV model [285], the tracking and scaling DHOST theories discussed in Ref. [333], Einstein Aether theory [334] and low-energy Hořava gravity [280, 279]. For a complete review about viable models after GW170817 see Ref. [335].

The constraints in eq. (5.2) assume $\alpha_{T}=0$ for any background. In other words, they have to be satisfied for any value of $\ddot{\phi}$. However, by requiring the scalar field to satisfy the equation of evolution for the scalar field “dynamically”, it is possible to obtain other constraint relations according to which a non trivial quintic Horndeski Lagrangian can be rescued. The latter however is ruled out by the effects of large scale inhomogeneities [336].

The range of application of the LIGO bound on MG models is still subject of debate. Let us mention that the applicability of the bound on $\alpha_{T}$ at any time is questionable since the source of GWs is at redshift $z\simeq 0.009$. Then such bound should be applied only to constrain the speed of MG models at recent time, $i.e.$ $z<10^{-2}$ [184, 335]. In addition, one of the issues put forward is that LIGO measurement relates to frequencies of $10-100$ Hz which correspond to energy scales several orders of magnitude larger than those describing DE. Typically the EFT description of low-energy phenomena breaks down at a cutoff $\sim 100$ Hz. In this regard the measurement of the speed of GWs can be considered dependent on the frequency at which it was measured $c_{t}(k_{LIGO})$ [290]. Then, the EFT may predict a sub-luminal propagation at low-energy since the speed of GWs is close to unity at LIGO scales, thus, $\bar{M}^{2}_{3}(k_{LIGO})=0$ yet $\bar{M}^{2}_{3}(k=0)\neq 0$ at cosmological scales where the transition mechanism would be provided by the partial UV completion of the theory [290]. An experimentum crucis for the GWs in cosmology will be the future LISA mission [242] with sensitivity near $10^{-3}$ Hz. Only at this frequency one can eventually aim at definitively ruling out all the theories with $c_{t}\neq 1$ and restrict the viable parameter space if the LIGO constraint will be confirmed.

The different scales at which MG phenomena are expected and LIGO measurement happened have been further discussed [337]. The former acts on scales $\sim H_{0}^{-1}$ while one can associate the lookback time of LIGO observation $10^{-4}H_{0}^{-1}$ to the latter. The corresponding wave number of GWs is $\sim 10^{19}H_{0}$. Defining $K_{grav}=k_{grav}/H_{0}$, one can use $K_{grav}$ to suppress modifications to gravity on small scales. The phase and group velocity of GWs for a wide range of models can be written as [337]

[TABLE]

where $M^{2}_{GW}(t)$ is the time dependent, dimensionless graviton mass and $K=k/aH$. In the limit $M_{GW}\ll K_{grav}\sim 10^{19}$ it follows that $v_{p}=v_{g}=\sqrt{1+\alpha_{T}}/a$ which implies $|\alpha_{T}|<10^{-15}$ as pointed out in previous works. However, the massive graviton could still have significant effects on cosmological scales.

In conclusion, investigation of gravitational models with $c_{t}^{2}\neq 1$ would lead to obtain independent constraints from cosmological observations in order to complement the astrophysical observations. Future space-based GW detectors such as LISA will enlarge the reach of multi-messenger GW astronomy. It will then be possible to test gravity to a much higher precision and eventually help in shedding light on the controversy about the applicability of LIGO bounds at cosmological scales.

6 Conclusion and outlook

The EFT provides a unifying framework for gauging general classes of DE/MG theories at large cosmological scales in terms of a variety of free functions of time. One can thereby make predictions and interpret observations directly in the space of class of theories and not within a single paradigm. For example, $\mu$ and $\Sigma$ are known to be a powerful phenomenological parameterization to accurately describe MG effects in the growth of structures and the lensing of light. Observational evidences of any deviations from GR in these functions are hard to connect with classes of MG models without their interpretation in terms of EFT functions yet possible for single specific models. EFT thus provides an appropriate framework to classify different aspects of MG according to their signatures [160]. It also helped in obtaining interesting trends in these phenomenological functions which are of large applicability, $e.g.$ the $\mu-\Sigma$ conjecture [160, 173, 136]. These informations can be translated into specific models since the EFT framework preserves the link with covariant theories allowing then to restrict the forms of the general functions characteristic of DE/MG Lagrangians. The use of numerical tools built employing the EFT framework, such as EFTCAMB [72, 73], hi_class [74], COOP [75] and EoS_class [76] made the explorations of DE/MG effects on observables straightforward. This initiated a systematic analysis of alternative models against cosmological data. A recollection of these can be found in Appendix C summarized in tables. Considerable progress was also made in identifying the stable parameter space of theories of gravity [85, 73, 163, 179]. General and theoretically rigorous conditions were derived in the EFT framework encompassing the most significative class of models [118]. Such conditions are now systematically enforced in numerical codes [85, 112, 113, 73, 48, 113].

Such a unifying description comes at a price: the functional form of the time dependent EFT functions is unknown. Observations generally do not have enough power to fix continuous functions of time but only numbers. One thus resorts to a phenomenological modeling of such functions. In other words, the unknown information contained in the structural functions is compressed into a finite set of parameters, but one faces the challenges of general parameterizations against oversimplifications. The risk would be to miss significant DE/MG signatures or eventually to give a false alert [212, 213]. The chosen parameterization should thus be universal enough to explore most of the space of stable theories and yet be effectively constrainable by observations. In other words, one must adopt functional behaviors useful for revealing the nature of cosmic acceleration and concurrently avoid including too many free parameters as they might loosen the constraining power of data [112]. Fortunately, cosmological observables seem not to be extremely sensitive to short time-scale variations and therefore smooth parameterizations are in general sufficient to describe the theory space in a satisfactory way [338]. A convenient way to fix the EFT functions might be to parameterize directly the stability conditions and derive, subsequently, the evolution of the EFT functions [184, 183, 185, 339]. Alternatively one can use more sophisticated data-driven analysis to reconstruct characteristic functions such as EFT functions, $\mu,\Sigma$ and $w_{\rm DE}$ and then derive specific model properties [107, 225]. Beyond this issue, the EFT framework, although not complete yet, already helped to acquire deeper knowledge about the nature of gravity force, and derive novel predictions at cosmological scales.

The attempt for a fundamental understanding of the nature of gravity will require new efforts in both theoretical and observational sides. The next generation of cosmological surveys such as Euclid [67], DESI [68], SKA [69] are specifically tailored to study the impact of DE on the distribution of clustered matter. They will deliver highly accurate data offering an unprecedented insight into gravity on cosmological scales. The EFT framework is likely to become a sound benchmark to interpret data in the context of MG and extract precious information about fundamental physics. Yet a number of improvements are required to make this approach as complete as possible. For instance, the original EFT description is not properly armed for a correct description of screening effects at intermediate non-linear scales. These are instead very important in the interpretation of data since substantial part of GC, CMB lensing and WL data originates from non-linear scales. A phenomenological way to mimic screening mechanisms has been recently introduced [146, 214] which highlighted the urgency of the inclusion for non-linearities in the data analysis as these have been proven to increase the constraints on the EFT functions [230]. In this regard, a number of pioneer work contributed to incorporate non-linear effects to the EFT framework [62, 64, 65]. Further developments will offer the possibility to include corrections to the power spectrum coming from non-linearities as well as high order correlation functions. In parallel, the interpretation of data in light of stability conditions opened in certain cases the possibility of multi-fields models to explain the region of the parameter space favored by data [113]. In this regard, an interesting field of investigation could be to generalize the EFT construction to include additional DoFs whose dynamics is relevant at late time following the example of the EFT of multi-field Inflation [340].

The EFT framework reveals to be an innovating and fascinating research field with very promising prospects opened for improvements and extensions. This not only stems from its flexibility and user friendly approach but also for the era within which it is developed. We believe that in coming years a major breakthrough in the understanding of our Universe will be possible. The window to test gravity is not limited to cosmological observations but now extends to multi-messenger probes which will allow to place constrains on MG/DE models to a much higher accuracy. This becomes even more exciting as our modeling of the Universe perfects, as new ways of making data whisper their secrets are developed from novel machine learning [341, 342] and gaussian process techniques, without forgetting that evidence for models beyond the standard model is starting to arise at the bayesian level [142].

Acknowledgments

We are grateful to E. Bellini, J. Beltrán Jiménez, C. Marinoni, F. Piazza and A. Silvestri for their detailed comments and precious feedback on the manuscript. We thank B. Hu, L. Lombriser, R. Maartens, M. Martinelli, F. Pace, G. Papadomanolakis, S. Peirone and D. Vernieri for useful discussions. We acknowledge the authors of [160, 136, 107, 175, 85, 73, 106, 163, 179, 173, 192, 193, 145, 201, 200, 204, 207, 205, 220, 112, 146, 227, 230, 228, 174, 135, 316, 317] for the permission to use their figures and we would like to thank also G. Brando and E. Linder for providing us an adapted version of their figure. The research of NF is supported by Fundação para a Ciência e a Tecnologia (FCT) through national funds (UID/FIS/04434/2019), by FEDER through COMPETE2020 (POCI-01-0145-FEDER-007672) and by FCT project “DarkRipple – Spacetime ripples in the dark gravitational Universe” with ref. number PTDC/FIS-OUT/29048/2017. LP is supported by the South African Radio Astronomy Observatory (SARAO) and the National Research Foundation (Grant No. 75415).

Appendix A Acronyms and symbols

This Appendix is dedicated to the acronyms and symbols used throughout the review. For the sake of clarity we have collected them respectively in Table 2 and Table 3.

Appendix B Alternative basis and pure EFT parameterizations

In this Appendix, we collect the different basis used in literature to identify the EFT functions in action (5). Table 4 provides the relations among them. The ones that have been largely used in this review are: the EFT basis in action (5), the EFTCAMB basis [86] and the $\alpha$-basis [128]. The EFTCAMB basis also includes the $\gamma_{6}$ function defined as [86]

[TABLE]

according to the action (5). It corresponds to eq. (2.6) in the $\alpha$-basis. This EFT function is necessary to parameterize Lorentz violating effects.

[FIGURE:]

We summarize the most common pure EFT parameterizations used to fix the functional form of the EFT functions. In the following, we denote an EFT function regardless of the basis considered by $\mathfrak{F}_{i}$, where $i$ spans on all the EFT functions in a given basis, $e.g.$ in the $\alpha$-basis $\mathfrak{F}_{i}\equiv\{M,\alpha_{M},\alpha_{B},\alpha_{K},\alpha_{T},\alpha_{H}\}$. Let us note that in the following we use $c_{i}$ to identify general constant parameters. This notation is not a common one, thus the name of the coefficients needs to be adapted to each paper considered in this review.

The pure EFT parameterizations are:

linear-de form: the EFT functions are chosen to be proportional to $\Omega_{\rm DE}(z)$ as follows

[TABLE]

where $c_{i}$ is the constant parameter associated to the coupling $\mathfrak{F}_{i}$. We note that the normalization with $\Omega_{\rm DE}(z=0)$ is not always used.

• 2.

de-N form: this parameterization assumes an expansion of the EFT functions in terms of $\Omega_{\rm DE}=1-\Omega_{\rm m}$. In particular it has been built in the formalism of [85] for Horndeski models. The EFT functions in this basis are defined as follows:

[TABLE]

where the free parameters $c_{i,0},c_{i,1},c_{i,2}$ are associated to the EFT functions $\mathfrak{F}_{i}$. The index $n$ is: $n=0$ for $\epsilon_{4}$, $n=1$ for $\mu$ and $\mu_{3}$, and $n=2$ for $\mu_{2}^{2}$. The correspondence with the basis used in this review is in Table 4. Note that the normalization with $1-\Omega_{{\rm m},0}$ is not always used. In the review we have defined de-1 form as the expansion up to first order ($c_{i,0}$, $c_{i,1}$ and $c_{2,i}=0$) and de-2 form the one including the second order ($c_{i,0}$, $c_{i,1}$ and $c_{2,i}$).

• 3.

6-parameters form: the EFT functions resemble eq. (194) with $\{c_{i,1},c_{i,2}\}=0$. For the running Planck mass function only, the coefficient $c_{M,0}$ is not a constant yet a function of $\Omega_{\rm m}$. Considering specifically a $w$CDM background it assumes the form [85]

[TABLE]

where $\alpha,\beta$ are constants. Thus one has two parameters plus three other parameters from the remaining EFT functions and finally the 6th parameter is $w_{\rm DE}=w_{0}$.

• 4.

scaling-a/linear scaling-a form: this parameterization considers the behavior of each EFT function to scale with the scale factor $a$ as follows

[TABLE]

where there are two free parameters $c_{i}$ (the amplitude) and $q_{i}$ (the slope) per coupling $\mathfrak{F}_{i}$. The linear scaling-a is obtained when the slope is zero ($q_{i}=0$).

• 5.

de-density form: the EFT functions are parametrized to be proportional to the DE density ($\rho_{\rm DE}$), i.e. $\mathfrak{F}_{i}\propto\rho_{\rm DE}$. In particular using a CPL parameterization for the background it implies:

[TABLE]

where $c_{i}$ are constants and $w_{0},w_{a}$ are the CPL parameters.

• 6.

e-fold form or $1+tanh$ form : this parameterization allows a transition for the EFT behavior from unity in the past to a constant value in the future as

[TABLE]

where there are 3 free parameters: $c_{i}$ giving the amplitude of the transition, $a_{t}$ the scale factor of the transition and $\tau$ its rapidity.

• 7.

hill form: this functional form follows the e-fold form previously presented but allows for the coupling to have a negative amplitude. It reads

[TABLE]

where the definition of the parameters follows the e-fold form.

• 8.

z-transition form: this behavior allows to switch on modifications to GR at a given redshift $z_{th}$ (a transition redshift) with a transition given by $\Delta z$. The expression is given by

[TABLE]

For example, the transition redshift $z_{th}$ has been chosen close to the redshift of neutrinos becoming non-relativistic and $\Delta z$ comparable to $z_{th}$ in Ref. [301].

• 9.

early/late time transition form: this behavior allows the EFT functions to shift from different early and late time values with a smooth transition. The expression is given by

[TABLE]

where $\mathfrak{F}_{i,early}$ and $\mathfrak{F}_{i,late}$ are early and late time values of $\mathfrak{F}_{i}$, $a_{T}$ is the scale factor at the time of the transition and $\Delta a$ is the transition sharpness.

• 10.

Padé expansion: the Padé expansion of order [N/M] is the rational function

[TABLE]

where the truncation order is given by $N$ and $M$, $\{c_{n},b_{m}\}$ are constants coefficients and $a_{0}$ is the point around which the EFT function is expanded. Choosing $a_{0}=0$ would select models exhibiting thawing behaviors and $a_{0}=1$ those having freezing behaviors. This modeling describes well EFT functions which show a transition from one value at small $a$ to another at large $a$. This parameterization is used in some of the Monte-Carlo exploration discussed in the review.

• 11.

Taylor expansion: it is a Taylor polynomial expansion defined as follows

[TABLE]

where $N$ is the order at which the expansion is truncated, $a_{0}$ is the point around which the EFT function is expanded (selecting $a_{0}=0$ or $a_{0}=1$ would give respectively thawing and freezing behaviors) and $c_{n}$ is a set of constant coefficients. This parameterization is usually used in the Monte-Carlo exploration and given $c_{n}$ coefficients with the same prior distributions would favor the lower order terms.

• 12.

Polynomial expansion: this choice for the parameterization of the EFT functions follows the Taylor expansion with the difference that the term $n!$ in the denominator is not present. Its form is given as follows

[TABLE]

where $N$ is the order at which the expansion is truncated, $a_{0}$ is the point around which the EFT function is expanded (selecting $a_{0}=0$ or $a_{0}=1$ would give respectively thawing and freezing behaviors) and $c_{n}$ are constant parameters. Differently from the Taylor expansion, in this case the high order terms are not suppressed because the $n!$ is absent. The polynomial expansion is also used in some of the Monte-Carlo exploration.

Appendix C Constraints summary

We summarize the observational constraints and forecasts on the pure EFT parameterizations discussed in Sections 4.1, 4.2 and 4.3 in Tables 5, 6, 7. They are organized according to the chosen parameterization making the comparison among different data sets easier.