# CMB lensing bispectrum as a probe of modified gravity theories

Toshiya Namikawa, François R. Bouchet, Atsushi Taruya

► To cite this version:

Toshiya Namikawa, François R. Bouchet, Atsushi Taruya. CMB lensing bispectrum as a probe of modified gravity theories. Physical Review D, 2018, 98 (4), pp.043530. 10.1103/PhysRevD.98.043530. hal-01815161

HAL Id: hal-01815161

https://hal.science/hal-01815161v1

Submitted on 2 Oct 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# CMB lensing bispectrum as a probe of modified gravity theories

Toshiya Namikawa

Leung Center for Cosmology and Particle Astrophysics,
National Taiwan University, Taipei, 10617, Taiwan

François R. Bouchet

Institut d'Astrophysique de Paris, UMR7095, CNRS & Sorbonne Université, Paris, France
and Yukawa Institute for Theoretical Physics, Kyoto University, Kyoto 606-8502, Japan

Atsushi Taruya

Center for Gravitational Physics, Yukawa Institute for Theoretical Physics,
Kyoto University, Kyoto 606-8502, Japan

and Kavli Institute for the Physics and Mathematics of the Universe (WPI),
The University of Tokyo Institutes for Advanced Study, The University of Tokyo, 5-1-5 Kashiwanoha,
Kashiwa, Chiba 277-8583, Japan

<figure><img src="image_2.png" /></figure>

(Received 2 June 2018; published 28 August 2018)

Cosmological structures grow differently in theories of gravity which are modified as compared to Einstein's general relativity (GR). Cosmic microwave background (CMB) fluctuation patterns at the last scattering surface are lensed by these structures along the photon path to the observer. The observed CMB pattern therefore keeps trace of the growth history of structures. We show that observations of the CMB lensing bispectrum offer an interesting way to constrain deviations from GR in a broad class of scalar-tensor theories of gravity called “beyond Horndeski”. We quantify how the constraints on generic parameters describing the deviations from GR depend on the effective multipole range of the analysis. Our results further indicate that an accurate nonlinear correction of the matter bispectrum in the modified gravity considered is necessary when the bispectrum is used to probe scales beyond a multipole $\ell_{max} \gtrsim 1500$. We also found that the results are insensitive to details of the implementation of the screening mechanism, at very small scales. We finally demonstrate the potential of the lensing bispectrum to provide a blind reconstruction of the redshift evolution of our modified gravity parameters by combining the analysis of CMB and low-z source lensing data.

DOI: 10.1103/PhysRevD.98.043530

## I. INTRODUCTION

Dark energy is an enduring mystery. One possibility is that the observational evidence pointing to its existence, i.e., late-time cosmic acceleration, is misinterpreted and should rather be considered as traces of the breakdown of general relativity (GR) at large scales. Many experiments are being conducted in order to test this hypothesis, in particular by seeking a possible deviation from GR in the growth rate of large scale structures of the Universe (e.g., [1–5]).

One fundamental difficulty in devising tests of gravity on cosmological scales is that the large scale structures of the Universe are delineated by visible objects—galaxies, clusters of galaxies—which are biased tracers of the dark matter distribution, while this is the latter which is the dominant source of gravity at these scales. Of course, there are well developed models relating the light distribution to the mass distribution and a perturbative description of the galaxy bias should work well at large scales, only introducing a small number of nuisance parameters (e.g., [6]).

Nevertheless, deviations from GR are bound to be very small, and controlling the accuracy of these models at the required level of precision is a daunting challenge which will be hard to meet convincingly. It is therefore all too natural to try finding observables which rely on a solid understanding of their physics. Gravitational lensing effects are obvious candidates, which are very promising. However, practical analyses of the weak lensing probed with luminous objects, referred to as cosmic shear, still have to cope with a number of difficulties, e.g., the alignment of objects with their surrounding structures or uncertainties in the redshift distribution of the lensed sources. In addition, most of the lensing effect comes from redshifts around half that of the most distant sources, typically $z \sim 1$. One interesting avenue to explore will be the lensing effect on the Lyman-$\alpha$ forest on the line of sight of distant quasars (e.g., [7–9]). Here, we instead look at the lensing effect on the CMB, whose physical origin is very well known. In that case, the main contribution to the
lensing signal comes from a large redshift range centered at $z \sim 2$, and this allows us to probe a redshift range hardly attainable by the lensing effects probed with luminous sources.

Multiple works have already considered the future forecast to constrain modified gravity theories from the CMB lensing measurements, focusing in particular on the (angular) power spectrum and/or cross-power spectra with other cosmological observables(e.g., [10–14]). But these two-point statistics are far from being exhaustive, and further information can be extracted with the higher-order statistics. Refs. [15,16] showed that the lensing bispectrum of the CMB is detectable at high statistical significance in near-term experiments. It is thus natural, as a second step, to explore the potential power of measurements of the CMB lensing bispectrum in assessing modified gravity theories. The aim of this paper is to present the first forecast study of this probe.

This paper is organized as follows. Section II recalls the basic theory of the lensing bi-spectrum and how it is altered in modified gravity theories. Section III presents our main results, while Sec. IV summarizes our work and discusses perspectives.

Throughout the paper, we adopt a fiducial cosmological model which is consistent with the latest Planck cosmology [17], i.e., a spatially flat $\Lambda$CDM cosmology with the baryon and matter density, $\Omega_b h^2 = 0.0223$ and $\Omega_m h^2 = 0.119$, the dark-energy density $\Omega_\Lambda = 0.689$, the amplitude of the primordial scalar power spectrum, $A_s = 2.13 \times 10^{-9}$, and its spectral index at $k = 0.05$ Mpc$^{-1}$, $n_s = 0.965$, with the reionization optical depth, $\tau = 0.0630$.

## II. LENSING BISPECTRUM

In this section, we first review briefly the formalism of the CMB lensing bispectrum in GR (see, e.g., [15,18]). We then describe the modification of the bispectrum in a general class of modified gravity theories (see, e.g., [19–21]).

### A. Lensing potential

The gravitational lensing effect on the CMB anisotropies is described as a remapping of the CMB fluctuations at the recombination by the so-called deflection angle, $d = \nabla\phi$,<fnref n="1" /> where $\phi$ is the lensing potential (e.g., [24,25]):

$$ \phi(\hat{n}) = -2 \int_{0}^{\chi_{*}} d\chi W(\chi, \chi_{*}) \Psi(\chi, \hat{n}). \quad (1) $$

Here $\chi_*$ is the comoving distance to the CMB last-scattering surface, and $\Psi$ is the Weyl potential. The lensing kernel, $W(\chi, \chi_*)$, is defined (for a flat cosmology) as

$$ W(\chi, \chi_*) = \frac{\chi_* - \chi}{\chi_* \chi_*} \quad (2) $$

Let us define the lensing convergence as $\kappa = -\nabla^2\phi/2$. By using Poisson equation, we obtain

$$ \kappa(\hat{n}) = \int_{0}^{\chi_{*}} d\chi \frac{3\Omega_{m,0}H_{0}^{2}\chi^{2}}{2a(\chi)} W(\chi, \chi_{*})\delta_{m}(\chi, \hat{n}), \quad (3) $$

