# Neural Network Acceleration of Large-scale Structure Theory Calculations

Joseph DeRose  Shi-Fan Chen  Martin White  Nickolas Kokron

###### Abstract

We make use of neural networks to accelerate the calculation of power spectra
required for the analysis of galaxy clustering and weak gravitational lensing
data. For modern perturbation theory codes, evaluation time for a single
cosmology and redshift can take on the order of two seconds. In combination
with the comparable time required to compute linear predictions using a
Boltzmann solver, these calculations are the bottleneck for many contemporary
large-scale structure analyses. In this work, we construct neural network-
based surrogate models for Lagrangian perturbation theory (LPT) predictions of
matter power spectra, real and redshift space galaxy power spectra, and
galaxy–matter cross power spectra that attain $\sim 0.1\%$ (at one sigma)
accuracy over a broad range of scales in a $w$CDM parameter space. The neural
network surrogates can be evaluated in approximately one millisecond, a factor
of 1000 times faster than the full Boltzmann code and LPT computations. In a
simulated full-shape redshift space galaxy power spectrum analysis, we
demonstrate that the posteriors obtained using our surrogates are accurate
compared to those obtained using the full LPT model. We make our surrogate
models public at https://github.com/sfschen/EmulateLSS, so that others may
take advantage of the speed gains they provide to enable rapid iteration on
analysis settings, something that is essential in complex contemporary large-
scale structure analyses.

##  1 Introduction

Recent years have seen the maturation of perturbation-theory based models of
large-scale structure statistics [1, 2, 3, 4, 5, 6] with applications from
redshift space distortions to synergistic combinations with $N$-body
simulations for real space clustering and lensing analyses [7, 8, 9, 8, 10].
In conjunction with the theory, a number of numerical techniques have been
devised to significantly speed up the evaluation of perturbation theory
integrals due to nonlinearities in gravitational clustering. With these
advances, state-of-the-art perturbation theory codes [11, 12, 13] can produce
fully general effective-theory predictions for galaxy clustering in any given
cosmology within (roughly) two seconds.

Despite these advances, however, direct evaluation of these theory predictions
remains a significant computational expense in the analysis of clustering data
via Markov-chain Monte Carlo (MCMC) samplers. For full $\Lambda$CDM parameter
fits in particular, the amount of time taken to generate linear power spectra
and cosmological evolution parameters through Boltzmann codes like CAMB or
CLASS [14, 15] typically exceeds even that spent running the perturbation
codes, making the combined computational cost of performing linear power
spectrum and perturbation theory calculations on the order of four seconds.

Recently, a number of authors have sought to circumvent this bottleneck by
constructing fast surrogates for linear and nonlinear matter power spectra
through neural networks [16] as well as 2LPT and perturbation-theory/$N$-body
hybrids for real space galaxy power spectra [9]. The idea of emulating such
predictions through interpolating a well-sampled set of pre-computed points is
not new [17, 18, 19], and due to the smoothness of the dependence of galaxy
clustering and weak lensing power spectra on cosmological parameters even
simpler approaches such as Taylor expanding about a fiducial cosmology have
been successfully implemented [12, 20, 21, 22, 23].

Neural networks can be used as “universal” function surrogates, whose accuracy
is limited primarily by the amount of data available to train them. In this
regard the combination of neural networks and perturbative models is
particularly attractive, because it is possible to generate $10^{5}-10^{6}$
models to train the network at modest computational cost. In a sense, we pay
the cost of evaluating models in a given family up front, and such models can
be used in many analyses. This technique is gaining currency in cosmology, and
similar approaches have been presented in refs. [24, 25, 26, 27, 28, 29, 30,
16].

In this work, we focus on constructing surrogate models for all power spectra
required for standard analyses of (real and redshift space) galaxy clustering,
weak lensing, and their cross-correlation. In particular, we construct
surrogates for real space galaxy auto power spectra, galaxy–matter cross power
spectra and matter power spectra, $P_{gg}(k)$, $P_{gm}(k)$ and $P_{mm}(k)$
respectively, as predicted by convolutional Lagrangian effective field theory
(CLEFT) [3, 4], and hybrid effective field theory (HEFT) [7, 8]. Additionally,
we construct surrogates for multipoles of the redshift space galaxy auto power
spectrum $P_{\ell}(k)$ as predicted by Lagrangian effective field theory [13,
6]. The surrogates that we construct have already seen their first use as
applied to the 3D power spectrum and correlation function of BOSS galaxies
[31] and the cross-correlation between luminous red galaxies (LRGs) targeted
by the Dark Energy Spectroscopic Instrument (DESI [32]) and CMB lensing from
Planck [33] in refs. [22, 34].

The structure of this paper is as follows. In section 2 we describe the models
that we build surrogates for, and the methodology that we use to produce our
training sets. In section 3 we describe the combination of principal
components and neural networks that we use as surrogate models, and in section
2 we show that these surrogates are accurate enough for precision cosmological
constraints, concluding with some final remarks in section 5.

##  2 Parameters and Training Sets

We have trained surrogates for a variety of models relevant for the prediction
of galaxy clustering and weak lensing statistics. The following subsections
briefly describe the models for which we have built surrogates. In all cases,
we sample our parameter space using a Korobov quasi-random sequence [35],
although the choice to use Korobov sequences instead of a different
experimental design is insignificant given the large number of samples that we
are able to produce in this work.

All training data is generated with a publicly available Cobaya [36]
likelihood111https://github.com/martinjameswhite/CobayaLSS in order to
facilitate the training of surrogates for additional statistics. One
particularly useful aspect of Cobaya that we make use of in this work is built
in fast-slow parameter splitting, which allows us to generate training
examples more quickly in directions of parameter space that are fast, such as
galaxy bias parameters.

Table 1: Parameters and priors Parameter | Prior | Statistic
---|---|---
Cosmology
$\Omega_{c}h^{2}$ | flat ($0.08,0.16$) | all
$\Omega_{b}h^{2}$ | flat ($0.019,0.024$) | all
$\log 10^{10}A_{s}$ | flat ($2,4$) | all
$n_{s}$ | fixed (0.967) | all
$h$ | flat ($0.55,0.91$) | all
$\sum m_{\nu}$ | fixed (0.06 eV) 222Fixed to zero for HEFT | all
Galaxy Bias
$b_{1}$ | flat ($0,3.0$) | $P_{\rm gm},P_{\rm gg},P_{\ell}$
$b_{2}$ | flat ($-5,5$) | $P_{\rm gm},P_{\rm gg},P_{\ell}$
$b_{s}$333Fixed to zero for CLEFT | flat ($-20,20$) | $P_{\rm gm},P_{\rm gg},P_{\ell}$
$b_{k}$444HEFT only | flat ($-5,5$) | $P_{\rm gm},P_{\rm gg}$
$\alpha_{0}$ | flat ($-200,200$) | $P_{\rm gm},P_{\rm gg},P_{\ell}$
$\alpha_{2}$ | flat ($-200,200$) | $P_{\rm gm},P_{\rm gg},P_{\ell}$
$R_{h}^{3}$ | flat ($-2\times 10^{4},2\times 10^{4}$) | $P_{\rm gg},P_{\ell}$
$R_{h}^{3}\sigma^{2}$ | flat ($-1\times 10^{5},1\times 10^{5}$) | $P_{\ell}$

###  2.1 Real Space Galaxy Power Spectra

In order to model projected galaxy clustering and the cross-correlation
between galaxy density and weak lensing statistics, we also train models for
galaxy auto power spectra, $P_{\rm gg}(k)$ and galaxy–matter cross power
spectra, $P_{\rm gm}(k)$. In particular, we build surrogate models for
Convolutional Lagrangian Effective Field Theory (CLEFT) [4] and Hybrid
Effective Field Theory (HEFT) [7, 8].

In CLEFT, the galaxy auto power spectrum and galaxy–matter cross spectrum are
given by:

$\displaystyle P_{\rm gg}$
$\displaystyle=\left(1-\frac{\alpha_{a}k^{2}}{2}\right)P_{\rm Z}+P_{\rm
1-loop}+b_{1}P_{\rm b_{1}}+b_{2}P_{\rm b_{2}}+b_{1}b_{2}P_{\rm
b_{1}b_{2}}+b_{1}^{2}P_{\rm b_{1}^{2}}+b_{2}^{2}P_{\rm b_{2}^{2}}+P_{\rm SN}$
(2.1) $\displaystyle P_{\rm gm}$
$\displaystyle=\left(1-\frac{\alpha_{\times}k^{2}}{2}\right)P_{\rm Z}+P_{\rm
1-loop}+\frac{b_{1}}{2}P_{\rm b_{1}}+\frac{b_{2}}{2}P_{\rm b_{2}}\,.$ (2.2)

Note that we have set the shear bias ($b_{s}$) to zero, consistent with the
findings in [37]. Here, $P_{\rm Z}$ and $P_{\rm 1-loop}$ are the Zeldovich and
1-loop matter contributions, $P_{\rm SN}$ is the shot-noise contribution,
$b_{1}$ and $b_{2}$ are the linear and quadratic bias parameters for the
galaxy sample, and $\alpha_{\times}$ and $\alpha_{a}$ are effective field
theory corrections, necessary to marginalize over small-scale physics not
included in this model [3]. Further details on how each contribution is
computed can be found in [4, 13, 34]: the power spectra can all be computed as
convolutions of powers of linear power spectra times known kernels. We compute
these linear spectra using the CAMB Boltzmann solver [14]. As we include a
non-zero neutrino mass, we make the approximation that galaxy power spectra
follow the CDM and baryon distribution [38], performing our perturbation
theory calculations using $P_{\rm cb,lin}(k)$, rather than the power spectrum
of all matter $P_{\rm mm,lin}(k)$.

We generate predictions at 50000 points in cosmological parameter space and 20
points in bias parameter space for each point in cosmology, for a total of
$10^{6}$ different points in cosmology and bias space. We use a distinct
Korobov sequence from which to sample our bias parameters, such that we use
the $20(i)$th through the $20(i+1)$th points in the bias parameter sequence at
the $i$th point in the sequence employed for the cosmological parameters. For
each point in cosmology and bias parameter space, we produce 20 outputs in
redshift, linearly spaced between $z=0$ and $z=2$. We use the CLEFT
implementation in velocileptors555https://github.com/sfschen/velocileptors
[13] to make our predictions.

Additionally, we construct surrogates for a hybrid effective field theory
(HEFT) model for $P_{\rm gg}(k)$ and $P_{\rm gm}(k)$, as implemented in
anzu666https://github.com/kokron/anzu[8]. The expressions for the power
spectra under the HEFT model are as follows:

$\displaystyle P_{\rm gg}(k)$
$\displaystyle=\sum_{X,Y}b_{X}b_{Y}P_{XY}(k)+P_{\rm SN}$ (2.3) $\displaystyle
P_{\rm gm}(k)$ $\displaystyle=\sum_{X}b_{X}P_{1X}(k)\,,$ (2.4)

where $X,Y\in\\{1,\delta,\delta^{2},s^{2},\nabla^{2}\\}$. Here, $P_{11}$ is
the nonlinear matter power spectrum, and $P_{1,\delta}$, $P_{1,\delta^{2}}$,
$P_{\delta,\delta^{2}}$, $P_{\delta,\delta}$, $P_{\delta^{2},\delta^{2}}$ are
the nonlinear completions of $P_{\rm b_{1}}$, $P_{\rm b_{2}}$, $P_{\rm
b_{1}b_{2}}$, $P_{\rm b_{1}^{2}}$ and $P_{\rm b_{2}^{2}}$, respectively,
asymptoting to the latter at low $k$. See [8] for additional implementation
details. Training data is generated by sampling over cosmological and bias
parameters in the same way as for CLEFT, although here we restrict the
boundaries of the parameter space such that they coincide with the minimum and
maximum parameter values sampled by the Aemulus simulations [39], fixing
$N_{\rm eff}=3.046$, as our HEFT model is not valid beyond those bounds.

###  2.2 Redshift Space Galaxy Power Spectra

Finally, we build surrogate models for multipoles of the redshift space galaxy
power spectrum as predicted by the Lagrangian perturbation theory model
described in [13, 6]. This model can be thought of as the redshift space
extension of the CLEFT model discussed above, replacing real space Lagrangian
displacements with their redshift space counterparts, and expanding the basis
of included counterterms and stochastic terms to account for the anisotropy
induced by galaxy peculiar velocities. The expression for the redshift space
power spectrum under this model is as follows:

$\displaystyle P_{s}(\mathbf{k})=P_{s}^{\rm
PT}(\mathbf{k})+\left(\alpha_{0}+\alpha_{2}\mu^{2}\right)k^{2}P_{\rm
Zel}(\mathbf{k})+R_{\rm h}^{3}\left(1+\sigma^{2}k^{2}\mu^{2}\right)\,,$ (2.5)

where $P_{s}^{\rm PT}(\mathbf{k})$ is the perturbative expression for the
redshift space power spectrum, containing the same bias parameters as eq. 2.2,
while the second and third terms in eq. 2.5 are effective field theory
counterterms and stochastic terms, respectively, included to remove
sensitivity to unmodeled small scale physics. $P_{\rm Zel}(\mathbf{k})$ is the
first order (Zel’dovich approximation) redshift space matter power spectrum.
The Alcock-Paczynski (AP) effect is included as described in Appendix C.1 of
[22], and $P_{s}(\mathbf{k})$ is summarized by computing its monopole,
quadrupole and hexadecapole moments, $P_{\ell=\\{0,2,4\\}}(k)$.

The training procedure for our $P_{\ell}(k)$ surrogates is slightly different
from that of the rest of the surrogates in this work. Modeling the AP effect
requires us to assume fiducial distances and expansion rates, and choices for
these often change from analysis to analysis making it difficult to build a
surrogate model for $P_{\ell}(k)$ that is appropriate for all analyses.
Additionally, nearly all RSD analyses to date have neglected the effect of the
evolution of $P_{\ell}(k)$ within a redshift bin. As such, we have opted to
build surrogates for a single redshift, $z=0.61$, assuming a fiducial
$\Lambda$CDM cosmology with $\Omega_{m}=0.3$ to incorporate the AP effect,
although we have made it simple to train models at different redshifts and
fiducial cosmologies.

###  2.3 Nonlinear Matter Power Spectrum

Finally, we also produce a surrogate for the nonlinear matter power spectrum
predicted by CAMB [14] and Halofit [40, 41]. The only free parameters of this
model are the cosmological parameters listed in Table 1. We generate training
examples at 50000 points in cosmological parameter space and 20 points in
redshift for each cosmology, linearly spaced between $z=0$ and $z=2$.

##  3 Surrogate Modeling Methodology

We use a combination of principal component analysis (PCA) and multi-layer
perceptrons (MLP), as our surrogate models, similar to the method used in
[42]. In particular, our surrogate model takes the form:

$\displaystyle\Gamma^{ab}(k,\mathbf{\Omega})\approx\sum_{i}^{N_{\rm
PC}^{ab}}\alpha_{i}^{ab}(\mathbf{\Omega})\mathrm{PC}_{i}^{ab}(k).$ (3.1)

where $\mathbf{\Omega}$ is a vector of cosmological and bias parameters,
$\Gamma^{ab}(k,\mathbf{\Omega})=\textrm{arcsinh}[\tilde{P}^{ab}(k,\mathbf{\Omega})]$,
and $\tilde{P}^{ab}(k,\mathbf{\Omega})$ is a whitened version of
$P^{ab}(k,\mathbf{\Omega})$, the power spectra that we are modeling. In
particular, the $k$ dependence of the spectra is modeled via PCA, and the
parameter dependence of the PCA coefficients, $\alpha_{i}^{ab}$, is modeled
with an MLP. In this section we briefly outline the PCA and MLP algorithms we
apply to the training sets described in section 2 with results described in
section 2.

###  3.1 Principal Component Analysis

