Cosmological constraints and phenomenology of a beyond-Horndeski model

TL;DR

This paper investigates a beyond-Horndeski dark energy model, deriving tight observational constraints on its parameters and comparing its fit to cosmological data against the standard Lambda-CDM model.

Contribution

It provides the first comprehensive observational constraints on the beyond-Horndeski parameter lpha_H and compares the model's fit to current cosmological data with lpha_H=0 and lpha_H.

Findings

Planck CMB data constrains lpha_H to 0^{-6}

The model suppresses the large-scale ISW effect compared to DM

The model fits the data slightly better than DM according to  statistics

Abstract

We study observational constraints on a specific dark energy model in the framework of Gleyzes-Langlois-Piazza-Vernizzi theories, which extends the Galileon ghost condensate (GGC) to the domain of beyond Horndeski theories. In this model, we show that the Planck cosmic microwave background (CMB) data, combined with datasets of baryon acoustic oscillations, supernovae type Ia, and redshift-space distortions, give the tight upper bound $|\alpha_{\rm H}^{(0)}| \le {\cal O}(10^{-6})$ on today's beyond-Horndeski (BH) parameter $\alpha_{\rm H}$. This is mostly attributed to the shift of CMB acoustic peaks induced by the early-time changes of cosmological background and perturbations arising from the dominance of $\alpha_{\rm H}$ in the dark energy density. In comparison to the $\Lambda$-cold-dark-matter ($\Lambda$CDM) model, our BH model suppresses the large-scale integrated-Sachs-Wolfe (ISW) tail of CMB temperature anisotropies due to the existence of cubic Galileons, and it modifies the small-scale CMB power spectrum because of the different background evolution. We find that the BH model considered fits the data better than $\Lambda$CDM according to the $\chi^2$ statistics, yet the deviance information criterion (DIC) slightly favors the latter. Given the fact that our BH model with $\alpha_{\rm H}=0$ (i.e., the GGC model) is favored over $\Lambda$CDM even by the DIC, there are no particular signatures for the departure from Horndeski theories in current observations.

Tables

Table 1. Table 1: List of starting values of the density parameters x i subscript 𝑥 𝑖 x_{i} at the redshift z s = 1.5 × 10 5 subscript 𝑧 𝑠 1.5 superscript 10 5 z_{s}=1.5\times 10^{5} and corresponding today’s values for three beyond-Horndeski (BH) models and the f Galileon Ghost Condensate (GGC) model with x 4 = 0 subscript 𝑥 4 0 x_{4}=0 . The BH1, BH2 and BH3 models differ in the starting values x i ( s ) superscript subscript 𝑥 𝑖 𝑠 x_{i}^{(s)} . All of them satisfy theoretically consistent conditions discussed in Sec.  IV . We study these models for the purpose of visualizing and quantifying the modifications from Λ Λ \Lambda CDM. The cosmological parameters (e.g., H 0 , Ω m , Ω r subscript 𝐻 0 subscript Ω 𝑚 subscript Ω 𝑟 H_{0},\Omega_{m},\Omega_{r} ) used for these models are the Planck 2015 best-fit values for Λ Λ \Lambda CDM Ade:2015xua .

Parameters	BH1	BH2	BH3	GGC
$x_1^{\left(s\right)}(\cdot {10}^{-16})$	$-1$	$-0.1$	$-0.01$	$-1$
$x_2^{\left(s\right)}(\cdot {10}^{-16})$	$5$	$0.05$	$0.0001$	$5$
$x_3^{\left(s\right)}(\cdot {10}^{-9})$	$1$	$1$	$0.1$	$10$
$x_4^{\left(s\right)}(\cdot {10}^{-6})$	$100$	$1$	$1$	$0$
$x_1^{\left(0\right)}$	$-1.37$	-1.03	$-0.73$	-1.23
$x_2^{\left(0\right)}$	$2.03$	1.02	$0.12$	1.63
$x_3^{\left(0\right)}$	$0.03$	0.69	$1.30$	0.29
$x_4^{\left(0\right)}$	$1\cdot {10}^{-5}$	$5\cdot {10}^{-6}$	$2\cdot {10}^{-4}$	$0$

Table 2. Table 2: Marginalized values of the model parameters x i ( 0 ) superscript subscript 𝑥 𝑖 0 x_{i}^{(0)} and their 95 % CL bounds, derived by Planck and PBRS datasets. In parenthesis, we also show the maximum likelihood values of these parameters.

Parameters	Planck	PBRS
$x_1^{\left(0\right)}$	$-{1.32}_{-0.12}^{+0.21}\left(-1.25\right)$	$-{1.35}_{-0.06}^{+0.01}\left(-1.25\right)$
$x_2^{\left(0\right)}$	${1.85}_{-0.69}^{+0.33}\left(1.62\right)$	${1.98}_{-0.29}^{+0.14}\left(1.68\right)$
$x_3^{\left(0\right)}$	${0.16}_{-0.18}^{+0.54}\left(0.34\right)$	${0.07}_{-0.1}^{+0.2}\left(0.27\right)$
$x_4^{\left(0\right)}(\cdot {10}^{-6})$	${0.7}_{-1.8}^{+2.2}\left(0.15\right)$	${0.3}_{-0.6}^{+0.7}\left(0.54\right)$

Table 3. Table 3: Marginalized values of H 0 subscript 𝐻 0 H_{0} , σ 8 ( 0 ) superscript subscript 𝜎 8 0 \sigma_{8}^{(0)} , and Ω m ( 0 ) superscript subscript Ω 𝑚 0 \Omega_{m}^{(0)} and their 95 % percent 95 95\% CL bounds in the BH and Λ Λ \Lambda CDM models, derived by Planck and PBRS datasets. The unit of H 0 subscript 𝐻 0 H_{0} is km sec -1  Mpc -1 . In parenthesis, we also show maximum likelihood values of these parameters.

Parameter	Model	Planck	PBRS
	BH	${68.7}_{-2.8}^{+3.2}\left(69.6\right)$	${68.0}_{-1.1}^{+1.1}\left(68.2\right)$
$H_0$
	$\Lambda$CDM	$67.9\pm 2.0\left(67.6\right)$	$68\pm 1\left(68\right)$
	BH	${0.849}_{-0.035}^{+0.037}\left(0.87\right)$	$0.84\pm 0.03\left(0.84\right)$
${\sigma }_8^{\left(0\right)}$
	$\Lambda$CDM	$0.841\pm 0.03\left(0.83\right)$	$0.84\pm 0.03\left(0.84\right)$
	BH	${0.300}_{-0.034}^{+0.033}\left(0.28\right)$	${0.306}_{-0.014}^{+0.014}$ (0.30)
${\Omega }_m^{\left(0\right)}$
	$\Lambda$CDM	$0.30\pm 0.03\left(0.31\right)$	$0.31\pm 0.01\left(0.31\right)$

Table 4. Table 4: Model comparisons in terms of Δ ​ χ eff 2 Δ superscript subscript 𝜒 eff 2 \Delta\chi_{\rm eff}^{2} and Δ ​ DIC Δ DIC \Delta{\rm DIC} . As the reference model, we use the value χ eff 2 superscript subscript 𝜒 eff 2 \chi_{\rm eff}^{2} in Λ Λ \Lambda CDM. The results for GGC are taken from Ref. PBFT .

Model	Dataset	$\Delta {\chi }_{\mathrm{eff}}^2$	$\Delta \mathrm{DIC}$
BH	Planck	-4.7	0.25
BH	PBRS	-1.8	0.1
GGC	Planck	-4.8	-2.5
GGC	PBRS	-2.8	-0.6

Equations

$$\mathcal{S}=\int {\mathrm{d}}^{4}x\sqrt{-g}\sum _{i=2}^4{L}_i+{\mathcal{S}}_M[g_{\mu \nu },{\chi }_M],$$

$$L_2$$

$$L_3$$

$$L_4$$

$${\alpha }_{\mathrm{H}}=-\frac{X^2{F}_4}{G_4-2X{G}_{4,X}+X^2{F}_4},$$

$$\mathrm{d}{s}^2=-\mathrm{d}{t}^2+a^2(t)\left({\delta }_{ij}+h_{ij}\right)\mathrm{d}{x}^i\mathrm{d}{x}^j,$$

$$c_t^2=\frac{G_4}{G_4-2X{G}_{4,X}+X^2{F}_4}.$$

$$F_4=\frac{2G_{4,X}}{X},$$

$$G_2$$

$$G_4$$

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

$$\mathrm{d}{s}^2$$

$$\delta K_{ij}$$

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