where $a$ is the scale factor and $\delta_m$ is the underlying density fluctuations of matter along the line-of-sight. The lensing potential (and thus the convergence) can be reconstructed from the observed CMB anisotropies by using the fact that a fixed lensing potential introduces statistical anisotropy into the observed CMB (e.g., [26–28]). The reconstructed lensing convergence map can then be used for cosmology by first transforming the map into its harmonic coefficients, $\kappa_{\ell m}$, and then by measuring various moments like the power spectrum and the bispectrum (which entails an accurate subtraction of the nonlensing contributions, see e.g., [29,30]). Recent CMB experiments have already detected the power spectrum of the lensing potential very precisely [31–35]. For instance, the Planck detection [33] has a $40\sigma$ significance. The detection and precise determination of the CMB lensing bispectrum are therefore obvious and important next step in CMB scientific analyses.

### B. Lensing bispectrum

The bispectrum of the lensing convergence, defined in harmonic space, has translational and rotational invariance. This is true as long as the statistical isotropy holds. It is thus sufficient to characterize it with a function $B_{\ell_1\ell_2\ell_3}$ of only three variables, weighted by the Wigner $3j$-symbols through

$$ \langle \kappa_{\ell_1 m_1} \kappa_{\ell_2 m_2} \kappa_{\ell_3 m_3} \rangle = \begin{pmatrix} \ell_1 & \ell_2 & \ell_3 \\ m_1 & m_2 & m_3 \end{pmatrix} B_{\ell_1 \ell_2 \ell_3}. \quad (4) $$

In what follows, we use the flat-sky approximation, and start by computing the bispectrum given by

$$ \langle \kappa_{\ell_1 \ell_2 \ell_3} \rangle = (2\pi)^2 \delta(\ell_1 + \ell_2 + \ell_3) B_{\ell_1 \ell_2 \ell_3}. \quad (5) $$

The full-sky bispectrum is then obtained from the flat-sky bispectrum through

$$ B_{\ell_1 \ell_2 \ell_3} = \begin{pmatrix} \ell_1 & \ell_2 & \ell_3 \\ 0 & 0 & 0 \end{pmatrix} \sqrt{\frac{(2\ell_1 + 1)(2\ell_2 + 1)(2\ell_3 + 1)}{4\pi}} \times B_{\ell_1 \ell_2 \ell_3}, \quad (6) $$

where the multipoles have to satisfy the triangle condition $|\ell_i - \ell_j| \le \ell_k \le \ell_i + \ell_j$. In order to evaluate this expression, it is enough at the scales of interest to use the following approximate form of the Wigner $3j$-symbol,

<footnote n="1">¹Here we ignore the curl mode of the deflection angle (e.g., [22]) because it contribution is negligible compared to that of the lensing potential in the standard $\Lambda$CDM cosmology (e.g., [23]).</footnote>
$$ \begin{pmatrix} \ell_1 & \ell_2 & \ell_3 \\ 0 & 0 & 0 \end{pmatrix} \approx (-1)^L \sqrt{\frac{e^3}{2\pi}} (L+1)^{-1/4} \prod_{i=1}^{3} (L-\ell_i+1)^{-1/4} \times \left( \frac{L-\ell_i+1/2}{L-\ell_i+1} \right)^{L-\ell_i+1/4} \quad (7) $$

for even $L$, where we define $L = (\ell_1 + \ell_2 + \ell_3)/2$. For odd $L$, the Wigner 3j-symbol becomes zero.

The CMB lensing bispectrum is sourced by

(1) the nonlinear evolution of the large-scale structure [15], and

(2) the so-called post-Born correction [16], i.e., the correction to the Born approximation for which the lensing effect is evaluated on the unperturbed geodesic.

We denote these contributions respectively by $B^{\text{LSS}}$ and $B^{\text{pb}}$, and present their explicit expressions below.

The CMB lensing bispectrum from the nonlinear growth of the density perturbations is given in the flat-sky limit by [15]

$$ B^{\text{LSS}}(\ell_1, \ell_2, \ell_3) = \int_0^{\chi_*} d\chi \left[ \frac{3\Omega_{m,0}H_0^2}{2a(\chi)} \right]^3 \times \chi^2 W^3(\chi, \chi_*) B_{\delta}(\mathbf{k}_1, \mathbf{k}_2, \mathbf{k}_3, \chi). \quad (8) $$

Here, $B_{\delta}$ is the matter bispectrum arising from the nonlinear growth of structure. In the weakly nonlinear regime, it can be obtained by using perturbation theory. The result at the tree-level order is of the general form

$$ B_{\delta}(\mathbf{k}_1, \mathbf{k}_2, \mathbf{k}_3, \chi) = 2F_2(\mathbf{k}_1, \mathbf{k}_2, z)P_{\delta}(k_1, z)P_{\delta}(k_2, z) + 2 \text{ perms}, \quad (9) $$

where $P_{\delta}(k, z)$ is the matter power spectrum at redshift $z(\chi)$, and the function $F_2$ is the second-order perturbation theory kernel (e.g., [36]). Writing $\mathbf{k}_1 \cdot \mathbf{k}_2 = k_1 k_2 \cos \theta$, it is given by

$$ F_2(\mathbf{k}_1, \mathbf{k}_2, z) = \frac{5}{7} a(k_1, z) a(k_2, z) + \frac{1}{2} \frac{k_1^2 + k_2^2}{k_1 k_2} b(k_1, z) b(k_2, z) \cos \theta + \frac{2}{7} c(k_1, z) c(k_2, z) \cos^2 \theta, \quad (10) $$

where $a(k, z)$, $b(k, z)$ and $c(k, z)$ are unity at the tree-level of perturbation theory. In the highly nonlinear regime, the deviation from tree-level prediction are significant, and a proper treatment of the nonlinear effects coming from the higher-order perturbations is needed. The scale- and time-dependent coefficients $a, b$,

and $c$ effectively characterize these, and their deviation from unity is calibrated with high-resolution N-body simulations. According to [37], they are given by

$$ a(k, z) = \frac{1 + \{\sigma_8(z)\}^{a_6} \sqrt{0.7Q(n_{\text{eff}})} (qa_1)^{n_{\text{eff}}+a_2}}{1 + (qa_1)^{n_{\text{eff}}+a_2}} \quad (11) $$

$$ b(k, z) = \frac{1 + 0.2a_3(n_{\text{eff}} + 3)(qa_7)^{n_{\text{eff}}+3+a_8}}{1 + (qa_7)^{n_{\text{eff}}+3.5+a_8}} \quad (12) $$

$$ c(k, z) = \frac{1 + [4.5a_4/(1.5 + (n_{\text{eff}} + 3)^4)](qa_5)^{n_{\text{eff}}+3+a_9}}{1 + (qa_5)^{n_{\text{eff}}+3.5+a_9}}, \quad (13) $$