The first step in our surrogate modeling methodology is to construct a basis
of principal components for each model. To do so, we generate a library of
model predictions covering a broad range of cosmological and bias parameters
and redshifts as described in section 2. The input data is transformed to have
zero mean and unit variance, and an arcsinh function is applied to reduce the
dynamic range. This data is then compressed into principal components.

Let $\mathbf{X}$ be the $N\times M$ array containing our training set, where
$N$ is the number of examples in our training set and $M$ is the
dimensionality of each training data vector. Then a basis of principal
components can be constructed by computing the eigenvectors of the covariance
matrix of each model:

$\displaystyle\mathbf{C}$ $\displaystyle=\mathbf{X}^{\rm T}\mathbf{X},$ (3.2)
$\displaystyle=\mathbf{W}\mathbf{\Lambda}\mathbf{W}^{\rm T},$

where the rows of $\mathbf{W}$ are the eigenvectors, i.e., the principal
components in question and $\mathbf{\Lambda}$ is a diagonal matrix of the
eigenvalues, which are equal to the variance of the model component described
by each eigenvector.

Depending on our application, we vary the number of principal components
employed in our surrogate model to achieve our desired accuracy. The
coefficients of these principal components are then predicted by a MLP. As our
training data is virtually noiseless (modulo numerical noise), we are not
penalized for choosing $N_{\rm PC}$ to be larger than necessary. As such, we
choose the number of principal components for each application so that

$1-{\sum_{i<N_{\rm
PC}}\Lambda_{i}^{1/2}/\mathrm{Tr}(\mathbf{\Lambda}^{1/2})}<10^{-4}\,.$ (3.3)

This ensures that the principal component representation of our data is not
close to the limiting factor in the accuracy of our surrogate models.

Figure 1: For each statistic we wish to model, we choose the number of
principal components employed in our surrogate model by determining the
minimum number of components required to satisfy the inequality in equation
3.3. This figure shows the fractional error induced by truncating our PCA
representation after $N_{\rm PC}$ principal components for the six statistics
that we focus on in this work. For the real space statistics in this figure,
we use the CLEFT model, but we find similar results for HEFT.

###  3.2 Network Architecture and Optimization

We use a fully connected MLP architecture with 4 hidden layers, each with 128
neurons and an activation function given by:

$a(\mathbf{x})=\left[\mathbf{\gamma}+\frac{1-\mathbf{\gamma}}{1+\exp(-\mathbf{\beta}\odot\mathbf{x})}\right]\odot\mathbf{x}$
(3.4)

where $\odot$ indicates element-wise multiplication, following [42]. The MLP
predicts $N_{\rm PC}$ values, which are then multiplied by their respective
principal component vectors and summed in order to compute a mean squared
error loss in the whitened training data space. We selected this architecture
by performing a hyper-parameter optimization over the number of layers and
neurons per layer. The optimization was performed using $P_{4}(k)$ as training
data, as this statistic consistently produced the largest residuals. The
resulting architecture performs well for all other statistics considered in
this work.

We use one MLP per statistic, optimizing the MLP parameters using adam [43]
and a mean squared error loss function. The inputs to our network are
cosmological and bias parameters, as well as redshift where appropriate,
transformed to have zero mean and unit variance. We begin the optimization
process using a learning rate of $\eta=10^{-2}$, and we decrease this by an
order of magnitude every 30 training epochs until we reach $\eta=10^{-6}$. We
use a batch size that decreases proportionally to our learning rate, starting
out by splitting our training data into $320$ batches, and halving the batch
size every time we decrease the learning rate.

##  4 Results

Having described how we construct our surrogate models, we now discuss their
accuracy, demonstrate their utility in a realistic analysis, and compare them
to a common alternative in Taylor series expansions.

Figure 2 shows the mean, 68th percentile and 95th percentile fractional
residuals for each of our surrogate models as a function of wave number. These
residuals are computed from a set of 10000 test spectra that are held out of
the training process. For all surrogates, the mean residuals are frequently
below $0.1\%$, with the $1\sigma$ residuals approximately $20-50\%$ larger
than the mean. The $2\sigma$ residuals approach and sometimes exceed $1\%$ at
high $k$, where the model predictions vary significantly more as a function of
the nonlinear parameters included in the models. The redshift space power
spectra are the most difficult to construct surrogates for, largely due to
their larger dynamic range driven by the broad priors we allow on the bias
parameters, counterterms, and stochastic terms.

Our surrogate models each take approximately one ms to evaluate, which is a
factor of more than 1000 times faster than a call to CAMB and velocileptors,
which takes approximately 4 seconds in combination. This means that the
computational cost of generating the training data for the surrogate models is
amortized with just one MCMC, assuming the chain requires approximately one
million model evaluations.

###  4.1 Fitting the PT Challenge

In order to demonstrate the accuracy of our surrogate models in a realistic
setting, we perform an MCMC analysis on the PT challenge simulations [44] and
compare posteriors obtained with surrogate models to those obtained with the
exact models. The PT challenge simulations are a suite of 10 $N$-body
simulations designed to test large-scale, full-shape RSD analyses. Each
simulation has a side length of $3840\,h^{-1}\,\rm Mpc$ and $3072^{3}$
particles. The simulations are populated with a halo occupation distribution
in order to roughly reproduce the large scale BOSS RSD power spectrum
multipoles. The multipoles measured from these simulations have a
significantly greater precision than any currently realized measurements, and
thus provide a useful testing ground for analytic models. In this analysis, we
focus on the $z=0.61$ measurements using the provided covariance appropriate
for the volume of the BOSS survey [45]777Measurements and covariance were
obtained from https://www2.yukawa.kyoto-u.ac.jp/
takahiro.nishimichi/data/PTchallenge/.

For the purposes of this test we fix $\sum m_{\nu}=0$, as opposed to what is
described in Sec. 2, as the PT challenge simulations assume massless
neutrinos. This requires retraining the surrogate models for $P_{\ell}(k)$,
but the change in accuracy when doing so is negligible. We assume priors for
bias parameters as listed in Table 1, but alter the priors on cosmological
parameters, setting flat priors of $0.2\leq\Omega_{m}\leq 0.4$, $50\leq
H_{0}\leq 80$ and $1.61\leq\log 10^{10}A_{s}\leq 3.91$, fixing $n_{s}=0.967$
and $\Omega_{b}=f_{b}\Omega_{m}$, where $f_{b}=0.1571$, consistent with the
information provided to us by the conveners of the PT challenge. We fit the
monopole, quadrupole, and hexadecapole to $k=0.12\,h\,\rm Mpc^{-1}$.

Figure 3 compares the cosmological constraints obtained using our surrogate
models for $P_{\ell}(k)$ to those obtained using the standard velocileptors
predictions. We have shifted the posteriors to be centered at 0 to avoid
unblinding the input PT challenge cosmology. The agreement between the two
posteriors is excellent, with shifts in the mean parameter constraints of
$\leq 0.15\sigma$, and errors agreeing at better than the $10\%$ level for all
parameters.

###  4.2 Comparison with Taylor Expansions

Finally, we compare our MLP surrogate model to a Taylor series expansion, a
common alternative for approximating EFT models for redshift space power
spectra [22, 20]. To set up the Taylor series at Nth order we compute
derivatives via finite difference on a $(2N+1)^{3}$ grid in
$\Omega_{m},h,\ln(10^{10}A_{s})$ using the FinDiff888
https://findiff.readthedocs.io/en/latest/ package. Specifically, since
perturbation theory predictions for the redshift space multipoles given a set
of cosmological parameters $\Theta$ are given by inhomogeneous quadratic
polynomials in the bias and EFT parameters we can write

$P_{\ell}(k|\Theta)=\sum_{ij}\beta_{i}\beta_{j}\mathcal{M}_{ij}(\Theta),\quad\beta=(1,b_{1},b_{2},b_{s},\alpha_{0},\alpha_{2},R_{h}^{3},\sigma^{2}R_{h}^{3}).$

Since the matrix $\mathcal{M}_{ij}$ is a smooth function of the cosmological
parameters $\Theta$ we can expand it as a Taylor series, such that the full
expression is