$$\mathcal{S}=$$

$$+H_0^2\frac{{\gamma }_1(a)}{2}{\left(\delta g^{00}\right)}^2-H_0\frac{{\gamma }_2(a)}{2}\delta g^{00}\delta K$$

$$-H_0^2\frac{{\gamma }_3(a)}{2}{\left(\delta K\right)}^2-H_0^2\frac{{\gamma }_4(a)}{2}\delta K_j^i\delta K_i^j$$

$$\displaystyle+\frac{\gamma_{5}(a)}{2}\delta g^{00}\delta{\cal R}\bigg{\}}+{\cal S}_{M}[g_{\mu\nu},\chi_{M}]\,,$$

$$3m_0^2{H}^2={\rho }_{\mathrm{DE}}+{\rho }_m+{\rho }_r,$$

$$-m_0^2\left(2\dot{H}+3H^2\right)=P_{\mathrm{DE}}+P_r,$$

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

$$P_{\mathrm{DE}}$$

$${\dot{\rho }}_{\mathrm{DE}}+3H\left({\rho }_{\mathrm{DE}}+P_{\mathrm{DE}}\right)=0.$$

$${\gamma }_3=0,{\gamma }_4=0,$$

$$\Omega =-\frac{2a_4{\dot{\varphi }}^4}{m_0^2},$$

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

$$P_{\mathrm{DE}}$$

$$x_1=-\frac{a_1{\dot{\varphi }}^2}{3m_0^2{H}^2},x_2=\frac{a_2{\dot{\varphi }}^4}{m_0^2{H}^2},$$

$$x_3=\frac{6a_3{\dot{\varphi }}^3}{m_0^{2}H},x_4=\frac{10a_4{\dot{\varphi }}^4}{m_0^2},$$

$${\Omega }_{\mathrm{DE}}=\frac{{\rho }_{\mathrm{DE}}}{3m_0^2{H}^2},{\Omega }_m=\frac{{\rho }_m}{3m_0^2{H}^2},{\Omega }_r=\frac{{\rho }_r}{3m_0^2{H}^2}.$$

$${\Omega }_m=1-{\Omega }_{\mathrm{DE}}-{\Omega }_r,$$

$${\Omega }_{\mathrm{DE}}=x_1+x_2+x_3+x_4.$$

$${\alpha }_{\mathrm{H}}=\frac{4x_4}{5-x_4},$$

$$x_1^{\prime }=2x_1({ϵ}_{\varphi }-h),$$

$$x_2^{\prime }=2x_2(2{ϵ}_{\varphi }-h),$$

$$x_3^{\prime }=x_3(3{ϵ}_{\varphi }-h),$$

$$x_4^{\prime }=4x_4{ϵ}_{\varphi },$$

Full text

I Introduction

Despite the tremendous progress of observational cosmology over the past two decades, the origin of today’s acceleration of the Universe has not been identified yet. The standard concordance scenario is the $\Lambda$CDM model, in which the cosmological constant $\Lambda$ is the source for cosmic acceleration. In addition to the difficulty of naturally explaining the origin of $\Lambda$ from the vacuum energy Weinberg ; Martin ; Tony , it is known that there are tensions between some datasets in the estimations of today’s value of the Hubble constant $H_{0}=100\,h$ km sec*-1* Mpc*-1* Riess:2011yx ; Delubac:2014aqe ; Riess:2016jrr ; Abbott:2017wau ; Abbott:2017smn and the amplitude $\sigma_{8}$ of the matter power spectrum on the scale of $8\,h^{-1}$ Mpc deJong:2015wca ; Kuijken:2015vca ; Hildebrandt:2016iqg ; Conti:2016gav ; Douspis:2018xlj . Such observational tensions along with the theoretical shortcoming of $\Lambda$CDM reinforce the idea to look for alternative models of dark energy dreview1 ; dreview2 ; dreview3 ; dreview4 ; dreview5 ; dreview6 ; Kase:2018aps .

Many extensions to the standard cosmological scenario include large-distance modifications of gravity due to an extra scalar degree of freedom (DOF), thus they are dubbed scalar-tensor theories Fujii:2003pa . Among those, the Horndeski class of theories Horndeski is the most general scheme with second-order equations of motion Def11 ; KYY ; Char11 . The latter feature ensures the absence of Ostrogradski instabilities, related to the existence of higher-order time derivatives.

It is possible to construct healthy theories beyond Horndeski gravity free from Ostrogradski instabilities. In Gleyzes-Langlois-Piazza-Vernizzi (GLPV) theories GLPV , for example, there are two extra Lagrangians beyond the Horndeski domain without increasing the extra propagating DOFs Hami1 ; Hami2 . GLPV theories have several peculiar properties: the propagation speeds of matter and the scalar field are mixed Gergely:2014rna ; Kase:2014yya ; DeFelice:2015isa ; DeFelice:2016ucp , a partial breaking of the Vainshtein mechanism occurs inside astrophysical bodies Kobayashi:2014ida ; Koyama:2015oma ; Sakstein:2015zoa ; Sakstein:2015aac ; Sakstein:2016ggl ; Babichev:2016jom , and a conical singularity can arise at the center of a spherically symmetric and static body DeFelice:2015sya ; Kase:2015gxi . We note that there exist also extensions of Horndeski theories containing higher-order spatial derivatives Gao:2014fra ; Kase:2014cwa ; Frusciante:2015maa (encompassing Horava gravity Horava ) and degenerate higher-order scalar-tensor theories with one scalar propagating DOF Langlois1 ; Langlois2 ; CKT ; Motohashi:2016ftl .

The detection of the gravitational wave (GW) signal GW170817 GW170817 accompanied by the gamma-ray burst event GRB170817A Goldstein shows that the speed of GWs $c_{t}$ is constrained to be in the range $-3\times 10^{-15}\leq c_{t}-1\leq 7\times 10^{-16}$ Abbott at the redshift $z\leq 0.009$, where we use the unit in which the speed of light $c$ is equivalent to 1. The Horndeski Lagrangian, which gives the exact value $c_{t}=1$ without the tuning among functions, is of the form $L=G_{4}(\phi)R+G_{2}(\phi,X)+G_{3}(\phi,X)\square\phi$, where $G_{4}$ is a function of the scalar field $\phi$, $R$ is the Ricci scalar, and $G_{2,3}$ depend on both $\phi$ and $X=\partial^{\mu}\phi\partial_{\mu}\phi$ GWcon1 ; GWcon2 ; GWcon3 ; GWcon4 ; GWcon5 ; GWcon6 . There are also dark energy models in which the GW speed consistent with the above observational bound of $c_{t}$ can be realized Battye:2018ssx ; Amendola:2018ltt ; Copeland:2018yuh . In GLPV theories with the $X$ dependence in $G_{4}$, it is also possible to realize $c_{t}=1$ by the existence of an additional quartic Lagrangian beyond the Horndeski domain Kase:2018iwp .

In addition to the bound on $c_{t}$, the absence of the decay of GWs into dark energy at LIGO/Virgo frequencies ($f\sim 100$ Hz) may imply that the parameter $\alpha_{\rm H}$ characterizing the deviation from Horndeski theories is constrained to be very tiny for the scalar sound speed $c_{s}$ different from 1, typically of order $|\alpha_{\rm H}|\lesssim 10^{-10}$ today Creminelli:2018xsv . If we literally use this bound, there is little room left for dark energy models in beyond-Horndeski theories Frusciante:2018tvu ; Kobayashi:2019hrl . If $c_{s}$ is equivalent to 1, the decay of GWs into dark energy is forbidden. However, it was argued in Ref. Creminelli:2018xsv that power-law divergent terms would appear, leading to the conclusion that the operator accompanying $\alpha_{\rm H}$ must be suppressed as well Creminelli:2018xsv .

We note that the LIGO/Virgo frequencies are close to those of the typical strong coupling scale or cut-off $\Lambda_{c}$ of dark energy models containing derivative field self-interactions deRham:2018red . Around this cut-off scale, we cannot exclude the possibility that some ultra-violet (UV) effects come into play to recover the propagation and property of GWs similar to those in General Relativity (GR). If this kind of UV completion occurs around the frequency $f\sim 100$ Hz, the aforementioned bounds on $c_{t}$ and $\alpha_{\rm H}$ are not applied to the effective field theory of dark energy exploited to describe the cosmological dynamics much below the energy scale $\Lambda_{c}$. Future space-based missions, such as LISA Audley:2017drz , are sensitive to much lower frequencies ($f\sim 10^{-3}$ Hz), so they will offer further valuable information on the properties of GWs with different frequencies.