with $Q(x) = (4 - 2^x)/(1 + 2^{x+1})$. Here, the variable $q$ is given by $q = k/k_{\text{NL}}$ with the nonlinear scale, $k_{\text{NL}}$, determined by $4\pi k_{\text{NL}}^3 P_{\text{m}}^{\text{lin}}(k_{\text{NL}}) = 1$. The quantity $\sigma_8(z)$ is the variance of the matter density fluctuations smoothed with a top-hat sphere of radius $8h^{-1}$ Mpc at redshift $z$. The logarithmic slope, $n_{\text{eff}} \equiv d \ln P_{\text{m}}^{\text{lin}}(k)/d \ln k$, is the effective spectral index of the linear power spectrum, $P_{\text{m}}^{\text{lin}}(k)$. The parameters, $a_i$, are determined by fitting results of N-body simulations, which yields [37]

$$ \begin{aligned} a_1 &= 0.250 & a_2 &= 3.50 & a_3 &= 2.00 \\ a_4 &= 1.00 & a_5 &= 2.00 & a_6 &= -0.200 \\ a_7 &= 1.00 & a_8 &= 0.00 & a_9 &= 0.00. \end{aligned} \quad (14) $$

Later on, Ref. [38] proposed an improved fit given by

$$ \begin{aligned} a_1 &= 0.484 & a_2 &= 3.74 & a_3 &= -0.849 \\ a_4 &= 0.392 & a_5 &= 1.01 & a_6 &= -0.575 \\ a_7 &= 0.128 & a_8 &= -0.722 & a_9 &= -0.926. \end{aligned} \quad (15) $$

In our baseline calculations, we use the parameters of Ref. [38] (hereafter “GM”) but we shall also use the earlier results of Ref. [37] (hereafter “SC”) as a means to assess the dependence of our results on the accuracy of these fitting formula.

The post-Born correction to the CMB lensing bispectrum, $B^{\text{pb}}$ is given as [16]

$$ B^{\text{pb}}(\ell_1, \ell_2, \ell_3) = 2 \frac{\ell_1 \cdot \ell_2}{\ell_1^2 \ell_2^2} \times [\ell_1 \cdot \ell_3 M(\ell_1, \ell_2) + \ell_2 \cdot \ell_3 M(\ell_2, \ell_1)] + \text{cyc perm}, \quad (16) $$

where
$$M(\ell, \ell') = \ell^4 \int_0^{\chi_*} d\chi \frac{[W(\chi, \chi_*)]^2}{\chi^2} \\ \times P_{\psi}\left(\frac{\ell}{\chi}, \chi\right) C_{\ell'}^{\kappa\kappa}(\chi, \chi_*). \quad (17)$$

$P_{\psi}(k, \chi)$ is the power spectrum of the Weyl potential at a comoving distance $\chi$, and

$$C_{\ell'}^{\kappa\kappa}(\chi', \chi_*) = \ell^4 \int_0^{\chi'} d\chi \frac{W(\chi, \chi')W(\chi, \chi_*)}{\chi^2} \\ \times P_{\psi}\left(\frac{\ell}{\chi}, \chi\right). \quad (18)$$

If these terms are known, or determined accurately, then the post-Born term is known as well.

All in all, the CMB lensing bispectrum is the sum of Eqs. (8) and (16).

### C. Effect of modified gravity on the bispectrum

In modified gravity theories, the perturbation theory kernel $F_2$ is altered and may be written, at the tree level in the quasistatic approximation, as (e.g., [20,21])

$$F_2(\mathbf{k}_1, \mathbf{k}_2, z) = \left( \kappa(z) - \frac{2}{7}\lambda(z) \right) \\ + \kappa(z) \frac{1}{2} \frac{k_1^2 + k_2^2}{k_1 k_2} \cos\theta \\ + \lambda(z) \frac{2}{7} \cos^2\theta. \quad (19)$$

The kernel above coincides with the GR case when $\lambda(z) = \kappa(z) = 1$. Let us note that here and in all the following we use $\kappa$ as a parameter characterizing the deviation from GR, as in many previous theoretical papers. This is then not to be confused with the lensing convergence which is denoted by the same symbol in other parts of the literature [and in (3)].

The Horndeski theory of gravity (e.g., Refs. [39–41]), is the most general nondegenerate scalar-tensor theory in 4D space-time that leads to second-order equations of motion. It may have $\lambda \neq 1$ in general, but $\kappa = 1$ is still preserved. An even wider class of theories imaginatively called “beyond Horndeski” theories, including GLPV [42,43] and DHOST [44,45], can explicitly break this latter condition, in close connection with the violation of Vainshtein mechanism to recover GR at nonlinear regime (see also [46]). Testing and constraining possible deviations of $\lambda$ and $\kappa$ from unity is thus very interesting, and could give important information on gravity at cosmological scales, rather independently of the growth rate of structure probed with galaxy redshift surveys. Further, no strong constraint on $\lambda$ has been obtained so far, and no theoretical upper/lower limits is known for $\kappa$ (see, e.g., Refs. [19–21] for

further discussion of the possible values of $\lambda$ and $\kappa$ under some specific models). The measurement of the bispectrum is therefore key to narrow down the constraints.

A particular subclass of modified gravity theories may have a specific redshift dependence of $\lambda$ and $\kappa$, and such a form will have to be used ultimately to get the tightest constraints on these specific theories. Here, we rather wish to look at the generic potential of the bispectrum probe. Following Ref. [20], we adopt the functional form of $\lambda$ and $\kappa$ as

$$\begin{aligned} \lambda(z) &= [\Omega_m(z)]^{\xi_\lambda}, \\ \kappa(z) &= [\Omega_m(z)]^{\xi_\kappa}, \end{aligned} \quad (20)$$

where $\Omega_m(z) = \Omega_{m,0}/(\Omega_{m,0} + a^3\Omega_\Lambda)$. This naturally embodies the expectation that a modified theory will converge to GR at high-$z$ and preserve the successful predictions of CMB anisotropies using GR. The form is also monotonic, in keeping with the idea of generic constraints. Our goal is therefore to assess how well the generic parameters $\xi_\lambda$ and $\xi_\kappa$ can be constrained by measurement of the CMB lensing bispectrum.

One should recall that the modification of the kernel $F_2$ given in Eq. (19) is only valid in the weakly nonlinear regime. In order to improve the constraints, we may want to use measurements at small scales, taking into account the nonlinear corrections introduced earlier in Sec. II B. This is, however, not trivial in the context of modified gravity, because a modification of gravity can also change the nonlinear corrections, and a proper account of these needs more elaborate work which we leave for future investigation. Here we rather adopt the fitting formula given in the GR case, in order to assess the potential power of the CMB lensing bispectrum in the intermediate regime. Therefore the kernel $F_2$ used in our analysis is given by

$$F_2(\mathbf{k}_1, \mathbf{k}_2, z) = \left( \kappa(z) - \frac{2}{7}\lambda(z) \right) a(k_1, z) a(k_2, z) \\ + \kappa(z) \frac{1}{2} \frac{k_1^2 + k_2^2}{k_1 k_2} b(k_1, z) b(k_2, z) \cos\theta \\ + \lambda(z) \frac{2}{7} c(k_1, z) c(k_2, z) \cos^2\theta. \quad (21)$$