$P_{\ell}(k|\Theta)=\sum_{ij}\beta_{i}\beta_{j}\sum_{M=0}^{N}\frac{1}{M!}(\Theta-\Theta_{0})_{a_{1}}\cdots(\Theta-\Theta_{0})_{a_{M}}\
\left.\partial_{a_{1}}\cdots\partial_{a_{M}}\mathcal{M}_{ij}\right|_{\Theta_{0}}.$
(4.1)

Note that in principle one can concatenate $\beta$ and $\Theta$ into one
“super” parameter vector to Taylor expand in, but that the mathematical
structure of perturbation theory limits the number of polynomial coefficients
in the bias parameters, similarly to how the smoothness of the multipoles in
$\Theta$ allows us to use a low-order Taylor series compared to a more generic
parametrization. Comparing a Taylor series so-constructed to the more generic
MLP surrogate described in this paper is the goal of this subsection.

Figure 4 compares the accuracy of our MLP model to that of the Taylor series
expansion with $N=4$. In particular, we compute the 68th percentile fractional
error when comparing to a test set spanning a $3\sigma$ parameter range around
the Planck 2018 best fit cosmology[46]. The MLP performs significantly better
in general, though the Taylor series often achieves errors significantly below
one percent, especially for the monopole and quadrupole, where measurements
are most precise. The dashed lines show the same quantities, but computed over
the full parameter range that our MLP is trained on. The MLP performs
similarly well in this expanded parameter range, while (unsurprisingly) the
Taylor series fares significantly worse with errors often exceeding one
percent. Reducing the error would require us to produce multiple, overlapping
Taylor expansions with different central points. This increases the complexity
of the method and introduces more hyperparameters to optimize.

Figure 2: Fractional residuals of our surrogates when comparing to a test set
drawn from the same prior bounds as the training set, but held out of the MLP
optimization procedure. We display the median absolute deviation (MAD), 68th
percentile and 95th percentile residuals as blue, orange and green lines
respectively for each statistic. The RSD statistics are the most difficult to
construct surrogates for, particularly $P_{2}$ and $P_{4}$, as they have a
large dynamic range due to the counterterms and stochastic terms used in the
EFT model that we employ. The residuals are made more significant by the zero
crossings of $P_{2}$ and $P_{4}$, where the denominator of the fractional
error becomes small. For the real space spectra, we display the residuals for
our CLEFT surrogates. Our HEFT surrogates achieve very comparable levels of
accuracy. Figure 3: Constraints on the three cosmological parameters left free
when fitting the PT challenge RSD multipoles. The black contours are obtained
with our surrogate models, while the red use velocileptors. Contours are
shifted to be centered at 0, so as to not reveal the input cosmology of the PT
challenge simulations. The agreement is excellent. Figure 4: Comparison of our
MLP surrogate model to a Taylor expansion for LPT predictions of
$P_{\ell}(k)$. 68th percentile fractional residuals of our MLP (blue) and
Taylor series (orange) surrogate models when comparing to our full test set
(dashed) or a test set restricted to a three sigma region in cosmological
parameter space around the best fit cosmology from [46] (solid). In general,
the MLP outperforms the Taylor series, but the surrogates achieve more
comparable levels of accuracy in the restricted test set, where accuracy is
most important given our prior knowledge of cosmological parameters. The
Taylor series surrogate is a factor of $\mathcal{O}(10^{3})$ times less
expensive to generate.

##  5 Conclusion

In this work we have combined principal component analysis (PCA) and multi-
layer perceptrons (MLP) to construct surrogate models for commonly used
ingredients in galaxy clustering and weak lensing analyses. In particular, we
have constructed surrogates for real space galaxy–galaxy, galaxy–matter, and
matter–matter power spectra, and multipoles of the redshift space galaxy auto
power spectrum as described in section 2. We use CLEFT [4] and HEFT [7, 8] to
predict the galaxy–galaxy and galaxy–matter power spectra, halofit to predict
the nonlinear matter power spectrum [40, 41], and the full LPT model discussed
in [6] to predict the redshift space galaxy power spectra. We train these
surrogates over a wide cosmology and bias parameter space to facilitate a
broad range of analyses. In section 2, we show that these surrogates achieve
mean fractional errors that are at or below $0.1\%$ for most scales. We
demonstrate the accuracy of our surrogate models by performing a simulated
analysis on redshift space power spectrum multipoles measured from the PT
challenge boxes [44] in section 4.1. Finally, we compare to a Taylor series
expansion of redshift space power spectrum multipoles, demonstrating that our
surrogate models achieve significantly better accuracy than these commonly
used alternatives. We make our surrogate models public at
github.com/sfschen/EmulateLSS, in the hope that they will be widely applicable
for ongoing cosmological analyses.

The computational expense of parameter estimation is rapidly becoming a
limiting factor in analyses of large-scale structure. This is largely due to
the increasing complexity of models driven by the stringent accuracy
requirements imposed by growing measurement precision. As such, it is
important to find ways to expedite the parameter estimation step of analyses.
In this work, we have focused on doing this by significantly speeding up the
prediction of power spectra by constructing PCA+MLP surrogate models. One
could also consider building surrogates for the final observables such as
angular power spectra, correlation functions, or window-convolved power
spectra, but doing so involves additional tracer and survey specific
information that would require building new surrogates for each analysis. This
may be beneficial in the case where projection or window convolution
operations make up a significant portion of a single likelihood evaluation,
e.g., if the Limber approximation is not applicable[47]. Acceleration of
theory calculations in the manner presented in this work will significantly
expedite complex contemporary cosmology analyses, saving countless hours of
waiting for MCMCs to converge, and empowering scientists to perform more
robust inference.

## Acknowledgments

J.D. is supported by the Chamberlain Fellowship at Lawrence Berkeley National
Laboratory. S.C and M.W. are supported by the DOE and the NSF. N.K. is
supported by the Gerald J. Lieberman fellowship. We acknowledge the use of
Cobaya [36, 48], GetDist [49] and CAMB [50] and thank their authors for making
these products public. This research has made use of NASA’s Astrophysics Data
System and the arXiv preprint server. This research used resources of the
National Energy Research Scientific Computing Center (NERSC), a U.S.
Department of Energy Office of Science User Facility operated under Contract
No. DE-AC02-05CH11231.