In GLPV theories, there are constraints on the parameter $\alpha_{\rm H}$ arising from the modifications to gravitational interactions inside astrophysical objects. For example, the consistency of the minimum mass for hydrogen burning in stars with the red dwarf of lowest mass shows that $|\alpha_{\rm H}|$ is at most of order 0.1 Sakstein:2015zoa ; Sakstein:2015aac ; Jain:2015edg ; Saltas:2018mxc . By using X-ray and lensing profiles of galaxy clusters, similar bounds on $\alpha_{\rm H}$ were obtained in Ref. Sakstein:2016ggl . From the orbital period of the Hulse-Taylor binary pulsar PSR B1913+1, the upper bound of $|\alpha_{\rm H}|$ is of order $10^{-3}$ Dima:2017pwp . Cosmological constraints on $\alpha_{\rm H}$ were derived by using particular parametric forms of dimensionless quantities appearing in the effective field theory of dark energy to describe their evolution. In this case, the constraints from CMB and large scale structure data on $|\alpha_{\rm H}|$ are of order $\mathcal{O}(1)$ Traykova:2019oyx .

In this paper, we place observational bounds on the beyond-Horndeski (BH) dark energy model proposed in Ref. Kase:2018iwp and study how the parameter $\alpha_{\rm H}$ is constrained from the cosmological datasets of CMB temperature anisotropies, baryon acoustic oscillations (BAO), supernovae type Ia (SN Ia), and redshift-space distortions (RSDs). In the limit $\alpha_{\rm H}\to 0$, the model reduces to the Galileon ghost condensate (GGC) in Horndeski theories. The recent analysis of Ref. PBFT reveals that the GGC model is observationally favored over $\Lambda$CDM according to several information criteria. We will investigate whether or not this property persists for the BH dark energy model ($\alpha_{\rm H}\neq 0$) of Ref. Kase:2018iwp . For the likelihood analysis, we will use the publicly available Effective-Field-Theory for CAMB (EFTCAMB) code111Web page: http://www.eftcamb.org Hu:2013twa ; Raveri:2014cka . In our investigation the gravitational theory is completely determined by a covariant action, while the analysis in Ref. Traykova:2019oyx follows a parameterized approach to GLPV theories. In this respect, the two cosmological models considered are completely different and the constraint on $\alpha_{\rm H}$ obtained in this paper cannot be straightforwardly compared to the results in Ref. Traykova:2019oyx .

The paper is organized as follows. In Sec II, we briefly review the basics of the BH dark energy model introduced in Ref. Kase:2018iwp . In Sec. III, we show how this model can be implemented in the EFT formulation and derive the background equations of motion together with theoretically consistent conditions. In Sec. IV, we discuss the evolution of cosmological perturbations in the presence of matter perfect fluids and investigate the impact of our model on observable quantities. In Sec. V, we present the Monte-Carlo-Markov-Chain (MCMC) constraints on model parameters and compute several information criteria to discuss whether the BH model is favored over the $\Lambda$CDM model. Finally, we conclude in Sec. VI.

II Dark energy model in GLPV theories

The dark energy model proposed in Ref. Kase:2018iwp belongs to the quartic-order GLPV theories given by the action

[TABLE]

where $g$ is the determinant of metric tensor $g_{\mu\nu}$, ${\cal S}_{M}$ is the matter action for all matter fields $\chi_{M}$, and the Lagrangians $L_{2,3,4}$ are defined by

[TABLE]

where $G_{2,3,4}$ and $F_{4}$ are functions of the scalar field $\phi$ and $X=\nabla^{\mu}\phi\nabla_{\mu}\phi$, $R$ is the Ricci scalar, and $\epsilon^{\mu\nu\rho\sigma}$ is the totally antisymmetric Levi-Civita tensor satisfying the normalization $\epsilon^{\mu\nu\rho\sigma}\epsilon_{\mu\nu\rho\sigma}=+4!$. We also define $G_{i,X}\equiv\partial G_{i}/\partial X$ and use the notations $\phi_{\mu}=\nabla_{\mu}\phi$ and $\phi_{\mu\nu}=\nabla_{\nu}\nabla_{\mu}\phi$ for the covariant derivative operator $\nabla_{\mu}$. We assume that the matter fields $\chi_{M}$ are minimally coupled to gravity.

The last term containing $F_{4}(\phi,X)$ in $L_{4}$ arises beyond the domain of Horndeski theories GLPV . The deviation from Horndeski theories can be quantified by the parameter

[TABLE]

which does not vanish for $F_{4}\neq 0$. The line element containing intrinsic tensor perturbations $h_{ij}$ on the flat Friedmann-Lemaître-Robertson-Walker (FLRW) space-time is given by

[TABLE]

where $a(t)$ is the time-dependent scale factor, and $h_{ij}$ satisfies the transverse and traceless conditions ($\nabla^{j}h_{ij}=0$ and ${h_{i}}^{i}=0$). The propagation speed squared of tensor perturbations is GLPV ; Gergely:2014rna ; Kase:2014yya

[TABLE]

In quartic-order Horndeski theories ($F_{4}=0$), the $X$ dependence in $G_{4}$ leads to the difference of $c_{t}^{2}$ from 1. In GLPV theories, it is possible to realize $c_{t}^{2}=1$ for the function

[TABLE]

under which $\alpha_{\rm H}=-2XG_{4,X}/G_{4}$.

In this paper, we will study observational constraints on the model proposed in Ref. Kase:2018iwp . This is characterized by the following functions

[TABLE]

where $m_{0}$ and $a_{1,2,3,4}$ are constants. This beyond-Horndeski model, hereafter BH, satisfies the condition (6), and hence $c_{t}^{2}=1$. When $a_{4}=0$, BH recovers the GGC model studied recently in Ref. PBFT . Taking the limits $a_{2}\to 0$ and $a_{3}\to 0$, GGC recovers the cubic covariant Galileon Nicolis ; Galileons and ghost condensate Arkani , respectively.

The BH model allows for the existence of self-accelerating de Sitter solutions finally approaching constant values of $X$. Before approaching the de Sitter attractor, the dark energy equation of state $w_{\rm DE}$ can exhibit a phantom behavior (i.e., $w_{\rm DE}<-1$) without the appearance of ghosts Kase:2018iwp . The cubic covariant Galileon gives rise to the tracker solution with $w_{\rm DE}=-2$ in the matter era DT10 , but this evolution is incompatible with the joint data analysis of CMB, BAO, and SN Ia NDT10 . On the other hand, in both BH and GGC, the $a_{2}X^{2}$ term works to prevent for approaching the tracker, so that $-2<w_{\rm DE}<-1$ in the matter era. This behavior of $w_{\rm DE}$ is consistent with the recent observational datasets of CMB, BAO, and SN Ia PBFT .

The BH model leads to the evolution of cosmological perturbations different from that in GR. The late-time modification to the cosmic growth rate arises mostly from the cubic Galileon term $3a_{3}X\square\phi$ Kase:2018iwp ; Bartolo . In GGC, the combined effect of $3a_{3}X\square\phi$ and $a_{2}X^{2}$ can suppress the power spectrum of large-scale CMB temperature anisotropies, so that the model shows a better compatibility with the Planck data with respect to the $\Lambda$CDM PBFT . It remains to be seen whether the similar property also holds for the BH model with $a_{4}\neq 0$, which we will address in this paper.

III Methodology

In this section, we discuss the evolution of the background and linear scalar perturbations in the BH model. We make use of the EFTCAMB/EFTCosmoMC codes Hu:2013twa ; Raveri:2014cka , in which the EFT of dark energy and modified gravity Gubitosi:2012hu ; Bloomfield:2013efa ; Bloomfield:2012ff ; Gleyzes:2013ooa ; Gleyzes:2014rba is implemented into CAMB/CosmoMC Lewis:1999bs ; Lewis:2002ah . The EFT framework enables one to deal with any dark energy and modified gravity model with one scalar propagating DOF $\phi$ in a unified and model-independent manner.

The EFT of dark energy is based on the 3+1 Arnowitt-Deser-Misner (ADM) decomposition of spacetime ADM given by the line element

[TABLE]

where $N$ is the lapse, $N^{i}$ is the shift vector, and $h_{ij}$ is the three-dimensional metric. A unit vector orthogonal to the constant time hyper-surface $\Sigma_{t}$ is given by $n_{\mu}=N\nabla_{\mu}t=(N,0,0,0)$. The extrinsic curvature is defined by $K_{ij}=h^{k}_{i}\nabla_{k}n_{j}$. The internal geometry of $\Sigma_{t}$ is quantified by the three-dimensional Ricci tensor ${\cal R}_{ij}={}^{(3)}R_{ij}$ associated with the metric $h_{ij}$.