We will then compare the forecast results based on Eq. (21) with those derived from the tree-level kernel at Eq. (19). The impact of the modification of this formula at small scales is also discussed in detail (see Sec. III C).

Let us conclude this theoretical section with a couple of comments. The first concerns the relation between the lensing potential and the density field, which was given in Eq. (3). This equation may be altered in modified gravity theories, and the expression given at Eq. (8) might not be relevant. Fortunately, Ref. [14] showed that the modified Weyl potential is given by a simple scaling of the matter
<figure><img src="image_3.png" /><figcaption>FIG. 1. The CMB lensing bispectrum in GR and modified theories of gravity ($\xi_{\lambda} = 1$ or $\xi_{\kappa} = 1$) with error bars expected from a future CMB experiment (S4) and in the cosmic-variance limit (CV), see specifications in Table I. The left and right panel show the case for an equilateral and folded configuration with, respectively, $\ell_1 = \ell_2 = \ell_3 \equiv \ell$ and $\ell_1 = 2\ell_2 = 2\ell_3 \equiv \ell$.</figcaption></figure>

density fluctuations by a factor $\mu$. This means that what is actually constrained from the lensing bispectrum are combinations of the parameters, i.e., $\mu^3\lambda$ and $\mu^3\kappa$. This being understood, in the following, we keep using the same notation, and denote these scaled parameters simply by $\lambda$ and $\kappa$.

The second point concerns the post-Born correction term, which could also receive corrections from the modification of gravity. However, the impact of such a modification is expected to be very small because the measurement of the lensing power spectrum severely limits a modification to the post-Born correction in the bispectrum. Indeed, the signal-to-noise ratio of the lensing power spectrum will be very high in future experiments [~$O(10^3)$] and the allowed modification to the power spectrum amplitude is therefore smaller than about 0.1%. The signal-to-noise ratio of the bispectrum is, on the other hand, much lower than that of the power spectrum [~$O(10)$], which is equivalent to a ~10% constraint on the modification of the bispectrum amplitude. This means that the prior information from the power spectrum limits the modification to the post-Born correction well below the measurement uncertainty of the bi-spectrum. For this reason, in our analysis, we ignore the effect of possible modification to the post-Born correction.

The third point is the validity of the nonlinear fitting formula beyond tree level.

### III. RESULTS

#### A. Expected constraints on parameters of modified gravity theories

We begin by inspecting the expected signal for lensing bi-spectrum for various experiments and configurations. In

Fig. 1, the bi-spectra for equilateral ($\ell_1 = \ell_2 = \ell_3 \equiv \ell$, left) and folded ($\ell_1 = 2\ell_2 = 2\ell_3 \equiv \ell$, right) configurations are shown, and the results are plotted as a function of $\ell$.<fnref n="2" /> Solid lines are the expected bispectrum signals in GR (black) and modified gravity models ($\xi_{\lambda} = 1$ for red, $\xi_{\kappa} = 1$ for blue), which are computed following the prescription in Sec. II C. The modification to gravity alters the scale-dependence of the bispectrum amplitude through the shape dependence of the kernel $F_2$, but only weakly, and the main effect appears to be a simple rescaling. In the folded bispectrum, however, the post-Born correction has non-negligible contributions to the bispectrum, and the total bispectrum has a bit complicated behavior. In the figure, the statistical errors are depicted by shaded areas, both for the upcoming experiment, CMB-S4 (green), and an hypothetical full-sky experiment limited only by the cosmic variance (red). Note that the quoted error is estimated from the Gaussian noise used in the Fisher matrix (see below).

Fig. 1 implies that with a realistic measurement of the lensing bi-spectrum, one can simultaneously constrain both modified gravity parameters. We note that these configurations at large scales ($\ell \lesssim 1000$) are a simple rescaling of the GR case, and as such might be hard to disentangle from other observational effects, a further reason to hope for small scale measurements. In addition, the folded

<footnote n="2">²The error bars are computed as follows. We first choose the width of the multipole bin and bin centers of the multipole, $\ell_1$, $\ell_2$ and $\ell_3$. Using the three multipole bin centers, the binned equilateral and folded bi-spectra satisfy $\ell_1 = \ell_2 = \ell_3$ and $\ell_1 = 2\ell_2 = 2\ell_3$, respectively. For each binned bi-spectrum, $B_{\ell_1, \ell_2, \ell_3}$, its signal-to-noise (i.e., the inverse of the error on its amplitude) is computed. The error bars plotted in the figure are obtained by multiplying $B_{\ell_1, \ell_2, \ell_3}$ to the errors of the bispectrum amplitude at each bin.</footnote>
<table><caption>TABLE I. Specifications for the CMB experiments considered in this paper: the noise level in the CMB polarization map, $\Delta_P$, in unit of $\mu$K arcmin, the angular resolution as expressed by the FWHM of a Gaussian beam, $\theta$, in unit of arcmin, and the fractional sky coverage, $f_{sky}$.</caption><thead><tr><th></th><th>Δ<sub>P</sub> [μK arcmin]</th><th>θ [arcmin]</th><th>f<sub>sky</sub></th></tr></thead><tbody><tr><td>CMB-S4</td><td>1</td><td>3</td><td>0.4</td></tr><tr><td>CV</td><td>0</td><td>...</td><td>1.0</td></tr></tbody></table>

bi-spectrum is much more sensitive to the modified gravity parameters than the equilateral bispectrum because the amplitude of the LSS bi-spectrum in the folded case is much larger than that in the equilateral case [15].

In order to see quantitatively how well the parameters $\xi_\lambda$ and $\xi_\kappa$ can be constrained, we follow Refs. [15,18] and define the Fisher matrix:

$$F_{ij} = \sum_{\ell_1 \le \ell_2 \le \ell_3}^{\ell_{\max}} f_{\text{sky}} \frac{B_{\ell_1 \ell_2 \ell_3, i} B_{\ell_1 \ell_2 \ell_3, j}}{\Delta_{\ell_1 \ell_2 \ell_3} C_{\ell_1} C_{\ell_2} C_{\ell_3}}, \quad (22)$$

where $B_{\ell_1 \ell_2 \ell_3, i}$ is the derivative of the lensing bispectrum with respect to the $i$th parameter. $\Delta_{\ell_1 \ell_2 \ell_3}$ is unity if all $\ell_i$ are different, 2 if two $\ell_i$ are equal, and 6 if all $\ell_i$ are equal. The lensing power spectra, $C_\ell$, includes the reconstruction noise of the lensing measurement.

The lensing reconstruction noise is computed by following the formula of Ref. [47] which is motivated by the