## References

  * [1] P. McDonald and A. Roy, Clustering of dark matter tracers: generalizing bias for the coming era of precision LSS, JCAP 2009 (August, 2009) 020, [arXiv:0902.0991].
  * [2] D. Baumann, A. Nicolis, L. Senatore, and M. Zaldarriaga, Cosmological non-linearities as an effective fluid, JCAP 2012 (July, 2012) 051, [arXiv:1004.2488].
  * [3] Z. Vlah, M. White, and A. Aviles, A Lagrangian effective field theory, JCAP 2015 (September, 2015) 014, [arXiv:1506.05264].
  * [4] Z. Vlah, E. Castorina, and M. White, The Gaussian streaming model and convolution Lagrangian effective field theory, JCAP 2016 (December, 2016) 007, [arXiv:1609.02908].
  * [5] M. Lewandowski, L. Senatore, F. Prada, C. Zhao, and C.-H. Chuang, EFT of large scale structures in redshift space, PRD 97 (March, 2018) 063526, [arXiv:1512.06831].
  * [6] S.-F. Chen, Z. Vlah, E. Castorina, and M. White, Redshift-Space Distortions in Lagrangian Perturbation Theory, arXiv e-prints (Dec., 2020) arXiv:2012.04636, [arXiv:2012.04636].
  * [7] C. Modi, S.-F. Chen, and M. White, Simulations and symmetries, MNRAS 492 (March, 2020) 5754–5763, [arXiv:1910.07097].
  * [8] N. Kokron, J. DeRose, S.-F. Chen, M. White, and R. H. Wechsler, The cosmology dependence of galaxy clustering and lensing from a hybrid N-body-perturbation theory model, MNRAS 505 (July, 2021) 1422–1440, [arXiv:2101.11014].
  * [9] M. Zennaro, R. E. Angulo, M. Pellejero-Ibáñez, J. Stücker, S. Contreras, and G. Aricò, The BACCO simulation project: biased tracers in real space, arXiv e-prints (Jan., 2021) arXiv:2101.12187, [arXiv:2101.12187].
  * [10] B. Hadzhiyska, C. García-García, D. Alonso, A. Nicola, and A. Slosar, Hefty enhancement of cosmological constraints from the DES Y1 data using a hybrid effective field theory approach to galaxy bias, JCAP 2021 (Sept., 2021) 020, [arXiv:2103.09820].
  * [11] A. Chudaykin, M. M. Ivanov, O. H. E. Philcox, and M. Simonović, Nonlinear perturbation theory extension of the Boltzmann code CLASS, PRD 102 (Sept., 2020) 063533, [arXiv:2004.10607].
  * [12] P. Zhang and G. D’Amico, PyBird: Python code for biased tracers in redshift space, Apr., 2021.
  * [13] S.-F. Chen, Z. Vlah, and M. White, Consistent modeling of velocity statistics and redshift-space distortions in one-loop perturbation theory, JCAP 2020 (July, 2020) 062, [arXiv:2005.00523].
  * [14] A. Lewis, A. Challinor, and A. Lasenby, Efficient Computation of Cosmic Microwave Background Anisotropies in Closed Friedmann-Robertson-Walker Models, ApJ 538 (Aug., 2000) 473–476, [astro-ph/9911177].
  * [15] D. Blas, J. Lesgourgues, and T. Tram, The Cosmic Linear Anisotropy Solving System (CLASS). Part II: Approximation schemes, JCAP 2011 (July, 2011) 034, [arXiv:1104.2933].
  * [16] A. Spurio Mancini, D. Piras, J. Alsing, B. Joachimi, and M. P. Hobson, $\it{CosmoPower}\,$: emulating cosmological power spectra for accelerated Bayesian inference from next-generation surveys, arXiv e-prints (June, 2021) arXiv:2106.03846, [arXiv:2106.03846].
  * [17] W. A. Fendt and B. D. Wandelt, Pico: Parameters for the Impatient Cosmologist, ApJ 654 (January, 2007) 2–11, [astro-ph/0606709].
  * [18] K. Heitmann, D. Higdon, C. Nakhleh, and S. Habib, Cosmic Calibration, ApJL 646 (July, 2006) L1–L4, [astro-ph/0606154].
  * [19] S. Agarwal, F. B. Abdalla, H. A. Feldman, O. Lahav, and S. A. Thomas, PkANN - II. A non-linear matter power spectrum interpolator developed using artificial neural networks, MNRAS 439 (April, 2014) 2102–2121, [arXiv:1312.2101].
  * [20] M. Cataneo, S. Foreman, and L. Senatore, Efficient exploration of cosmology dependence in the EFT of LSS, JCAP 2017 (Apr., 2017) 026, [arXiv:1606.03633].
  * [21] T. Colas, G. d’Amico, L. Senatore, P. Zhang, and F. Beutler, Efficient cosmological analysis of the SDSS/BOSS data from the Effective Field Theory of Large-Scale Structure, JCAP 2020 (June, 2020) 001, [arXiv:1909.07951].
  * [22] S.-F. Chen, Z. Vlah, and M. White, A new analysis of the BOSS survey, including full-shape information and post-reconstruction BAO, arXiv e-prints (October, 2021) arXiv:2110.05530, [arXiv:2110.05530].
  * [23] K. Osato, T. Nishimichi, A. Taruya, and F. Bernardeau, Implementing spectra response function approaches for fast calculation of power spectra and bispectra, arXiv e-prints (July, 2021) arXiv:2107.04275, [arXiv:2107.04275].
  * [24] T. Auld, M. Bridges, and M. P. Hobson, COSMONET: fast cosmological parameter estimation in non-flat models using neural networks, MNRAS 387 (July, 2008) 1575–1582, [astro-ph/0703445].
  * [25] S. Agarwal, F. B. Abdalla, H. A. Feldman, O. Lahav, and S. A. Thomas, PkANN - II. A non-linear matter power spectrum interpolator developed using artificial neural networks, MNRAS 439 (Apr., 2014) 2102–2121, [arXiv:1312.2101].
  * [26] J. Albers, C. Fidler, J. Lesgourgues, N. Schöneberg, and J. Torrado, CosmicNet. Part I. Physics-driven implementation of neural networks within Einstein-Boltzmann Solvers, JCAP 2019 (Sept., 2019) 028, [arXiv:1907.05764].
  * [27] A. Manrique-Yus and E. Sellentin, Euclid-era cosmology for everyone: neural net assisted MCMC sampling for the joint 3 × 2 likelihood, MNRAS 491 (Jan., 2020) 2655–2663, [arXiv:1907.05881].
  * [28] M. F. Kasim, D. Watson-Parris, L. Deaconu, S. Oliver, P. Hatfield, D. H. Froula, G. Gregori, M. Jarvis, S. Khatiwala, J. Korenaga, J. Topp-Mugglestone, E. Viezzer, and S. M. Vinko, Building high accuracy emulators for scientific simulations with deep neural architecture search, arXiv e-prints (Jan., 2020) arXiv:2001.08055, [arXiv:2001.08055].
  * [29] G. Aricò, R. E. Angulo, and M. Zennaro, Accelerating Large-Scale-Structure data analyses by emulating Boltzmann solvers and Lagrangian Perturbation Theory, arXiv e-prints (Apr., 2021) arXiv:2104.14568, [arXiv:2104.14568].
  * [30] R. E. Angulo, M. Zennaro, S. Contreras, G. Aricò, M. Pellejero-Ibañez, and J. Stücker, The BACCO Simulation Project: Exploiting the full power of large-scale structure for cosmology, arXiv e-prints (Apr., 2020) arXiv:2004.06245, [arXiv:2004.06245].
  * [31] K. S. Dawson, D. J. Schlegel, C. P. Ahn, S. F. Anderson, É. Aubourg, S. Bailey, R. H. Barkhouser, J. E. Bautista, A. Beifiori, A. A. Berlind, V. Bhardwaj, D. Bizyaev, C. H. Blake, M. R. Blanton, M. Blomqvist, A. S. Bolton, A. Borde, J. Bovy, W. N. Brandt, H. Brewington, J. Brinkmann, P. J. Brown, J. R. Brownstein, K. Bundy, N. G. Busca, W. Carithers, A. R. Carnero, M. A. Carr, Y. Chen, J. Comparat, N. Connolly, F. Cope, R. A. C. Croft, A. J. Cuesta, L. N. da Costa, J. R. A. Davenport, T. Delubac, R. de Putter, S. Dhital, A. Ealet, G. L. Ebelke, D. J. Eisenstein, S. Escoffier, X. Fan, N. Filiz Ak, H. Finley, A. Font-Ribera, R. Génova-Santos, J. E. Gunn, H. Guo, D. Haggard, P. B. Hall, J.-C. Hamilton, B. Harris, D. W. Harris, S. Ho, D. W. Hogg, D. Holder, K. Honscheid, J. Huehnerhoff, B. Jordan, W. P. Jordan, G. Kauffmann, E. A. Kazin, D. Kirkby, M. A. Klaene, J.-P. Kneib, J.-M. Le Goff, K.-G. Lee, D. C. Long, C. P. Loomis, B. Lundgren, R. H. Lupton, M. A. G. Maia, M. Makler, E. Malanushenko, V. Malanushenko, R. Mandelbaum, M. Manera, C. Maraston, D. Margala, K. L. Masters, C. K. McBride, P. McDonald, I. D. McGreer, R. G. McMahon, O. Mena, J. Miralda-Escudé, A. D. Montero-Dorta, F. Montesano, D. Muna, A. D. Myers, T. Naugle, R. C. Nichol, P. Noterdaeme, S. E. Nuza, M. D. Olmstead, A. Oravetz, D. J. Oravetz, R. Owen, N. Padmanabhan, N. Palanque-Delabrouille, K. Pan, J. K. Parejko, I. Pâris, W. J. Percival, I. Pérez-Fournon, I. Pérez-Ràfols, P. Petitjean, R. Pfaffenberger, J. Pforr, M. M. Pieri, F. Prada, A. M. Price-Whelan, M. J. Raddick, R. Rebolo, J. Rich, G. T. Richards, C. M. Rockosi, N. A. Roe, A. J. Ross, N. P. Ross, G. Rossi, J. A. Rubiño-Martin, L. Samushia, A. G. Sánchez, C. Sayres, S. J. Schmidt, D. P. Schneider, C. G. Scóccola, H.-J. Seo, A. Shelden, E. Sheldon, Y. Shen, Y. Shu, A. Slosar, S. A. Smee, S. A. Snedden, F. Stauffer, O. Steele, M. A. Strauss, A. Streblyanska, N. Suzuki, M. E. C. Swanson, T. Tal, M. Tanaka, D. Thomas, J. L. Tinker, R. Tojeiro, C. A. Tremonti, M. Vargas Magaña, L. Verde, M. Viel, D. A. Wake, M. Watson, B. A. Weaver, D. H. Weinberg, B. J. Weiner, A. A. West, M. White, W. M. Wood-Vasey, C. Yeche, I. Zehavi, G.-B. Zhao, and Z. Zheng, The Baryon Oscillation Spectroscopic Survey of SDSS-III, AJ 145 (Jan., 2013) 10, [arXiv:1208.0022].
  * [32] DESI Collaboration, A. Aghamousa, J. Aguilar, S. Ahlen, S. Alam, L. E. Allen, C. Allende Prieto, J. Annis, S. Bailey, C. Balland, O. Ballester, C. Baltay, L. Beaufore, C. Bebek, T. C. Beers, E. F. Bell, J. L. Bernal, R. Besuner, F. Beutler, C. Blake, H. Bleuler, M. Blomqvist, R. Blum, A. S. Bolton, C. Briceno, D. Brooks, J. R. Brownstein, E. Buckley-Geer, A. Burden, E. Burtin, N. G. Busca, R. N. Cahn, Y.-C. Cai, L. Cardiel-Sas, R. G. Carlberg, P.-H. Carton, R. Casas, F. J. Castander, J. L. Cervantes-Cota, T. M. Claybaugh, M. Close, C. T. Coker, S. Cole, J. Comparat, A. P. Cooper, M. C. Cousinou, M. Crocce, J.-G. Cuby, D. P. Cunningham, T. M. Davis, K. S. Dawson, A. de la Macorra, J. De Vicente, T. Delubac, M. Derwent, A. Dey, G. Dhungana, Z. Ding, P. Doel, Y. T. Duan, A. Ealet, J. Edelstein, S. Eftekharzadeh, D. J. Eisenstein, A. Elliott, S. Escoffier, M. Evatt, P. Fagrelius, X. Fan, K. Fanning, A. Farahi, J. Farihi, G. Favole, Y. Feng, E. Fernandez, J. R. Findlay, D. P. Finkbeiner, M. J. Fitzpatrick, B. Flaugher, S. Flender, A. Font-Ribera, J. E. Forero-Romero, P. Fosalba, C. S. Frenk, M. Fumagalli, B. T. Gaensicke, G. Gallo, J. Garcia-Bellido, E. Gaztanaga, N. Pietro Gentile Fusillo, T. Gerard, I. Gershkovich, T. Giannantonio, D. Gillet, G. Gonzalez-de-Rivera, V. Gonzalez-Perez, S. Gott, O. Graur, G. Gutierrez, J. Guy, S. Habib, H. Heetderks, I. Heetderks, K. Heitmann, W. A. Hellwing, D. A. Herrera, S. Ho, S. Holland, K. Honscheid, E. Huff, T. A. Hutchinson, D. Huterer, H. S. Hwang, J. M. Illa Laguna, Y. Ishikawa, D. Jacobs, N. Jeffrey, P. Jelinsky, E. Jennings, L. Jiang, J. Jimenez, J. Johnson, R. Joyce, E. Jullo, S. Juneau, S. Kama, A. Karcher, S. Karkar, R. Kehoe, N. Kennamer, S. Kent, M. Kilbinger, A. G. Kim, D. Kirkby, T. Kisner, E. Kitanidis, J.-P. Kneib, S. Koposov, E. Kovacs, K. Koyama, A. Kremin, R. Kron, L. Kronig, A. Kueter-Young, C. G. Lacey, R. Lafever, O. Lahav, A. Lambert, M. Lampton, M. Landriau, D. Lang, T. R. Lauer, J.-M. Le Goff, L. Le Guillou, A. Le Van Suu, J. H. Lee, S.-J. Lee, D. Leitner, M. Lesser, M. E. Levi, B. L’Huillier, B. Li, M. Liang, H. Lin, E. Linder, S. R. Loebman, Z. Lukić, J. Ma, N. MacCrann, C. Magneville, L. Makarem, M. Manera, C. J. Manser, R. Marshall, P. Martini, R. Massey, T. Matheson, J. McCauley, P. McDonald, I. D. McGreer, A. Meisner, N. Metcalfe, T. N. Miller, R. Miquel, J. Moustakas, A. Myers, M. Naik, J. A. Newman, R. C. Nichol, A. Nicola, L. Nicolati da Costa, J. Nie, G. Niz, P. Norberg, B. Nord, D. Norman, P. Nugent, T. O’Brien, M. Oh, K. A. G. Olsen, C. Padilla, H. Padmanabhan, N. Padmanabhan, N. Palanque-Delabrouille, A. Palmese, D. Pappalardo, I. Pâris, C. Park, A. Patej, J. A. Peacock, H. V. Peiris, X. Peng, W. J. Percival, S. Perruchot, M. M. Pieri, R. Pogge, J. E. Pollack, C. Poppett, F. Prada, A. Prakash, R. G. Probst, D. Rabinowitz, A. Raichoor, C. H. Ree, A. Refregier, X. Regal, B. Reid, K. Reil, M. Rezaie, C. M. Rockosi, N. Roe, S. Ronayette, A. Roodman, A. J. Ross, N. P. Ross, G. Rossi, E. Rozo, V. Ruhlmann-Kleider, E. S. Rykoff, C. Sabiu, L. Samushia, E. Sanchez, J. Sanchez, D. J. Schlegel, M. Schneider, M. Schubnell, A. Secroun, U. Seljak, H.-J. Seo, S. Serrano, A. Shafieloo, H. Shan, R. Sharples, M. J. Sholl, W. V. Shourt, J. H. Silber, D. R. Silva, M. M. Sirk, A. Slosar, A. Smith, G. F. Smoot, D. Som, Y.-S. Song, D. Sprayberry, R. Staten, A. Stefanik, G. Tarle, S. Sien Tie, J. L. Tinker, R. Tojeiro, F. Valdes, O. Valenzuela, M. Valluri, M. Vargas-Magana, L. Verde, A. R. Walker, J. Wang, Y. Wang, B. A. Weaver, C. Weaverdyck, R. H. Wechsler, D. H. Weinberg, M. White, Q. Yang, C. Yeche, T. Zhang, G.-B. Zhao, Y. Zheng, X. Zhou, Z. Zhou, Y. Zhu, H. Zou, and Y. Zu, The DESI Experiment Part I: Science,Targeting, and Survey Design, arXiv e-prints (Oct., 2016) arXiv:1611.00036, [arXiv:1611.00036].
  * [33] Planck Collaboration, N. Aghanim, Y. Akrami, M. Ashdown, J. Aumont, C. Baccigalupi, M. Ballardini, A. J. Banday, R. B. Barreiro, N. Bartolo, S. Basak, K. Benabed, J.-P. Bernard, M. Bersanelli, P. Bielewicz, J. J. Bock, J. R. Bond, J. Borrill, F. R. Bouchet, F. Boulanger, M. Bucher, C. Burigana, E. Calabrese, J.-F. Cardoso, J. Carron, A. Challinor, H. C. Chiang, L. P. L. Colombo, C. Combet, B. P. Crill, F. Cuttaia, P. de Bernardis, G. de Zotti, J. Delabrouille, E. Di Valentino, J. M. Diego, O. Doré, M. Douspis, A. Ducout, X. Dupac, G. Efstathiou, F. Elsner, T. A. Enßlin, H. K. Eriksen, Y. Fantaye, R. Fernandez-Cobos, F. Forastieri, M. Frailis, A. A. Fraisse, E. Franceschi, A. Frolov, S. Galeotta, S. Galli, K. Ganga, R. T. Génova-Santos, M. Gerbino, T. Ghosh, J. González-Nuevo, K. M. Górski, S. Gratton, A. Gruppuso, J. E. Gudmundsson, J. Hamann, W. Handley, F. K. Hansen, D. Herranz, E. Hivon, Z. Huang, A. H. Jaffe, W. C. Jones, A. Karakci, E. Keihänen, R. Keskitalo, K. Kiiveri, J. Kim, L. Knox, N. Krachmalnicoff, M. Kunz, H. Kurki-Suonio, G. Lagache, J.-M. Lamarre, A. Lasenby, M. Lattanzi, C. R. Lawrence, M. Le Jeune, F. Levrier, A. Lewis, M. Liguori, P. B. Lilje, V. Lindholm, M. López-Caniego, P. M. Lubin, Y.-Z. Ma, J. F. Macías-Pérez, G. Maggio, D. Maino, N. Mandolesi, A. Mangilli, A. Marcos-Caballero, M. Maris, P. G. Martin, E. Martínez-González, S. Matarrese, N. Mauri, J. D. McEwen, A. Melchiorri, A. Mennella, M. Migliaccio, M.-A. Miville-Deschênes, D. Molinari, A. Moneti, L. Montier, G. Morgante, A. Moss, P. Natoli, L. Pagano, D. Paoletti, B. Partridge, G. Patanchon, F. Perrotta, V. Pettorino, F. Piacentini, L. Polastri, G. Polenta, J.-L. Puget, J. P. Rachen, M. Reinecke, M. Remazeilles, A. Renzi, G. Rocha, C. Rosset, G. Roudier, J. A. Rubiño-Martín, B. Ruiz-Granados, L. Salvati, M. Sandri, M. Savelainen, D. Scott, C. Sirignano, R. Sunyaev, A.-S. Suur-Uski, J. A. Tauber, D. Tavagnacco, M. Tenti, L. Toffolatti, M. Tomasi, T. Trombetti, J. Valiviita, B. Van Tent, P. Vielva, F. Villa, N. Vittorio, B. D. Wandelt, I. K. Wehus, M. White, S. D. M. White, A. Zacchei, and A. Zonca, Planck 2018 results. VIII. Gravitational lensing, ArXiv e-prints (July, 2018) [arXiv:1807.06210].
  * [34] M. White, R. Zhou, J. DeRose, S. Ferraro, S.-F. Chen, N. Kokron, S. Bailey, D. Brooks, J. Garcia-Bellido, J. Guy, K. Honscheid, R. Kehoe, A. Kremin, M. Levi, N. Palanque-Delabrouille, C. Poppett, D. Schlegel, and G. Tarle, Cosmological constraints from the tomographic cross-correlation of DESI Luminous Red Galaxies and Planck CMB lensing, arXiv e-prints (Nov., 2021) arXiv:2111.09898, [arXiv:2111.09898].
  * [35] A. Korobov, The approximate computation of multiple integrals, in Dokl. Akad. Nauk SSSR, vol. 124, pp. 1207–1210, 1959.
  * [36] J. Torrado and A. Lewis, Cobaya: Bayesian analysis in cosmology, Oct., 2019.
  * [37] C. Modi, M. White, and Z. Vlah, Modeling CMB lensing cross correlations with CLEFT, JCAP 2017 (August, 2017) 009, [arXiv:1706.03173].
  * [38] E. Castorina, C. Carbone, J. Bel, E. Sefusatti, and K. Dolag, DEMNUni: the clustering of large-scale structures in the presence of massive neutrinos, JCAP 2015 (July, 2015) 043, [arXiv:1505.07148].
  * [39] J. DeRose, R. H. Wechsler, J. L. Tinker, M. R. Becker, Y.-Y. Mao, T. McClintock, S. McLaughlin, E. Rozo, and Z. Zhai, The AEMULUS Project. I. Numerical Simulations for Precision Cosmology, ApJ 875 (Apr., 2019) 69, [arXiv:1804.05865].
  * [40] 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, Stable clustering, the halo model and non-linear cosmological power spectra, MNRAS 341 (June, 2003) 1311–1332, [astro-ph/0207664].
  * [41] R. Takahashi, M. Sato, T. Nishimichi, A. Taruya, and M. Oguri, Revising the Halofit Model for the Nonlinear Matter Power Spectrum, ApJ 761 (December, 2012) 152, [arXiv:1208.2701].
  * [42] J. Alsing, H. Peiris, J. Leja, C. Hahn, R. Tojeiro, D. Mortlock, B. Leistedt, B. D. Johnson, and C. Conroy, SPECULATOR: Emulating Stellar Population Synthesis for Fast and Accurate Galaxy Spectra and Photometry, ApJS 249 (July, 2020) 5, [arXiv:1911.11778].
  * [43] D. P. Kingma and J. Ba, Adam: A Method for Stochastic Optimization, arXiv e-prints (December, 2014) arXiv:1412.6980, [arXiv:1412.6980].
  * [44] T. Nishimichi, G. D’Amico, M. M. Ivanov, L. Senatore, M. Simonović, M. Takada, M. Zaldarriaga, and P. Zhang, Blinded challenge for precision cosmology with large-scale structure: Results from effective field theory for the redshift-space galaxy power spectrum, PRD 102 (December, 2020) 123541, [arXiv:2003.08277].
  * [45] K. S. Dawson, J.-P. Kneib, W. J. Percival, S. Alam, F. D. Albareti, S. F. Anderson, E. Armengaud, É. Aubourg, S. Bailey, J. E. Bautista, A. A. Berlind, M. A. Bershady, F. Beutler, D. Bizyaev, M. R. Blanton, M. Blomqvist, A. S. Bolton, J. Bovy, W. N. Brandt, J. Brinkmann, J. R. Brownstein, E. Burtin, N. G. Busca, Z. Cai, C.-H. Chuang, N. Clerc, J. Comparat, F. Cope, R. A. C. Croft, I. Cruz-Gonzalez, L. N. da Costa, M.-C. Cousinou, J. Darling, A. de la Macorra, S. de la Torre, T. Delubac, H. du Mas des Bourboux, T. Dwelly, A. Ealet, D. J. Eisenstein, M. Eracleous, S. Escoffier, X. Fan, A. Finoguenov, A. Font-Ribera, P. Frinchaboy, P. Gaulme, A. Georgakakis, P. Green, H. Guo, J. Guy, S. Ho, D. Holder, J. Huehnerhoff, T. Hutchinson, Y. Jing, E. Jullo, V. Kamble, K. Kinemuchi, D. Kirkby, F.-S. Kitaura, M. A. Klaene, R. R. Laher, D. Lang, P. Laurent, J.-M. Le Goff, C. Li, Y. Liang, M. Lima, Q. Lin, W. Lin, Y.-T. Lin, D. C. Long, B. Lundgren, N. MacDonald, M. A. Geimba Maia, E. Malanushenko, V. Malanushenko, V. Mariappan, C. K. McBride, I. D. McGreer, B. Ménard, A. Merloni, A. Meza, A. D. Montero-Dorta, D. Muna, A. D. Myers, K. Nandra, T. Naugle, J. A. Newman, P. Noterdaeme, P. Nugent, R. Ogando, M. D. Olmstead, A. Oravetz, D. J. Oravetz, N. Padmanabhan, N. Palanque-Delabrouille, K. Pan, J. K. Parejko, I. Pâris, J. A. Peacock, P. Petitjean, M. M. Pieri, A. Pisani, F. Prada, A. Prakash, A. Raichoor, B. Reid, J. Rich, J. Ridl, S. Rodriguez-Torres, A. Carnero Rosell, A. J. Ross, G. Rossi, J. Ruan, M. Salvato, C. Sayres, D. P. Schneider, D. J. Schlegel, U. Seljak, H.-J. Seo, B. Sesar, S. Shandera, Y. Shu, A. Slosar, F. Sobreira, A. Streblyanska, N. Suzuki, D. Taylor, C. Tao, J. L. Tinker, R. Tojeiro, M. Vargas-Magaña, Y. Wang, B. A. Weaver, D. H. Weinberg, M. White, W. M. Wood-Vasey, C. Yeche, Z. Zhai, C. Zhao, G.-b. Zhao, Z. Zheng, G. Ben Zhu, and H. Zou, The SDSS-IV Extended Baryon Oscillation Spectroscopic Survey: Overview and Early Data, AJ 151 (Feb., 2016) 44, [arXiv:1508.04473].
  * [46] Planck Collaboration, N. Aghanim, Y. Akrami, M. Ashdown, J. Aumont, C. Baccigalupi, M. Ballardini, A. J. Banday, R. B. Barreiro, N. Bartolo, S. Basak, R. Battye, K. Benabed, J. P. Bernard, M. Bersanelli, P. Bielewicz, J. J. Bock, J. R. Bond, J. Borrill, F. R. Bouchet, F. Boulanger, M. Bucher, C. Burigana, R. C. Butler, E. Calabrese, J. F. Cardoso, J. Carron, A. Challinor, H. C. Chiang, J. Chluba, L. P. L. Colombo, C. Combet, D. Contreras, B. P. Crill, F. Cuttaia, P. de Bernardis, G. de Zotti, J. Delabrouille, J. M. Delouis, E. Di Valentino, J. M. Diego, O. Doré, M. Douspis, A. Ducout, X. Dupac, S. Dusini, G. Efstathiou, F. Elsner, T. A. Enßlin, H. K. Eriksen, Y. Fantaye, M. Farhang, J. Fergusson, R. Fernandez-Cobos, F. Finelli, F. Forastieri, M. Frailis, A. A. Fraisse, E. Franceschi, A. Frolov, S. Galeotta, S. Galli, K. Ganga, R. T. Génova-Santos, M. Gerbino, T. Ghosh, J. González-Nuevo, K. M. Górski, S. Gratton, A. Gruppuso, J. E. Gudmundsson, J. Hamann, W. Handley, F. K. Hansen, D. Herranz, S. R. Hildebrandt, E. Hivon, Z. Huang, A. H. Jaffe, W. C. Jones, A. Karakci, E. Keihänen, R. Keskitalo, K. Kiiveri, J. Kim, T. S. Kisner, L. Knox, N. Krachmalnicoff, M. Kunz, H. Kurki-Suonio, G. Lagache, J. M. Lamarre, A. Lasenby, M. Lattanzi, C. R. Lawrence, M. Le Jeune, P. Lemos, J. Lesgourgues, F. Levrier, A. Lewis, M. Liguori, P. B. Lilje, M. Lilley, V. Lindholm, M. López-Caniego, P. M. Lubin, Y. Z. Ma, J. F. Macías-Pérez, G. Maggio, D. Maino, N. Mandolesi, A. Mangilli, A. Marcos-Caballero, M. Maris, P. G. Martin, M. Martinelli, E. Martínez-González, S. Matarrese, N. Mauri, J. D. McEwen, P. R. Meinhold, A. Melchiorri, A. Mennella, M. Migliaccio, M. Millea, S. Mitra, M. A. Miville-Deschênes, D. Molinari, L. Montier, G. Morgante, A. Moss, P. Natoli, H. U. Nørgaard-Nielsen, L. Pagano, D. Paoletti, B. Partridge, G. Patanchon, H. V. Peiris, F. Perrotta, V. Pettorino, F. Piacentini, L. Polastri, G. Polenta, J. L. Puget, J. P. Rachen, M. Reinecke, M. Remazeilles, A. Renzi, G. Rocha, C. Rosset, G. Roudier, J. A. Rubiño-Martín, B. Ruiz-Granados, L. Salvati, M. Sandri, M. Savelainen, D. Scott, E. P. S. Shellard, C. Sirignano, G. Sirri, L. D. Spencer, R. Sunyaev, A. S. Suur-Uski, J. A. Tauber, D. Tavagnacco, M. Tenti, L. Toffolatti, M. Tomasi, T. Trombetti, L. Valenziano, J. Valiviita, B. Van Tent, L. Vibert, P. Vielva, F. Villa, N. Vittorio, B. D. Wandelt, I. K. Wehus, M. White, S. D. M. White, A. Zacchei, and A. Zonca, Planck 2018 results. VI. Cosmological parameters, A&A 641 (September, 2020) A6, [arXiv:1807.06209].
  * [47] E. Krause, X. Fang, S. Pandey, L. F. Secco, O. Alves, H. Huang, J. Blazek, J. Prat, J. Zuntz, T. F. Eifler, N. MacCrann, J. DeRose, M. Crocce, A. Porredon, B. Jain, M. A. Troxel, S. Dodelson, D. Huterer, A. R. Liddle, C. D. Leonard, A. Amon, A. Chen, J. Elvin-Poole, A. Ferté, J. Muir, Y. Park, S. Samuroff, A. Brandao-Souza, N. Weaverdyck, G. Zacharegkas, R. Rosenfeld, A. Campos, P. Chintalapati, A. Choi, E. Di Valentino, C. Doux, K. Herner, P. Lemos, J. Mena-Fernández, Y. Omori, M. Paterno, M. Rodriguez-Monroy, P. Rogozenski, R. P. Rollins, A. Troja, I. Tutusaus, R. H. Wechsler, T. M. C. Abbott, M. Aguena, S. Allam, F. Andrade-Oliveira, J. Annis, D. Bacon, E. Baxter, K. Bechtol, G. M. Bernstein, D. Brooks, E. Buckley-Geer, D. L. Burke, A. Carnero Rosell, M. Carrasco Kind, J. Carretero, F. J. Castander, R. Cawthon, C. Chang, M. Costanzi, L. N. da Costa, M. E. S. Pereira, J. De Vicente, S. Desai, H. T. Diehl, P. Doel, S. Everett, A. E. Evrard, I. Ferrero, B. Flaugher, P. Fosalba, J. Frieman, J. García-Bellido, E. Gaztanaga, D. W. Gerdes, T. Giannantonio, D. Gruen, R. A. Gruendl, J. Gschwend, G. Gutierrez, W. G. Hartley, S. R. Hinton, D. L. Hollowood, K. Honscheid, B. Hoyle, E. M. Huff, D. J. James, K. Kuehn, N. Kuropatkin, O. Lahav, M. Lima, M. A. G. Maia, J. L. Marshall, P. Martini, P. Melchior, F. Menanteau, R. Miquel, J. J. Mohr, R. Morgan, J. Myles, A. Palmese, F. Paz-Chinchón, D. Petravick, A. Pieres, A. A. Plazas Malagón, E. Sanchez, V. Scarpine, M. Schubnell, S. Serrano, I. Sevilla-Noarbe, M. Smith, M. Soares-Santos, E. Suchyta, G. Tarle, D. Thomas, C. To, T. N. Varga, and J. Weller, Dark Energy Survey Year 3 Results: Multi-Probe Modeling Strategy and Validation, arXiv e-prints (May, 2021) arXiv:2105.13548, [arXiv:2105.13548].
  * [48] J. Torrado and A. Lewis, Cobaya: code for Bayesian analysis of hierarchical physical models, JCAP 2021 (May, 2021) 057, [arXiv:2005.05290].
  * [49] A. Lewis, GetDist: a Python package for analysing Monte Carlo samples, arXiv:1910.13970.
  * [50] A. Lewis, A. Challinor, and A. Lasenby, Efficient computation of CMB anisotropies in closed FRW models, ApJ 538 (2000) 473–476, [astro-ph/9911177].