On the flat FLRW background, we consider the line element containing three scalar metric perturbations $\delta N$, $\psi$, and $\zeta$, as

[TABLE]

where $\partial_{i}\equiv\partial/\partial x^{i}$. We also choose the the unitary gauge in which the perturbation $\delta\phi$ of the scalar field $\phi$ vanishes. Then, the perturbations of extrinsic and intrinsic curvatures are expressed as Gleyzes:2013ooa ; Gergely:2014rna ; Kase:2014cwa ; Gleyzes:2014rba

[TABLE]

where $\partial^{2}\equiv\delta^{kl}\partial_{k}\partial_{l}$, and $H=\dot{a}/a$ is the Hubble expansion rate, and a dot represents a derivative with respect to $t$. The perturbations of traces $K\equiv{K^{i}}_{i}$ and ${\cal R}\equiv{{\cal R}^{i}}_{i}$ are denoted as $\delta K$ and $\delta{\cal R}$, respectively, with $\delta g^{00}=2\delta N$.

In the ADM language, the Lagrangian of GLPV theories depends on the scalar quantities $N$, $K$, $K_{ij}K^{ij}$, ${\cal R}$, $K_{ij}{\cal R}^{ij}$, and $t$ Gleyzes:2013ooa . Expanding the corresponding action up to second order in scalar perturbations of those quantities, it follows that

[TABLE]

where $m_{0}$ is a constant having a dimension of mass, and $\Omega$, $\Lambda$, $c$, $\gamma_{i}$ are called EFT functions that depend on the background scale factor $a$. The explicit relations between those EFT functions and the functions $G_{2,3,4},F_{4}$ in the action (1) are given in Ref. Frusciante:2016xoj .

The first three variables $\Omega$, $\Lambda$, $c$ determine both the background evolution and linear perturbations, whereas the functions $\gamma_{i}$ solely appear at the level of linear perturbations. For the matter action ${\cal S}_{M}$, we take dark matter and baryons (background density $\rho_{m}$ and vanishing pressure) and radiation (background density $\rho_{r}$ and pressure $P_{r}=\rho_{r}/3$) into account. Then, the background equations are expressed as Gubitosi:2012hu ; Bloomfield:2012ff

[TABLE]

where

[TABLE]

The density $\rho_{\rm DE}$ and pressure $P_{\rm DE}$ of dark energy obey the continuity equation

[TABLE]

In GLPV theories, there is the specific relation $\gamma_{3}=-\gamma_{4}$. If we restrict the theories to those satisfying $c_{t}^{2}=1$, it follows that $\gamma_{4}=0$. Then, the model given by the functions (II) corresponds to the coefficients

[TABLE]

so that we are left with three functions $\gamma_{1},\gamma_{2},\gamma_{5}$ at the level of linear perturbations.

To study the cosmological evolution of our model in EFTCAMB, we first solve the background equations of motion and then map to the EFT functions according to the procedure given in Refs. Gubitosi:2012hu ; Bloomfield:2012ff ; Bloomfield:2013efa ; Gleyzes:2013ooa ; Gleyzes:2014rba ; Frusciante:2015maa ; Frusciante:2016xoj .

III.1 Background equations in the BH model

For the model (II), the background equations are given by Eqs. (13) and (14), with

[TABLE]

and

[TABLE]

The parameters $c$ and $\Lambda$ in Eqs. (15) and (16) can be expressed in terms of quantities on the right hand sides of Eqs. (20) and (21). Following Ref. Kase:2018iwp , we define the dimensionless variables (density parameters):

[TABLE]

and

[TABLE]

From Eq. (13), we have

[TABLE]

where the dark energy density parameter is given by

[TABLE]

In terms of $x_{4}$, the deviation parameter (3) from Horndeski theories is expressed as

[TABLE]

and hence $\alpha_{\rm H}$ is of the same order as $x_{4}$ for $|x_{4}|\leq 1$.

The variables $x_{1,2,3,4}$ and $\Omega_{r}$ are known by solving the ordinary differential equations

[TABLE]

where a prime denotes the derivative with respect to ${\cal N}=\ln(a)$. On using Eqs. (13) and (14), it follows that

[TABLE]

with

[TABLE]

For a given set of initial conditions $x_{1,2,3,4}$ and $\Omega_{r}$, we can solve Eqs. (27)-(31) to determine the evolution of density parameters as well as $\phi$ and $H$. Practically, we start to solve the above dynamical system at redshift $z_{s}=1.5\times 10^{5}$ and iteratively scan over initial conditions leading to the viable cosmology satisfying the constraint (24) today ($z=0$). Additionally, evaluating Eq. (25) at present time, we can eliminate one model parameter, for example $x_{2}^{(0)}$, as $x_{2}^{(0)}=\Omega_{\rm DE}^{(0)}-x_{1}^{(0)}-x_{3}^{(0)}-x_{4}^{(0)}$, where “(0)” represents today’s quantities.

III.2 Mapping

To study the evolution of scalar perturbations and observational constraints on dark energy models in EFTCAMB, it is convenient to use the mapping between EFT functions and model parameters in BH. In Sec. III.1, we already discussed the mapping of the background quantities $\Omega,\Lambda$ and $c$. The functions $\gamma_{1,2,5}$, which are associated with scalar perturbations, are given by

[TABLE]

The expressions of these EFT functions allow us to draw already some insight about the contributions of each $x_{i}$ to the dynamics of linear perturbations. In general, the variable $\gamma_{1}$ cannot be well constrained by data being its contribution to the observables below the cosmic variance Frusciante:2018jzw . The main modification to the evolution of perturbations compared to GR arises from $\gamma_{2}$ and $\gamma_{5}$, which are mostly affected by $x_{3}$ and $x_{4}$. The variables $x_{1}$ and $x_{2}$ contribute to the perturbation dynamics through the Hubble expansion rate $H$ in $\gamma_{2}$.

III.3 Viability constraints

There are theoretically consistent conditions under which the perturbations are not plagued by the appearance of ghosts and Laplacian instabilities in the small-scale limit. For the BH model (II), the conditions for the absence of ghosts in tensor and scalar sectors are given, respectively, by Kase:2018iwp

[TABLE]

where $q_{s}$ is defined in Eq. (32). Then, we have the following constraints

[TABLE]

The BH model has the property $c_{t}^{2}=1$, so there is no Laplacian instability for tensor perturbations. We note that the reduced Planck mass $M_{\rm pl}$ is related to $m_{0}$ according to the relation $M_{\rm pl}^{2}=m_{0}^{2}(1+\Omega_{0})$ in the local environment with screened fifth forces Kobayashi:2014ida , where $\Omega_{0}$ is today’s value of $\Omega$. Then, the Newton gravitational constant $G_{\rm N}$ is given by

[TABLE]

which is positive under the absence of tensor ghosts.

For scalar perturbations, there are three propagation speed squares $c_{s}^{2}$, $\tilde{c}_{r}^{2}$, and $\tilde{c}_{m}^{2}$ associated with the scalar field $\phi$, radiation, and nonrelativistic matter, respectively. In Horndeski theories, they are not coupled to each other, so that the propagation speed squares of radiation and nonrelativistic matter are given, respectively, by $c_{r}^{2}=1/3$ and $c_{m}^{2}=+0$. In GLPV theories, they are generally mixed with each other, apart from $\tilde{c}_{m}^{2}$ (which has the value $\tilde{c}_{m}^{2}=+0$) GLPV ; Gergely:2014rna ; Kase:2014yya ; DeFelice:2015isa ; DeFelice:2016ucp . Then, the Laplacian instabilities of scalar perturbations can be avoided under the two conditions

[TABLE]

where

[TABLE]

When $|\alpha_{\rm H}|\ll 1$ we have $c_{s}^{2}\simeq c_{\rm H}^{2}-\beta_{\rm H}$ and $\tilde{c}_{r}^{2}\simeq c_{r}^{2}=1/3$, so the second stability condition (41) is satisfied.