<figure><img src="image_4.png" /><figcaption>FIG. 2. Left: Expected constraints on $\xi_\lambda$ as a function of $\ell_{\max}$ and experimental specifications, with and without the nonlinear correction (see text). The constraints are derived by marginalizing only $\xi_\lambda$ as a free parameter. Right: Same as left but for $\xi_\kappa$. Note that the curves labeled “tree” correspond to Eq. (19) and those labeled “GM” to Eq. (21).</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>FIG. 3. The expected joint constraints on the two parameters, $\xi_\lambda$ and $\xi_\kappa$, assuming the CMB-S4 experiment (left) and in the cosmic-variance case (right). Ellipses of different colors correspond to changes in the maximum multipole of the bispectrum used to constrain the parameters ($\ell_{\max}$ = 500, 1000, 1500). We also show with dashes the case when including a correction of the nonlinear evolution of the large-scale structure in the bispectrum (denoted by “GM”).</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>FIG. 4. Expected constraints on $\xi_{\lambda}$ (left) and $\xi_{\kappa}$ (right) as a function of $\ell_{\text{max}}$ with and without the fitting model of the nonlinear correction and experimental specifications. Compared to Fig. 2, we assume that the two parameters are simultaneously constrained.</figcaption></figure>

maximum-likelihood approach to reconstructing the lensing potential [27,48]. The specifications of the CMB instrumental noise we use are summarized in Table I.

The constraining power at various scales is investigated by varying the maximum multipole $\ell_{\text{max}}$ of the summation of Eq. (22). Derivatives are numerically computed based on the symmetric difference quotient. The fiducial values of the parameters are GR ones, i.e., $(\xi_{\lambda}, \xi_{\kappa}) = (0, 0)$. The other cosmological parameters are fixed in the analysis because these parameters are severely constrained by other observables such as the primary CMB power spectrum and lensing potential power spectrum. The linear matter power spectrum is computed with CAMB [49]. We use the nonlinear correction to the matter power spectrum of Refs. [50,51].

Figure 2 shows the expected size of the statistical error on the parameters $\xi_{\lambda}$ (left) and $\xi_{\kappa}$ (right) as a function of the maximum multipole considered, $\ell_{\text{max}}$, when only one parameter at a time, either $\xi_{\lambda}$ or $\xi_{\kappa}$, is free to vary while keeping the other fixed. These plots show quantitatively how fast the constraints improve when smaller scales are included. Let us also note that the $1\sigma$ constraints do not vary simply as a power law of $\ell_{\text{max}}$.<fnref n="3" /> One also sees that the constraint on $\xi_{\kappa}$ is much better than that on $\xi_{\lambda}$. This is mostly due to the factor 2/7 in front of $\lambda$; in addition, the terms proportional to $\lambda$ vanish for the folded configuration ($\ell_1 = 2\ell_2 = 2\ell_3$ and its permutations) which has dominant contributions to the large-scale structure bispectrum [15]. In any case, future CMB experiments will allow

exploring a relevant part of the parameter space which is so far nearly unconstrained.

Figure 3 shows the error contours of the joint constraint on the parameters $\xi_{\lambda}$ and $\xi_{\kappa}$. The degeneracy between the two parameters is clearly seen. It does not change with $\ell_{\text{max}}$, nor is it really broken by the nonlinear corrections. This degeneracy comes from the first term of the kernel $F_2$ in Eq. (21).

Figure 4 shows the dependence of the statistical error of $\xi_{\lambda}$ (left) and $\xi_{\kappa}$ (right) on $\ell_{\text{max}}$, when the two parameters are simultaneously constrained. By comparing Fig. 4 with Fig. 2, due to the parameter degeneracy, the $1\sigma$ uncertainty in each constrained parameter becomes much larger compared to that in the absence of the degeneracy.

It is worth stressing that the correction due to the nonlinear growth beyond tree level significantly increases the total signal-to-noise of the bispectrum [15]. Correspondingly, the constraints on the parameters become significantly tighter.

## B. Impact of nonlinear loop correction on the lensing bispectrum

Our results imply that the study of the nonlinear growth in modified gravity theories is very important to best extract information in the lensing bispectrum in the era of CMB-S4 and beyond. We have used the fitting formula of GM for the matter bispectrum beyond the tree-level prediction [see Eq. (15)]. Unfortunately, its validity for deriving constraints on the modified gravity is not yet well studied. However, Ref. [52] recently uses the same parametrization of the fitting formula as we proposed in Eq. (15), and shows that the formula in the nDGP model, which falls into the Horndeski theories, is in good agreement with the simulation results. Although the fitting formula gives a reasonable approximation in at least the nDGP model as one of representative Horndeski theories, a complete treatment of nonlinear correction in modified gravity theories is still

<footnote n="3">³In the case of the power spectrum, $C_{\ell}$, if the derivative of $C_{\ell}$ with respect to a parameter, $p$, is proportional to $C_{\ell}$, the Fisher matrix for $p$ in the cosmic variance limit is simply given by $F_{pp} = \sum_{1 \le l_1 \le l_2 \le l_3} (l_1 + 1/2) \approx \ell_{\text{max}}^2 (l_{\text{max}} \gg 1)$. Then the $1\sigma$ constraint on $p$ varies as $1/\ell_{\text{max}}$. Nevertheless, Eq. (22) is not simply written as a function of $\ell_{\text{max}}$ and the $\ell_{\text{max}}$ dependence of $\sigma$ becomes more complicated.</footnote>
<figure><img src="image_7.png" /><figcaption>FIG. 5. Expected constraints using the nonlinear loop correction of Ref. [38] compared with that of Ref. [37].</figcaption></figure>

highly involved and well beyond the scope of this paper. Here, we get a feel of the sensitivity of our results to the specifics of the nonlinear correction by contrasting those obtained with GM with those derived with the SC fitting parameters [see Eq. (14)].

The expected $1\sigma$ error contour in the two dimensional parameter space using the SC fitting parameters is shown in Fig. 5. The expected constraints on $\xi_\lambda$ and $\xi_\kappa$ do not significantly deviate from those obtained in the GM case. This suggests than our results are actually quite robust against the details of the nonlinear correction.

Still, in addition to the effect on the uncertainties of the parameter constrains, the inaccuracy in the fitting formula could lead to a bias in the best fit parameters. Assuming that the GM fitting formula correctly captures the matter bispectrum beyond the tree level, we evaluate the expected bias in the parameter estimations by using on purpose an “inaccurate” fitting formula as follows. Specifically, we regard the SC fitting formula as the inaccurate model, and compute the expected parameter bias as [53,54]

$$b_i = \sum_j \tilde{F}_{ij}^{-1} \Delta_j, \qquad (23)$$

<table><caption>TABLE II. Expected bias in the parameters, $\xi_\lambda$ and $\xi_\kappa$, if the nonlinear loop correction of Ref. [29] is inaccurate while that of Ref. [30] correctly captures the loop correction. We also show between parentheses the expected $1\sigma$ constraints with the GM fitting parameters.</caption><thead><tr><th></th><th>ℓ<sub>max</sub></th><th>Δξ<sub>λ</sub> (σ(ξ<sub>λ</sub>))</th><th>Δξ<sub>κ</sub> (σ(ξ<sub>κ</sub>))</th></tr></thead><tbody><tr><td rowspan="3">CMB-S4</td><td>500</td><td>0.26 (2.8)</td><td>0.048 (0.57)</td></tr><tr><td>1000</td><td>0.45 (1.3)</td><td>0.053 (0.28)</td></tr><tr><td>1500</td><td>0.49 (0.92)</td><td>0.033 (0.20)</td></tr><tr><td rowspan="3">CV</td><td>500</td><td>0.26 (1.6)</td><td>0.048 (0.34)</td></tr><tr><td>1000</td><td>0.48 (0.65)</td><td>0.053 (0.14)</td></tr><tr><td>1500</td><td>0.53 (0.35)</td><td>0.026 (0.079)</td></tr></tbody></table>