There are also constraints on today’s parameter $\alpha_{\rm H}^{(0)}$ (or equivalently, $x_{4}^{(0)}$) from massive astrophysical objects Sakstein:2015zoa ; Sakstein:2015aac ; Dima:2017pwp . Among those constraints, the orbital period of Hulse-Taylor binary pulsar gives the tightest bound $-0.0031\leq x_{4}^{(0)}\leq 0.0094$ Dima:2017pwp ; Kase:2018iwp . If we literally use the bound arising from the absence of the GW decay into dark energy at LIGO/Virgo frequencies, the parameter $\alpha_{\rm H}^{(0)}$ should be less than the order of $10^{-10}$ Creminelli:2018xsv . As we mentioned in Introduction, it is still a matter of debate whether the EFT of dark energy is valid around the frequency $f\sim 100$ Hz deRham:2018red . In this paper, we will not impose such a bound and independently test how the cosmological observations place the upper limit of $x_{4}^{(0)}$.

In Fig. 1, we show the physically viable parameter space (blue colored region) for the initial conditions $x_{1}^{(s)}$, $x_{2}^{(s)}$, $x_{3}^{(s)}$, $x_{4}^{(s)}$ (at redshift $z_{s}=1.5\times 10^{5}$) and today’s values $x_{1}^{(0)},x_{2}^{(0)},x_{3}^{(0)}$ (at redshift $z=0$). We find that $x_{1}^{(0)}$ is negative, while $x_{2}^{(0)}$ and $x_{3}^{(0)}$ are positive. We note that the ghost condensate model Arkani has a de Sitter solution satisfying $x_{1}<0$ and $x_{2}>0$. The Galileon term $x_{3}$ modifies the cosmological dynamics of ghost condensate, but there is also a de Sitter attractor characterized by $x_{1}<0$, $x_{2}>0$, and $x_{3}>0$ Kase:2018iwp . As we see in Fig. 1, the parameter $x_{3}^{(0)}$ is not well constrained from the theoretically viable conditions alone.

The parameter space of the variable $x_{4}^{(0)}$ is not shown in Fig. 1, but it is in the range $|x_{4}^{(0)}|\ll 1$ to satisfy all the theoretically consistent conditions. As $x_{4}^{(0)}$ approaches the order 1, the scalar perturbation is typically prone to the Laplacian instability associated with the negative value of $c_{s}^{2}$ Kase:2018iwp .

The above results will be used to set theoretical priors for the MCMC analysis.

IV Cosmological perturbations

In this section, we discuss the evolution of scalar cosmological perturbations in the BH model for the perturbed line element given by Eq. (9). We introduce the two gauge-invariant gravitational potentials:

[TABLE]

For the matter sector, we consider scalar perturbations of the matter-energy momentum tensor $T^{\mu}_{\nu}$ arising from the action ${\cal S}_{M}$, as $\delta T^{0}_{0}=-\delta\rho$, $\delta T^{0}_{i}=\partial_{i}\delta q$, and $\delta T^{i}_{j}=\delta P\delta^{i}_{j}$. The density perturbation $\delta\rho$, the momentum perturbation $\delta q$, and the pressure perturbation $\delta P$ are expressed in terms of the sum of each matter component, as $\delta\rho=\sum_{i}\delta\rho_{i}$, $\delta q=\sum_{i}\delta q_{i}$, and $\delta P=\sum_{i}\delta P_{i}$, where $i=m,r$. We introduce the gauge-invariant density contrast:

[TABLE]

where $\rho_{i}$ is the background density of each component. In the BH model, the full linear perturbation equations of motion were derived in Ref. Kase:2018iwp .

In Fourier space with the comoving wavenumber $k$, we relate the gravitational potentials in Eq. (43) with the total matter density contrast $\Delta=\sum_{i}\Delta_{i}$, as Amendola:2007rr ; Bertschinger:2008zb ; Pogosian:2010tj

[TABLE]

where $G_{\rm N}$ is the Newton gravitational constant given by Eq. (39), and $\rho=\sum_{i}\rho_{i}$ is the total background matter density. The dimensionless quantities $\mu$ and $\Sigma$ correspond to the effective gravitational couplings felt by matter and light, respectively. For nonrelativistic matter, the density contrast $\Delta_{m}$ obeys Kase:2018iwp

[TABLE]

where ${\cal B}\equiv\zeta+H\delta q_{m}/\rho_{m}$. This means that the matter density contrast grows due to the gravitational instability through the modified Poisson Eq. (45). In GR, both $\mu$ and $\Sigma$ are equivalent to 1, but in the BH model, they are different from 1. Hence the growth of structures and gravitational potentials is subject to modifications.

For the perturbations deep inside the sound horizon ($c_{s}^{2}k^{2}/a^{2}\gg H^{2}$), the common procedure is to resort to a quasi-static approximation for the estimations of $\mu$ and $\Sigma$ Pola00 ; Tsujikawa:2007gd ; DKT . This amounts to picking up the terms containing $k^{2}/a^{2}$ and $\Delta_{m}$ in the perturbation equations of motion. In Horndeski theories, it is possible to obtain the closed form expressions of $\Psi,\Phi,\zeta$ DKT ; Kase:2018aps . In GLPV theories, the additional time derivatives $\alpha_{\rm H}\dot{\psi}$ and $\alpha_{\rm H}\dot{\zeta}$ appear even under the quasi-static approximation DeFelice:2015isa ; Tsujikawa:2015mga , so the perturbation equations are not closed. If $|\alpha_{\rm H}|$ is very much smaller than 1 and $x_{4}$ is subdominant to $x_{1,2,3}$, we may ignore the contributions of the term $x_{4}$ to the perturbation equations. In this case, we can estimate $\mu$ and $\Sigma$ in the BH model, as Kase:2018iwp

[TABLE]

Since $\mu$ and $\Sigma$ are identical to each other, it follows that $\Psi\simeq\Phi$. Under the theoretically consistent conditions (36), (37), and (40), we also have $\mu\simeq\Sigma>1$ and hence the gravitational interaction is stronger than that in GR. Let us note that in the following we will not rely on this approximation and we will solve the complete linear perturbation equations.

To understand the evolution of perturbations, we consider four different cases (BH1, BH2, BH3, GGC) listed in Table 1. The difference between these models is characterized by the different choices of initial conditions $x_{i}^{(s)}$ at the redshift $z_{s}=1.5\times 10^{5}$. Among them, BH1 has the largest initial value $x_{4}^{(s)}$, while $x_{4}$ is always 0 in GGC (which belongs to Horndeski theories). In Fig. 2, we plot the evolution of $x_{i}$ from the past to today for these four different cases. In BH1, the variable $x_{4}$ dominates over other variables $x_{1,2,3}$ for $a\lesssim 10^{-2}$, but it becomes subdominant at low redshifts with today’s value of order $10^{-5}$. Comparing BH1 with BH3, we observe that the initial largeness of $x_{4}$ does not necessarily imply the large present-day value $x_{4}^{(0)}$. At low redshifts, $x_{4}$ is typically less than the order $10^{-3}$ to avoid $c_{s}^{2}<0$ with the amplitude smaller than $x_{1,2,3}$, in which case the analytic estimation (48) can be trustable. Indeed, for all the models given in Table 1, we numerically checked that the quasi-static approximation holds with sub-percent precision for the wavenumbers $k>0.01$ Mpc*-1* (as confirmed in Horndeski theories in Ref. Peirone:2017ywi ; Frusciante:2018jzw ).

In the top panel of Fig. 3, we plot the evolution of $\Psi$ normalized by its initial value $\Psi^{(s)}$ for the four models in Table 1 and for the $\Lambda$CDM. In the bottom panel, we depict the percentage difference of $\Psi$ for the chosen models with respect to $\Lambda$CDM. At the late epoch, the deviations from $\Lambda$CDM show up with the enhanced gravitational potential (around $a\sim 0.2$ for the BH2, BH3, GGC models). The largest deviation arises for BH3, in which case the difference is more than 75 % today. As estimated from Eq. (48), the modified evolution of $\Psi$ is mostly attributed to the cubic Galileon term $x_{3}$. For larger today’s values of $x_{3}^{(0)}$, the difference of $\Psi$ from $\Lambda$CDM tends to be more significant with the larger deviation of $\mu$ from 1. In Fig. 3, we observe that the deviation from $\Lambda$CDM increases with the order of BH1, GGC, BH2, BH3, by reflecting their increasing values of $x_{3}^{(0)}$ given in Table 1.

In BH1, there is the suppression of $|\Psi|$ in comparison to $\Lambda$CDM at high redshifts ($a\lesssim 10^{-2}$). This property arises from the dominance of $x_{4}$ over $x_{1,2,3}$ at early times, in which case the relative density abundances between dark energy and matter fluids are modified. Besides this effect, the non-negligible early-time contribution of $x_{4}$ to scalar perturbations gives rise to a scale-dependent evolution of gravitational potentials, which manifests itself in the $k$-dependent variation of $\mu(a,k)$ and $\Sigma(a,k)$. In Fig. 4, we plot the evolution of $\Psi$ in BH1 for three different values of $k$. For perturbations on smaller scales, the deviation from $\Lambda$CDM tends to be more significant. In models BH2, BH3, and GGC, the early-time evolution of $\Psi$ is similar to that in $\Lambda$CDM, but they exhibit large deviations from $\Lambda$CDM at late times.

At low redshifts, the lensing gravitational potential $\phi_{\rm len}=(\Psi+\Phi)/2$ evolves in a similar way to $\Psi$, by reflecting the property $\mu\simeq\Sigma$ for $x_{4}^{(0)}\ll 1$. The lensing angular power spectrum can be computed by using the line of sight integration method, with the convention Lewis:2006fu

[TABLE]

where $\mathcal{P}(k)=\Delta^{2}_{\mathcal{R}}(k)$ is the primordial power spectrum of curvature perturbations, and $j_{\ell}$ is the spherical Bessel function. The source $S_{\phi}$ is expressed in terms of the transfer function

[TABLE]

with $T_{\phi}(k,\tau)=k\phi_{\rm len}$, $\chi$ is the comoving distance with $\chi_{\ast}$ corresponding to that to the last scattering surface, $\tau_{0}$ is today’s conformal time $\tau=\int a^{-1}{\rm d}t$ satisfying the relation $\chi=\tau_{0}-\tau$. In Fig. 5, we show the lensing power spectra $D_{\ell}^{\phi\phi}=\ell(\ell+1)C_{\ell}^{\phi\phi}/(2\pi)$ and relative differences in units of the cosmic variance for four models listed in Table 1. Since $\Sigma>1$ at low redshifts in BH and GGC models, this works to enhance $D_{\ell}^{\phi\phi}$ compared to $\Lambda$CDM. We note that the amplitude of matter density contrast $\delta_{m}$ in these models also gets larger than that in $\Lambda$CDM by reflecting the fact that $\mu>1$. In Fig. 5, we observe that, apart from BH1 in which $\Sigma$ is close to 1, the lensing power spectra in other three cases are subject to the enhancement with respect to $\Lambda$CDM. Since today’s values of $\mu$ and $\Sigma$ increase for larger $x_{3}^{(0)}$, the deviation from $\Lambda$CDM tends to be more significant with the order of GGC, BH2, and BH3.

Let us proceed to the discussion of the impact of BH and GGC models on the CMB temperature anisotropies. The CMB temperature-temperature (TT) angular spectrum can be expressed as Seljak:1996is

[TABLE]

where

[TABLE]

with $\tilde{\mu}$ being the angular separation, and $S_{\rm T}(k,\tau)$ is the radiation transfer function. The contribution to $S_{\rm T}(k,\tau)$ arising from the integrated-Sachs-Wolfe (ISW) effect is of the form

[TABLE]

where $\kappa$ is the optical depth. Besides the early ISW effect which occurs during the transition from the radiation to matter eras by the time variation of $\Psi+\Phi$, the presence of dark energy induces the late-time ISW effect. In the $\Lambda$CDM model, the gravitational potential $-(\Psi+\Phi)$, which is positive, decreases by today with at least more than 30 % relative to its initial value (see Fig. 3). As we observe in Fig. 6 we have $\dot{\Psi}+\dot{\Phi}>0$ in this case, so the ISW effect gives rise to the positive contribution to Eq. (51). In Fig. 7, we plot the CMB TT power spectra $D_{\ell}^{\rm TT}=\ell(\ell+1)C_{\ell}^{\rm TT}/(2\pi)$ for the models listed in Table 1 and $\Lambda$CDM. In BH1 the parameter $\Sigma$ is close to 1 at low redshifts due to the smallness of $x_{3}^{(0)}$, so the late-time ISW effect works in the similar way to the GR case. Hence the TT power spectrum in BH1 for the multipoles $\ell\lesssim 30$ is similar to that in $\Lambda$CDM.

In the GGC model of Fig. 7, we observe that the large-scale ISW tail is suppressed relative to that in $\Lambda$CDM. This reflects the fact that the larger deviation of $\Sigma$ from 1 leads to the time derivative $\dot{\Psi}+\dot{\Phi}$ closer to 0, see Fig. 6. Hence the late-time ISW effect is not significant, which results in the suppression of $D_{\ell}^{\rm TT}$ with respect to $\Lambda$CDM. In Ref. PBFT this fact was first recognized in the GGC model, which exhibits a better fit to the Planck CMB data. As the deviation of $\Sigma$ from 1 increases further, the sign of $\dot{\Psi}+\dot{\Phi}$ changes to be negative (see Fig. 6). The BH2 model can be regarded as such a marginal case in which the large-scale ISW tail is nearly flat. In BH3, the increase of $\Sigma$ at low redshifts is so significant that the largely negative ISW contribution to Eq. (51) leads to the enhanced low-$\ell$ TT power spectrum relative to $\Lambda$CDM.

The modified evolution of the Hubble expansion rate from $\Lambda$CDM generally leads to the shift of CMB acoustic peaks at high-$\ell$. In Fig. 8, we observe that the largest deviation of $H(a)$ at high redshifts occurs for BH1 by the dominance of $x_{4}$ over $x_{1,2,3}$. This leads to the shift of acoustic peaks toward lower multipoles (see Fig. 7). We also find that BH3 is subject to non-negligible shifts of high-$\ell$ peaks due to the large modification of $H(a)$ at low redshifts, in which case the peaks shift toward higher multipoles. Moreover, there is the large enhancement of ISW tails for BH3, so it should be tightly constrained from the CMB data. We note that the shift of CMB acoustic peaks is further constrained by the datasets of BAO and SN Ia. For BH2 and GGC the changes of peak positions are small in comparion to BH1 and BH3, but still they are in the range testable by the CMB data. Moreover, the large-scale ISW tail is subject to the suppression relative to $\Lambda$CDM in BH2 and GGC.

In BH1, we also notice a change in the amplitude of acoustic peaks occurring dominantly at high $\ell$. This is known to be present in models with early-time modifications of gravity Xiang-Raveri ; Benevento2018 . The modification of gravitational potentials affects the evolution of radiation perturbations (monopole and dipole) through the radiation driving effect Hu:1994uz ; Xiang-Raveri , thus resulting in the changes in amplitude and phase of acoustic peaks at high $\ell$.

The modified time variations of $\Psi$ and $\Phi$ around the recombination epoch also give a contribution to the early ISW effect. This is important on scales around the first acoustic peak, corresponding to the wavenumber $k\simeq 0.016$ Mpc*-1* for our choice of model parameters. To have a more qualitative feeling of this contribution, we have estimated the impact of the early ISW effect on $D_{\ell}^{\rm TT}$ by using the approximate ISW integral presented in Ref. Hu:1994uz :

[TABLE]

where $\tau_{*}$ is the conformal time at the last scattering. Then, we find a negative difference of about 4.9 % between BH1 and $\Lambda$CDM. This is in perfect agreement with the change in amplitude of the first acoustic peak shown in Fig. 7. Thus, the BH models in which $x_{4}$ is the dominant contribution to the dark energy dynamics at early times can be severely constrained from the CMB data.

We stress that, in the late Universe, $x_{4}$ is typically suppressed compared to $x_{1,2,3}$ for the viable cosmological background, so the main impact on the evolution of perturbations comes from the cubic Galileon term $x_{3}$. The analytic estimation (48) is sufficiently trustable for studying the evolution of gravitational potentials and matter perturbations at low redshifts. However, we solve the full perturbation equations of motion for the MCMC analysis without resorting to the quasi-static approximation.

V Observational constraints

We place observational bounds on the BH model by performing the MCMC simulation with different combinations of datasets at high and low redshifts.

V.1 Datasets

For the MCMC likelihood analysis, based on the EFTCosmoMC code, we use the Planck 2015 (Aghanim:2015xee, ; Ade:2015xua, ) data of CMB temperature and polarization on large angular scales, for multipoles $\ell<29$ (low-$\ell$ TEB likelihood) and the CMB temperature on smaller angular scales (PLIK TT Likelihood). We also consider the BAO measurements from the 6dF Galaxy Survey (BAO1, ) and from the SDSS DR7 Main Galaxy Sample (BAO2, ). Moreover, we include the combined BAO and RSD datasets from the SDSS DR12 consensus release (Alam:2016hwk, ) and the JLA SN Ia sample (Betoule:2014frx, ). We will refer to the full combined datasets as “PBRS”.