where we define

$$\Delta_j = \sum_{\ell_1 \le \ell_2 \le \ell_3} f_{\text{sky}} \frac{(B_{\ell_1 \ell_2 \ell_3} - \tilde{B}_{\ell_1 \ell_2 \ell_3}) \tilde{B}_{\ell_1 \ell_2 \ell_3, j}}{\Delta_{\ell_1 \ell_2 \ell_3} C_{\ell_1} C_{\ell_2} C_{\ell_3}}, \qquad (24)$$

The quantities, $\tilde{F}_{ij}$ and $\tilde{B}_{\ell_1 \ell_2 \ell_3}$, are the Fisher matrix and lensing bispectrum computed with the incorrect theoretical model of the fitting formula, respectively.

Table II shows the results of the parameter bias for several experimental specifications and maximum multipole. As the maximum multipole of the bispectrum increases, the bias also increases. In the CMB-S4 case, the bias is still well within the $1\sigma$ expected constraints. It will thus not require heavy theoretical work to solidify the constraints on modified gravity theories at loop level. On the other hand, the bias for $\ell_{\text{max}} \gtrsim 1500$ in the CV limit is actually fairly significant compared to the $1\sigma$ constraint.

### C. Impact of the screening mechanism on the lensing bispectrum

At small scales, modified theories of gravity have to become close to GR in order to evade the observational constraints based on solar-system measurements. This is insured by a screening mechanism. Following the treatment of Ref. [55], we model this screening mechanism through the following modification of the parameters in $F_2(k_m, k_n, z)$,

$$\lambda(z) \rightarrow 1 + (\lambda(z) - 1) f(k_m, k_n, z), \qquad (25)$$

$$\kappa(z) \rightarrow 1 + (\kappa(z) - 1) f(k_m, k_n, z), \qquad (26)$$

where $(m, n)$ is (1,2), (2,3) or (1,3) and

$$f(k_m, k_n, z) = \exp \left[ - \frac{k_m^2 + k_n^2}{k_V^2(z)} \right]. \qquad (27)$$
<figure><img src="image_8.png" /><figcaption>FIG. 6. The expected constraints as a function of the screening scale, $v = k_v/k_{NL}$, divided by those without the screening mechanism.</figcaption></figure>

Here, we further introduce the parameter, $v$, defined through $k_V(z) = vk_{NL}(z)$. The value of $v$ would be typically $v \gg 1$, but could be $v \le 1$. Figure 6 shows the impact of the screening mechanism, i.e., the constraints on $\xi_\lambda$ (left) and $\xi_\kappa$ (right) divided by those without the screening mechanism, plotted as a function of $v$. For $v \ge 10$, the impact of screening mechanism on the parameter constraints is less significant, and it slightly worsens the constraints by 1% (10%) at $\ell_{max} = 500(1500)$. On the other hand, if $v \lesssim 1$, the change is significant, and the statistical errors get increased approximately by a factor of more than 2. Our result indicates that to constrain the models with $v \lesssim 1$, a proper account of the screening effect is very crucial. But, this may in turn suggest a possibility of the simultaneous constraints on the scale of the screening mechanism and parameters $\xi_{\lambda,\kappa}$.

## D. Comparison with low-z lensing bispectrum

So far we have focused on the constraints using the CMB lensing bispectrum as a clean cosmological probe. Indeed the lensing bispectrum directly measures the underlying gravitational potential and, unlike the galaxy bispectrum discussed by many previous works (e.g., [20,46]), it is immune to many observational and theoretical uncertainties, e.g., galaxy biases. The bispectrum of the cosmic shear is an alternative observational probe of modified gravity theories (e.g., [56]) which is based on the same physics. The resulting constraints on the modified gravity theories are expected to be tighter than those obtained from the CMB lensing bispectrum (at least for monotonous deviations). However, the cosmic shear has significant uncertainties in, e.g., the theoretical modeling of the nonlinear growth of the large-scale structure, baryon physics, intrinsic alignment, and observational difficulties such as photo-z, PFS, and calibration.

Table III shows the expected bias and $1\sigma$ constraints on the parameters, $\xi_\lambda$ and $\xi_\kappa$, using the cosmic shear at $z_s = 1$

and 2, to be compared with the CMB results of Table II. Compared to the latter case, the bias is larger than the $1\sigma$ constraint even at $\ell_{max} = 500$. The lensing bi-spectrum of the low-z galaxies is highly sensitive to the modeling of the nonlinear correction, and the constraints derived from the cosmic shear bispectrum would be easily biased by limitation of the accuracy of the nonlinear correction. In this respect, the constraints derived from the CMB lensing bi-spectrum will serve as an important cross-check of the results obtained from the cosmic shear bispectrum. We shall now see that CMB and shear lensing constraints are actually fairly complementary.

Here we discuss which redshifts the lensing bispectra of CMB and galaxies are sensitive to. To that effect, we consider a time dependence of the parameter, $\lambda(z)$, varying according to a top-hat function selecting a redshift bin:

$$ \lambda(z) = \begin{cases} 1 + \epsilon & (z_i \le z \le z_i + 0.1) \\ 1 & (\text{otherwise}) \end{cases}. \quad (28) $$

We compute the expected constraints on $\epsilon$ with varying $z_i$. Note that the results of the Fisher analysis do not depend on the fiducial value of $\epsilon$. We also change the source comoving distance, $\chi_*$, to include cases of galaxy weak lensing.

Figure 7 shows the expected constraints on the parameter, $\lambda(z)$, at each redshift bin. We vary the maximum

<table><caption>TABLE III. Same as Table II but for the lensing bispectrum probed with low-z sources, in the cosmic-variance limit.</caption><thead><tr><th></th><th>ℓ<sub>max</sub></th><th>Δξ<sub>λ</sub> (σ(ξ<sub>λ</sub>))</th><th>Δξ<sub>κ</sub> (σ(ξ<sub>κ</sub>))</th></tr></thead><tbody><tr><td rowspan="3">z<sub>s</sub> = 1</td><td>500</td><td>0.068 (0.11)</td><td>0.027 (0.025)</td></tr><tr><td>1000</td><td>0.076 (0.044)</td><td>-0.030 (0.011)</td></tr><tr><td>1500</td><td>0.069 (0.027)</td><td>-0.074 (0.0068)</td></tr><tr><td rowspan="3">z<sub>s</sub> = 2</td><td>500</td><td>0.18 (0.22)</td><td>0.063 (0.049)</td></tr><tr><td>1000</td><td>0.24 (0.086)</td><td>0.034 (0.021)</td></tr><tr><td>1500</td><td>0.27 (0.051)</td><td>-0.0081 (0.013)</td></tr></tbody></table>
<figure><img src="image_9.png" /><figcaption>FIG. 7. Expected $1\sigma$ uncertainty on $\lambda$ at each redshift [see Eq. (28)] using multipoles up to $\ell_{max} = 500$ (left) and 1500 (right) for various source redshifts, $z_s$. We do not include noise contributions in the covariance (i.e., this is the CV limit). Of course the impact of the nonlinear correction beyond the tree level increases with decreasing redshift.</figcaption></figure>