Finally, we impose the flat priors on the model parameters: $x_{1}^{(s)}\in[-10,10]\times 10^{-16}$, $x_{3}^{(s)}\in[-10,10]\times 10^{-9}$, and $x_{4}^{(s)}\in[0,10]\times 10^{-6}$. Even by increasing the prior volume by one order of magnitude, we confirmed that the likelihood results are not subject to the priors choice.

V.2 Constrained parameter space

In this section, we show observational constraints on model parameters in the BH model. We use the datasets presented in Sec. V.1 with two combinations: (i) Planck and (ii) PBRS. For reference, we also present the results of the $\Lambda$CDM model.

In Table 2, we show the marginalized values of today’s four density parameters $x_{i}^{(0)}$ with 95 % confidence level (CL) limits. In Fig. 9, we plot the observationally allowed regions derived by two combinations of datasets with the $68\%$ and $95\%$ CL boundaries. The best-fit values of $x_{1}^{(0)}$ and $x_{2}^{(0)}$ constrained by the Planck data are not affected much by including the datasets of BAO, SN Ia, and RSDs. In the observationally allowed region we have $x_{1}^{(0)}<0$ and $x_{2}^{(0)}>0$, but there are neither ghosts nor Laplacian instabilities in the constrained parameter space (as in the ghost condensate model Arkani ).

With the Planck data alone, the 95 % CL upper bound on $x_{3}^{(0)}$ is close to 1, but the PBRS datasets give the tighter limit $x_{3}^{(0)}\leq 0.27$ at 95 % CL. The maximum likelihood value of $x_{3}^{(0)}$ derived with the Planck data is 0.34, which is similar to the corresponding value 0.27 constrained with PBRS. The non-vanishing best-fit value of $x_{3}^{(0)}$ is attributed to the facts that, relative to $\Lambda$CDM, (i) the Galileon term can suppress the large-scale ISW tale, and (ii) the modified background evolution gives rise to the TT power spectrum showing a better fit to the Planck CMB data at high-$\ell$. In Fig. 10, these properties can be seen in the best-fit TT power spectrum of the BH model. Increasing $x_{3}^{(0)}$ further eventually leads to the enhancement of the ISW tale in comparison to $\Lambda$CDM. As we see in BH3 of Fig. 7, the models with large $x_{3}^{(0)}$ do not fit the TT power spectrum well at high-$\ell$ either. Such models are disfavored from the CMB data (as in the case of covariant Galileons Renk ; Peirone ), so that $x_{3}^{(0)}$ is bounded from above. The RSD data at low redshifts can be also consistent with the intermediate values of $x_{3}^{(0)}$ constrained from CMB.

In Fig. 11, we show the evolution of $w_{\rm DE}$ for the best-fit BH model. As discussed in Ref. Kase:2018iwp , the existence of $x_{2}$ besides $x_{3}$ prevents the approach to a tracker solution characterized by $w_{\rm DE}=-2$ during the matter-dominated epoch. The best-fit background solution first enters the region $-2<w_{\rm DE}<-1$ in the matter era and finally approaches a de Sitter attractor characterized by $w_{\rm DE}=-1$. Thus, the BH and GGC models with $x_{2}\neq 0$ alleviate the observational incompatibility problem of tracker solutions of covariant Galileons NDT10 . For the best-fit BH model, there is the deviation of $w_{\rm DE}$ from $-1$ with the value $w_{\rm DE}\approx-1.1$ at the redshift $1<z<3$, so the model is different from $\Lambda$CDM even at the background level.

From the PBRS datasets, today’s value of $x_{4}$ is constrained to be

[TABLE]

so that $|x_{4}^{(0)}|$ is at most of order $10^{-6}$. With the Planck data alone, the upper bound of $|x_{4}^{(0)}|$ is also of the same order. This means that the upper limit of $x_{4}^{(0)}$ is mostly determined by the CMB data. As we discussed in Sec. IV, the CMB TT power spectrum is sensitive to the dominance of $x_{4}$ over $x_{1,2,3}$ in the early cosmological epoch. Then, today’s value of $x_{4}$ is also tightly constrained as Eq. (55), which translates to the bound

[TABLE]

Apart from the constraint arising from the GW decay to dark energy Creminelli:2018xsv , the above upper limit on $\alpha_{\rm H}^{(0)}$ is the most stringent bound derived from cosmological observations so far.

In Table 3, we present the values of $H_{0}$, $\sigma_{8}^{(0)}$, and $\Omega_{m}^{(0)}$ constrained from the datasets of Planck and PBRS for the BH and $\Lambda$CDM models. The bounds on $H_{0}$, $\sigma_{8}^{(0)}$, and $\Omega_{m}^{(0)}$ derived with the PBRS datasets are similar to those in $\Lambda$CDM. In Fig. 12, we also plot the two-dimensional observational contours for these parameters constrained by the Planck data. The direct measurements of $H_{0}$ at low redshifts Riess:2018uxu give the bound $H_{0}>70$ km sec*-1* Mpc*-1*, whereas the Planck data tend to favor lower values of $H_{0}$. Thus, as in the case of $\Lambda$CDM, the BH model does not alleviate the tension of $H_{0}$ between the Planck data and its local measurements. The similar property also holds for $\sigma_{8}^{(0)}$, where the Planck data favor higher values of $\sigma_{8}^{(0)}$ than those constrained in low-redshift measurements. We can also put further bounds on $\sigma_{8}^{(0)}$ by using the datasets of weak lensing measurements, such as KiDS Hildebrandt:2016iqg ; Kuijken:2015vca ; deJong:2015wca . For this purpose, we need to take non-linear effects into account in the MCMC analysis, which is beyond the scope of the current paper.

V.3 Model Selection

The BH model has three more parameters compared to those in $\Lambda$CDM. This means that the former has more freedom to fit the model better with the data. In order to study whether the former is statistically favored over the latter, we compute the Deviance Information Criterion (DIC) RSSB :

[TABLE]

where $\chi^{2}_{\rm eff}(\hat{\theta})=-2\ln\mathcal{L}(\hat{\theta})$, and $\hat{\theta}$ is a vector associated with model parameters maximizing the likelihood function $\mathcal{L}$. The quantity $p_{\rm D}$ is defined by $p_{\rm D}=\bar{\chi}^{2}_{\rm eff}({\theta})-\chi^{2}_{\rm eff}(\hat{\theta})$, where the bar represents an average over the posterior distribution. From its definition, the DIC accounts for the goodness of fit, $\chi^{2}_{\rm eff}(\hat{\theta})$, and the Bayesian complexity of the model, $p_{\rm D}$. The complex models with more free parameters give larger $p_{\rm D}$. To compare the BH model with the $\Lambda$CDM model, we calculate

[TABLE]

If $\Delta{\rm DIC}$ is negative, then BH is favored over $\Lambda$CDM. For positive $\Delta{\rm DIC}$, the situation is reversed.

In Table 4, we present the relative differences of $\Delta\chi_{\rm eff}^{2}$ and $\Delta{\rm DIC}$ in BH and GGC models, as compared to $\Lambda$CDM. Since $\Delta\chi_{\rm eff}^{2}$ are always negative, these models provide the better fit to the data relative to $\Lambda$CDM. In particular, we find that $\Delta\chi_{\rm eff}^{2}$ constrained by the Planck data alone are smaller than those derived with the PBRS datasets. This preference of BH over $\Lambda$CDM by the Planck data arises from combined effects of the suppressed large-scale ISW tale caused by the Galileon term and the modified high-$\ell$ TT power spectrum induced by the different background evolution relative to $\Lambda$CDM (as shown in Fig. 10). The former contributes by $\sim 20\,\%$ to a better $\chi^{2}_{\rm eff}$, while the latter to the remaining $\sim 80\,\%$. We note that a further lowering of the ISW tail is limited by the shift of acoustic peaks at high-$\ell$. Such modifications are also subject to further constraints from the datasets of BAO and SN Ia, but the values of $\Delta\chi_{\rm eff}^{2}$ constrained with the PBRS datasets are still negative in both BH and GGC models.

According to the DIC, the BH model is slightly disfavored over $\Lambda$CDM with the PBRS datasets. The GGC model, which has one parameter less than those in BH, is favored over $\Lambda$CDM with both Planck and PBRS datasets. This implies that the existence of an additional parameter $x_{4}$ does not contribute to provide better fits to the data. Indeed, today’s value of $x_{4}$ is severely constrained as Eq. (55) mostly from the CMB data. At the same time, this implies that there are no observational signatures for the deviation $\alpha_{\rm H}$ from Horndeski theories. It is interesting to note that the GGC model, which belongs to a sub-class of Horndeski theories, is statistically favored over $\Lambda$CDM even with two additional parameters, but this property does not persist in the BH model due to the extra beyond-Horndeski term $\alpha_{\rm H}$ modifying the cosmic expansion and growth histories.

VI Conclusion

We studied observational constraints on the BH model given by the action (1) with the functions (II). This model belongs to a sub-class of GLPV theories with the tensor propagation speed squared $c_{t}^{2}$ equivalent to 1. The deviation from Horndeski theories is weighed by the dimensionless parameter $\alpha_{\rm H}=4x_{4}/(5-x_{4})$, where $x_{4}$ is defined in Eq. (22). The BH model also has the $a_{2}X^{2}$ and $3a_{3}X\square\phi$ terms in the Lagrangian, which allow the possibility for approaching a de Sitter attractor from the region $-2<w_{\rm DE}<-1$ without reaching a tracker solution ($w_{\rm DE}=-2$).

Compared to the standard $\Lambda$CDM model, the beyond-Horndeski term $x_{4}$ can change the background cosmological dynamics in the early Universe. Since the Hubble expansion rate $H$ is modified by the non-vanishing $x_{4}$ term, this leads to the shift of acoustic peaks of CMB temperature anisotropies at high-$\ell$, see BH1 in Fig. 7. Moreover, as we observe in Fig. 4, the early-time dominance of $x_{4}$ over $x_{1,2,3}$ leads to the modified evolution of gravitational potentials $\Psi$ and $\Phi$ in comparison to $\Lambda$CDM, whose effect is more significant for small-scale perturbations. This modification also affects the evolution of radiation perturbations and the early-time ISW effect. As a result, the amplitude of CMB acoustic peaks is changed by the $x_{4}$ term. These modifications allow us to put bounds on the deviation from Horndeski theories.

The cubic Galileon existing in the BH model leads to the modified growth of matter perturbations and gravitational potentials at low redshifts. Provided that $x_{4}$ is subdominant to $x_{1,2,3}$, the dimensional quantities $\mu$ and $\Sigma$, which characterize the gravitational interactions with matter and light respectively, are given by Eq. (48) under the quasi-static approximation deep inside the sound horizon. Thus, the Galileon term $x_{3}$ enhances the linear growth of perturbations without the gravitational slip ($\mu\simeq\Sigma>1$). This enhancement can be seen in the lensing power spectrum $D_{\ell}^{\phi\phi}$ plotted in Fig. 5.

For the CMB temperature anisotropies, the late-time modified growth of perturbations caused by the cubic Galileon manifests itself in the large-scale ISW tale. The ISW effect is attributed to the variation of the lensing gravitational potential $\Psi+\Phi$ related to the quantity $\Sigma$. Unlike the $\Lambda$CDM model in which the time derivative $\dot{\Psi}+\dot{\Phi}$ is positive, the Galileon term $x_{3}$ allows the possibility for realizing $\dot{\Psi}+\dot{\Phi}$ closer to 0. In this case, the large-scale TT power spectrum is lower than that in $\Lambda$CDM, see GGC and BH2 in Fig. 7. Moreover, the modified background evolution at low redshifts induced by the Galileon leads to the shift of small-scale CMB acoustic peaks toward higher multipoles. If the contribution of $x_{3}$ to the total dark energy density is increased further, the ISW tale is subject to the significant enhancement compared to $\Lambda$CDM, together with the large shift of high-$\ell$ CMB acoustic peaks (see BH3 in Fig. 7). These large modifications to the TT power spectrum also arise for covariant Galileons without the $x_{2}$ term, whose behavior is disfavored from the CMB data Renk ; Peirone . In the BH model, the existence of $x_{2}$ besides $x_{3}$ can give rise to the moderately modified TT power spectrum being compatible with the data.

We put observational constraints on free parameters in the BH model by running the MCMC simulation with the datasets of CMB, BAO, SN Ia, and RSDs. With the Planck CMB data, we showed that today’s value of $x_{4}$ is constrained to be smaller than the order $10^{-6}$. Inclusion of other datasets does not modify the order of upper limit of $x_{4}^{(0)}$, and hence $|\alpha_{\rm H}^{(0)}|\leq{\cal O}(10^{-6})$. Apart from the bound arising from the GW decay to dark energy, this is the tightest bound on $|\alpha_{\rm H}^{(0)}|$ derived so far from cosmological observations.

The other dark energy density parameters $x_{1}^{(0)},x_{2}^{(0)},x_{3}^{(0)}$ are constrained to be in a similar way to those derived in Ref. PBFT . The best-fit value of $x_{3}^{(0)}$ is smaller than $|x_{1}^{(0)}|$ and $x_{2}^{(0)}$ by one order of magnitude. This intermediate value of $x_{3}^{(0)}$ leads to the CMB TT power spectrum with modifications at both large and small scales, in such a way that the BH model can be observationally favored over $\Lambda$CDM. The evolution of matter perturbations at low redshifts is not subject to the large modification by this intermediate value of $x_{3}^{(0)}$ in comparison to $\Lambda$CDM, so the BH model is also compatible with the RSD data. The best-fit background expansion history corresponds to the case in which $w_{\rm DE}$ finally approaches $-1$ from the phantom region $-2<w_{\rm DE}<-1$, whose behavior is consistent with the datasets of SN Ia and BAO. We also showed that, as in the $\Lambda$CDM model, the tensions in $H_{0}$ and $\sigma_{8}^{(0)}$ between CMB and low-redshift measurements are not alleviated for the datasets used in our analysis. Future investigations including non-linear effects and additional probes from weak lensing measurements will allow us to shed light on the possibility for alleviating such tensions in the BH model.

To make comparison between BH and $\Lambda$CDM models, we computed the DIC defined by Eq. (57) penalizing complex models with more free parameters. In BH, there are three additional parameters than those in $\Lambda$CDM. We found that the effective $\chi_{\rm eff}^{2}$ statistics in BH is smaller than that in $\Lambda$CDM for two combinations of datasets (Planck and PBRS). This is mostly due to both the suppressed ISW tail in BH and the shifts of high-$\ell$ acoustic peaks of the CMB TT power spectrum. These combined effects allow the BH model to fit the Planck data better. According to the DIC, however, there is a slight preference of $\Lambda$CDM over BH with both Planck and PBRS datasets. The beyond-Horndeski term $x_{4}$ generally works to prevent better fits to the data. The GGC model, which corresponds to $x_{4}=0$ with one parameter less than those in BH, is statistically favored over $\Lambda$CDM even with the DIC PBFT . This means that, at least in the BH model, there is no preference for the departure from Horndeski theories in cosmological observations.

We have thus shown that the deviation from Horndeski theories is severely constrained by the current observational data, especially from CMB. In spite of this restriction, the best-fit BH model gives the effective $\chi_{\rm eff}^{2}$ statistics smaller than that in $\Lambda$CDM. Moreover, the GGC model with $\alpha_{\rm H}=0$ leads to the smaller DIC relative to $\Lambda$CDM, even with two additional parameters. Thus, the BH and GGC models can be compelling and viable candidates for dark energy. Further investigations may be performed in several directions. In this work we considered massless neutrinos, but we plan to extend the analysis to include massive neutrinos and inquire about any degeneracy which can arise between such fluid components and modified gravitational interactions. Moreover, it is of interest to investigate cross-correlations between the ISW signal and galaxy distributions, which can be used to place further constraints on BH and GGC models.

Acknowledgements.

We thank Matteo Martinelli for support in the numerical implementation. We are grateful to N. Bartolo, A. De Felice, R. Kase, M. Liguori, S. Nakamura, M. Raveri and A. Silvestri for useful discussions and comments. SP acknowledge support from the NWO and the Dutch Ministry of Education, Culture and Science (OCW), and also from the D-ITP consortium, a program of the NWO that is funded by the OCW. GB acknowledges financial sup- port from Fondazione Ing. Aldo Gini. The research of NF is supported by Fundação para a Ciência e a Tecnologia (FCT) through national funds (UID/FIS/04434/2013), 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. NF, SP and GB acknowledge the COST Action (CANTATA/CA15117), supported by COST (European Cooperation in Science and Technology). ST is supported by the Grant-in-Aid for Scientific Research Fund of the JSPS No. 19K03854 and MEXT KAKENHI Grant-in-Aid for Scientific Research on Innovative Areas “Cosmic Acceleration” (No. 15H05890).