multipole, fitting function and source redshift. At high-$z$ ($z \gtrsim 1$), the constraints from the CMB lensing bi-spectrum are much tighter than those from the lensing bi-spectrum of low-$z$ sources. For both $\ell_{max}$, the best constraints on the $z$-bin of $\lambda$ are obtained from $z_s = 1$ at $z \lesssim 0.7$, from $z_s = 2$ for $0.7 \lesssim z \lesssim 1.3$, and at higher $z$, the CMB constraints are the tightest. In addition, the plot indicates by how much the constraints depend on the somewhat uncertain nonlinear correction beyond tree-level.

While lensing analyses will initially be done independ-
ently for sources at various redshift, in the long run, the
combination of measurements of the CMB and low-$z$
lensing bispectra will allow a full blind reconstruction of
the evolution of the modified gravity parameters.

IV. SUMMARY AND DISCUSSION

Besides their academic interest, modified gravity theo-
ries have been proposed as a means to interpret differently
the observational evidences for dark energy. Given the great
success otherwise of Einstein theory at all scales, from the
solar system tests of GR to the largest observational scales
probed by CMB anisotropies, deviations have to be
tenuous, and any search for such modifications are likely
difficult and not immune to a host of troublesome system-
atics effects.

Since the physics of CMB anisotropies is now well
understood, the CMB lensing bispectrum will offer a very
clean probe of modified gravity theories in the near future,
compared to other cosmological probes such as galaxy
clustering and optical weak lensing. In this paper, we
quantitatively evaluated the expected constraints on a
generic two-parameter model of the modified gravity
theories assuming the specifications of near future and
ultimate CMB experiments, which was not explored
previously. This intermediate step allows avoiding a
detailed comparison of the observation with a specific

class of modified gravity models (in much the same way
that the slope of the primordial curvature spectrum is a
good agnostic point of contact between specific theories
and observational constraints).

We quantified the information coming from various
scales and showed the impact of uncertainties in the
theoretical description of the non-linear evolution of large
scale structures. While ultimate, cosmic variance limited
experiments will, in the long run, require further theoretical
advances, we showed that our calculations are already
adequate for the next decade experiments. And we further
checked that details of the needed screening mechanism (to
satisfy Solar System constraints) are at such very small
scales that they have little impact at the scales which can be
realistically probed with such CMB observations.

Finally we point out the exciting prospect of blindly
reconstructing the redshift evolution of the distortion
parameters which we used to parametrize generically the
beyond Horndeski class of modified gravity models (e.g.,
[57]). This will be achieved by performing a joint analysis
of the bispectrum of CMB and source lensing, each
providing the best constraints within different redshift
range, in addition to allowing welcome cross-checks in
the overlap region.

ACKNOWLEDGMENTS

T. N. and A. T. thank Daisuke Yamauchi for references
and comments on modified gravity effects on the lensing.
This work was partly supported by the International
Research Unit of Advanced Future Studies at Kyoto
University. T. N. acknowledges the support from the
Ministry of Science and Technology (MOST), Taiwan,
R. O. C. through the MOST research project Grants
(No. 107-2112-M-002-002-MY3). A. T. acknowledges
financial support from MEXT/JSPS KAKENHI Grants
No. JP15H05899 and No. JP16H03977.
[1] K. S. Dawson, D. J. Schlegel et al., The baryon oscillation spectroscopic survey of SDSS-III, Astron. J. **145**, 10 (2013).

[2] K. S. Dawson, J.-P. Kneib, W. J. Percival et al., The SDSS-IV extended baryon oscillation spectroscopic survey: Overview and early data, Astron. J. **151**, 44 (2016).

[3] The DESI Collaboration, The DESI Experiment, a white-paper for Snowmass 2013, arXiv:1308.0847.

[4] T. Abbott et al. (Dark Energy Survey Collaboration), The Dark Energy Survey Collaboration, The dark energy survey: More than dark energy—an overview Mon. Not. R. Astron. Soc., **460**, 1270 (2016).

[5] L. Amendola et al., Euclid Theory Working Group, Cosmology and fundamental physics with the Euclid satellite, Living Rev. Relativity **16**, 6 (2013).

[6] V. Desjacques, D. Jeong, and F. Schmidt, Large-scale galaxy bias, Phys. Rep. **733**, 1 (2018).

[7] M. Bartelmann and A. Loeb, Gravitational lensing of quasi-stellar objects by their damped LY alpha absorbers, Astrophys. J. **457**, 529 (1996).

[8] J. Lesgourgues, M. Viel, M. G. Haehnelt, and R. Massey, A combined analysis of 3D weak lensing, Lyman-α forest and WMAP year three data, J. Cosmol. Astropart. Phys. **11** (2007) 008.

[9] R. B. Metcalf, R. A. C. Croft, and A. Romeo, Noise estimates for measurements of weak lensing from the Ly α forest, Mon. Not. R. Astron. Soc. **477**, 2841 (2018).

[10] V. Acquaviva and C. Baccigalupi, Dark energy records in lensed cosmic microwave background, Phys. Rev. D **74**, 103510 (2006).

[11] E. Calabrese, A. Cooray, M. Martinelli, A. Melchiorri, L. Pagano, A. Slosar, and G. F. Smoot, CMB lensing constraints on dark energy and modified gravity scenarios, Phys. Rev. D **80**, 103516 (2009).

[12] L. Amendola, G. Ballesteros, and V. Pettorino, Effects of modified gravity on B-mode polarization, Phys. Rev. D **90**, 043009 (2014).

[13] D. Munshi, B. Hu, A. Renzi, A. Heavens, and P. Coles, Probing modified gravity theories with ISW and CMB lensing, Mon. Not. R. Astron. Soc. **442**, 821 (2014).

[14] G. D'Amico, Z. Huang, M. Mancarella, and F. Vernizzi, Weakening gravity on redshift-survey scales with kinetic matter mixing, J. Cosmol. Astropart. Phys. 02 (2017) 014.

[15] T. Namikawa, CMB lensing bispectrum from nonlinear growth of the large scale structure, Phys. Rev. D **93**, 121301 (2016).

[16] G. Pratten and A. Lewis, Impact of Post-Born lensing on the CMB, J. Cosmol. Astropart. Phys. 08 (2016) 047.

[17] Planck Collaboration, Planck 2015 results. XIII. Cosmological parameters, Astron. Astrophys. **594**, A13 (2016).

[18] M. Takada and B. Jain, Cosmological parameters from lensing power spectrum and bispectrum tomography, Mon. Not. R. Astron. Soc. **348**, 897 (2004).

[19] Y. Takushima, A. Terukina, and K. Yamamoto, Bispectrum of cosmological density perturbations in the most general second-order scalar-tensor theory, Phys. Rev. D **89**, 104007 (2014).

[20] D. Yamauchi, S. Yokoyama, and H. Tashiro, Constraining modified theory of gravity with galaxy bispectrum, Phys. Rev. D **96**, 123516 (2017).

[21] S. Hirano, T. Kobayashi, and S. Yokoyama, Matter bispectrum beyond horndeski, Phys. Rev. D **97**, 103517 (2018).

[22] T. Namikawa, D. Yamauchi, and A. Taruya, Full-sky lensing reconstruction of gradient and curl modes from CMB maps, J. Cosmol. Astropart. Phys. 01 (2012) 007.

[23] S. Saga, D. Yamauchi, and K. Ichiki, Weak lensing induced by second-order vector mode, Phys. Rev. D **92**, 063533 (2015).

[24] A. Lewis and A. Challinor, Weak gravitational lensing of the CMB, Phys. Rep. **429**, 1 (2006).

[25] D. Hanson, A. Challinor, and A. Lewis, Weak lensing of the CMB, Gen. Relativ. Gravit. **42**, 2197 (2010).

[26] W. Hu, Mapping the dark matter through the cosmic microwave background damping tail, Astrophys. J. **557**, L79 (2001).

[27] C. M. Hirata and U. Seljak, Analyzing weak lensing of the cosmic microwave background using the likelihood function, Phys. Rev. D **67**, 043001 (2003).

[28] A. Lewis and A. Challinor, Weak gravitational lensing of the CMB, Phys. Rep. **429**, 1 (2006).

[29] D. Hanson, A. Challinor, G. Efstathiou, and P. Bielewicz, CMB temperature lensing power reconstruction, Phys. Rev. D **83**, 043005 (2011).

[30] T. Namikawa, D. Hanson, and R. Takahashi, Bias-hardened CMB lensing, Mon. Not. R. Astron. Soc. **431**, 609 (2013).

[31] B. D. Sherwin et al. (ACT Collaboration), The Atacama Cosmology Telescope: Two-season ACTPol lensing power spectrum, Phys. Rev. D **95**, 123529 (2017).

[32] BICEP2/Keck Array Collaboration, BICEP2/Keck Array VIII: Measurement of gravitational lensing from large-scale B-mode polarization, Astrophys. J. **833**, 228 (2016).

[33] Planck Collaboration, Planck 2015 results. XV. Gravitational lensing, Astron. Astrophys. **594**, A15 (2015).

[34] POLARBEAR Collaboration, Measurement of the Cosmic Microwave Background Polarization Lensing Power Spectrum with the POLARBEAR experiment, Phys. Rev. Lett. **113**, 021301 (2014).

[35] A. van Engelen et al. (SPT Collaboration), A measurement of gravitational lensing of the microwave background using South Pole Telescope data, Astrophys. J. **756**, 142 (2012).

[36] F. Bernardeau, S. Colombi, E. Gaztanaga, and R. Scoccimarro, Large scale structure of the universe and cosmological perturbation theory, Phys. Rep. **367**, 1 (2002).

[37] R. Scoccimarro and H. M. Couchman, A fitting formula for the nonlinear evolution of the bispectrum, Mon. Not. R. Astron. Soc. **325**, 1312 (2001).

[38] H. Gil-Marin, C. Wagner, F. Fragkoudi, R. Jimenez, and L. Verde, An improved fitting formula for the dark matter bispectrum, J. Cosmol. Astropart. Phys. 02 (2012) 047.

[39] G. W. Horndeski, Second-order scalar-tensor field equations in a four-dimensional space, Int. J. Theor. Phys. **10**, 363 (1974).

[40] T. Kobayashi, M. Yamaguchi, and J. Yokoyama, Generalized G-inflation—Inflation with the most general second-order field equations—, Prog. Theor. Phys. **126**, 511 (2011).

[41] C. Deffayet, X. Gao, D. A. Steer, and G. Zahariade, From k-essence to generalized Galileons, Phys. Rev. D **84**, 064039 (2011).
[42] J. Gleyzes, D. Langlois, F. Piazza, and F. Vernizzi, Exploring gravitational theories beyond Horndeski, J. Cosmol. Astropart. Phys. 02 (2015) 018.

[43] J. Gleyzes, D. Langlois, F. Piazza, and F. Vernizzi, New Class of Consistent Scalar-Tensor Theories, Phys. Rev. Lett. **114**, 211101 (2015).

[44] D. Langlois and K. Noui, Degenerate higher derivative theories beyond Horndeski: evading the Ostrogradski instability, J. Cosmol. Astropart. Phys. 02 (2016) 034.

[45] D. Langlois and K. Noui, Hamiltonian analysis of higher derivative scalar-tensor theories, J. Cosmol. Astropart. Phys. 07 (2016) 016.

[46] M. Zumalacregui and J. Garca-Bellido, Transforming gravity: From derivative couplings to matter to second-order scalar-tensor theories beyond the Horndeski Lagrangian, Phys. Rev. D **89**, 064046 (2014).

[47] K. M. Smith, D. Hanson, M. LoVerde, C. M. Hirata, and O. Zahn, Delensing CMB polarization with external datasets, J. Cosmol. Astropart. Phys. 06 (2012) 014.

[48] C. M. Hirata and U. Seljak, Reconstruction of lensing from the cosmic microwave background polarization, Phys. Rev. D **68**, 083002 (2003).

[49] A. Lewis, A. Challinor, and A. Lasenby, Efficient computation of CMB anisotropies in closed FRW models, Astrophys. J. **538**, 473 (2000).

[50] R. E. Smith, J. A. Peacock, A. Jenkins, S. D. M. White, C. S. Frenk, F. R. Pearce, P. A. Thomas, G. Efstathiou, and H. M. P. Couchman (Virgo Consortium Collaboration), Stable clustering, the halo model and nonlinear cosmological power spectra, Mon. Not. R. Astron. Soc. **341**, 1311 (2003).

[51] R. Takahashi, M. Sato, T. Nishimichi, A. Taruya, and M. Oguri, Revising the Halofit model for the nonlinear matter power spectrum, Astrophys. J. **761**, 152 (2012).

[52] B. Bose and A. Taruya, The one-loop matter bispectrum as a probe of gravity and dark energy, arXiv:1808.01120.

[53] B. Joachimi and P. Schneider, The removal of shear-ellipticity correlations from the cosmic shear signal: Influence of photometric redshift errors on the nulling technique, Astron. Astrophys. **507**, 105 (2009).

[54] T. Namikawa, T. Okamura, and A. Taruya, Magnification effect on the detection of primordial non-Gaussianity from photometric surveys, Phys. Rev. D **83**, 123514 (2011).

[55] M. Fasiello and Z. Vlah, Screening in perturbative approaches to LSS, Phys. Rev. B **773**, 236 (2017).

[56] B. R. Dinda, Weak lensing probe of cubic Galileon model, J. Cosmol. Astropart. Phys. 06 (2018) 017.

[57] J. Kennedy, L. Lombriser, and A. Taylor, Reconstructing Horndeski theories from phenomenological modified gravity and dark energy models on cosmological scales, arXiv:1804.04582.