
Future Baryon Acoustic Oscillation surveys aim at observing galaxy clustering
over a wide range of redshift and galaxy populations at great precision,
reaching tenths of a percent, in order to detect any deviation of dark energy
from the $\LCDM$ model. We utilize a set of paired quasi\Nb\, FastPM
simulations that were designed to mitigate the sample variance effect on the
BAO feature and evaluated the BAO systematics as precisely as $\sim 0.01\%$. We
report anisotropic BAO scale shifts before and after density field
reconstruction in the presence of redshiftspace distortions over a wide range
of redshift, galaxy/halo biases, and shot noise levels. We test different
reconstruction schemes and different smoothing filter scales, and introduce
physicallymotivated BAO fitting models. For the first time, we derive a
Galileaninvariant infrared resummed model for halos in real and redshift
space. We test these models from the perspective of robust BAO measurements and
nonBAO information such as growth rate and nonlinear bias. We find that
prereconstruction BAO scale has moderate fittingmodel dependence at the level
of $0.1\%0.2\%$ for matter while the dependence is substantially reduced to
less than $0.07\%$ for halos. We find that postreconstruction BAO shifts are
generally reduced to below $0.1\%$ in the presence of galaxy/halo bias and show
much smaller fitting model dependence. Different reconstruction conventions can
potentially make a much larger difference on the lineofsight BAO scale, upto
$0.3\%$. Meanwhile, the precision (error) of the BAO measurements is quite
consistent regardless of the choice of the fitting model or reconstruction
convention.

As a major source of cosmological information, galaxy clustering is
susceptible to longwavelength density and tidal fluctuations. These long modes
modulate the growth and expansion rate of local structures, shifting them in
both amplitude and scale. These effects are often named the growth and dilation
effects, respectively. In particular the dilation shifts the baryon acoustic
oscillation (BAO) peak and breaks the assumption of the AlcockPaczynski (AP)
test. This cannot be removed with reconstruction techniques because the effect
originates from long modes outside the survey. In redshift space, the long
modes generate a largescale radial peculiar velocity that affects the
redshiftspace distortion (RSD) signal. We compute the redshiftspace response
functions of the galaxy power spectrum to long density and tidal modes at
leading order in perturbation theory, including both the growth and dilation
terms. We validate these response functions against measurements from simulated
galaxy mock catalogs. As one application, long density and tidal modes beyond
the scale of a survey correlate various observables leading to an excess error
known as the supersample covariance, and thus weaken their constraining power.
We quantify the supersample effect on BAO, AP, and RSD measurements, and study
its impact on current and future surveys.

Motivated by recent developments in perturbative calculations of the
nonlinear evolution of largescale structure, we present an iterative algorithm
to reconstruct the initial conditions in a given volume starting from the dark
matter distribution in real space. In our algorithm, objects are first moved
back iteratively along estimated potential gradients, with a progressively
reduced smoothing scale, until a nearly uniform catalog is obtained. The linear
initial density is then estimated as the divergence of the cumulative
displacement, with an optional secondorder correction. This algorithm should
undo nonlinear effects up to oneloop order, including the higherorder
infrared resummation piece. We test the method using dark matter simulations in
real space. At redshift $z=0$, we find that after eight iterations the
reconstructed density is more than $95\%$ correlated with the initial density
at $k\le 0.35\; h\mathrm{Mpc}^{1}$. The reconstruction also reduces the power
in the difference between reconstructed and initial fields by more than 2
orders of magnitude at $k\le 0.2\; h\mathrm{Mpc}^{1}$, and it extends the
range of scales where the full broadband shape of the power spectrum matches
linear theory by a factor of 23. As a specific application, we consider
measurements of the baryonic acoustic oscillation (BAO) scale that can be
improved by reducing the degradation effects of largescale flows. In our
idealized dark matter simulations, the method improves the BAO signaltonoise
ratio by a factor of 2.7 at $z=0$ and by a factor of 2.5 at $z=0.6$, improving
standard BAO reconstruction by $70\%$ at $z=0$ and $30\%$ at $z=0.6$, and
matching the optimal BAO signal and signaltonoise ratio of the linear density
in the same volume. For BAO, the iterative nature of the reconstruction is the
most important aspect.

We study the matter bispectrum of largescale structure by comparing the
predictions of different perturbative and phenomenological models with the full
threedimensional bispectrum from $N$body simulations estimated using modal
methods. We show that among the perturbative approaches, effective field theory
succeeds in extending the range of validity furthest on intermediate scales, at
the cost of free additional parameters. By studying the halo model, we show
that although it is satisfactory in the deeply nonlinear regime, it predicts a
deficit of power on intermediate scales, worsening at redshifts $z>0$. By
comparison with the $N$body bispectrum on those scales, we show that there is
a significant squeezed component underestimated in the halo model. On the basis
of these results, we propose a new threeshape model, based on the treelevel,
squeezed and constant bispectrum shapes we identified in the halo model; after
calibration this fits the simulations on all scales and redshifts of interest.
We extend this model further to primordial nonGaussianity of the local and
equilateral types by showing that the same shapes can be used to describe the
additional nonGaussian component in the matter bispectrum. This method
provides a HALOFITlike prototype of the bispectrum that could be used to
describe and test parameter dependencies and should be relevant for the
bispectrum of weak gravitational lensing and wider applications.

CMB and lensing reconstruction power spectra are powerful probes of
cosmology. However they are correlated, since the CMB power spectra are lensed
and the lensing reconstruction is constructed using CMB multipoles. We perform
a full analysis of the auto and crosscovariances, including polarization
power spectra and minimum variance lensing estimators, and compare with
simulations of idealized future CMBS4 observations. Covariances sourced by
fluctuations in the unlensed CMB and instrumental noise can largely be removed
by using a realizationdependent subtraction of lensing reconstruction noise,
leaving a relatively simple covariance model that is dominated by
lensinginduced terms and well described by a small number of principal
components. The correlations between the CMB and lensing power spectra will be
detectable at the level of $\sim 5\sigma$ for a CMBS4 mission, and neglecting
those could underestimate some parameter error bars by several tens of percent.
However we found that the inclusion of external priors or data sets to estimate
parameter error bars can make the impact of the correlations almost negligible.

Modeling the largescale structure of the universe on nonlinear scales has
the potential to substantially increase the science return of upcoming surveys
by increasing the number of modes available for model comparisons. One way to
achieve this is to model nonlinear scales perturbatively. Unfortunately, this
involves highdimensional loop integrals that are cumbersome to evaluate.
Trying to simplify this, we show how twoloop (nexttonexttoleading order)
corrections to the density power spectrum can be reduced to lowdimensional,
radial integrals. Many of those can be evaluated with a onedimensional Fast
Fourier Transform, which is significantly faster than the fivedimensional
MonteCarlo integrals that are needed otherwise. The general idea of this
FFTPT method is to switch between Fourier and position space to avoid
convolutions and integrate over orientations, leaving only radial integrals.
This reformulation is independent of the underlying shape of the initial linear
density power spectrum and should easily accommodate features such as those
from baryonic acoustic oscillations. We also discuss how to account for halo
bias and redshift space distortions.

Forming a three dimensional view of the Universe is a longstanding goal of
astronomical observations, and one that becomes increasingly difficult at high
redshift. In this paper we discuss how tomography of the intergalactic medium
(IGM) at $z\simeq 2.5$ can be used to estimate the redshifts of massive
galaxies in a large volume of the Universe based on spectra of galaxies in
their background. Our method is based on the fact that hierarchical structure
formation leads to a strong dependence of the halo density on largescale
environment. A map of the latter can thus be used to refine our knowledge of
the redshifts of halos and the galaxies and AGN which they host. We show that
tomographic maps of the IGM at a resolution of $2.5\,h^{1}$Mpc can determine
the redshifts of more than 90 per cent of massive galaxies with redshift
uncertainty $\Delta z/(1+z)=0.01$. Higher resolution maps allow such redshift
estimation for lower mass galaxies and halos.

Given the importance of future large scale structure surveys for delivering
new cosmological information, it is crucial to reliably predict their
observables. The Effective Field Theory of Large Scale Structures (EFTofLSS)
provides a manifestly convergent perturbative scheme to compute the clustering
of dark matter in the weakly nonlinear regime in an expansion in $k/k_{\rm
NL}$, where $k$ is the wavenumber of interest and $k_{\rm NL}$ is the
wavenumber associated to the nonlinear scale. It has been recently shown that
the EFTofLSS matches to $1\%$ level the dark matter power spectrum at redshift
zero up to $k\simeq 0.3 h\,$Mpc$^{1}$ and $k\simeq 0.6 h\,$Mpc$^{1}$ at one
and two loops respectively, using only one counterterm that is fit to data.
Similar results have been obtained for the momentum power spectrum at one loop.
This is a remarkable improvement with respect to former analytical techniques.
Here we study the prediction for the equaltime dark matter bispectrum at one
loop. We find that at this order it is sufficient to consider the same
counterterm that was measured in the power spectrum. Without any remaining free
parameter, and in a cosmology for which $k_{\rm NL}$ is smaller than in the
previously considered cases ($\sigma_8=0.9$), we find that the prediction from
the EFTofLSS agrees very well with $N$body simulations up to $k\simeq 0.25
h\,$Mpc$^{1}$, given the accuracy of the measurements, which is of order a few
percent at the highest $k$'s of interest. While the fit is very good on average
up to $k\simeq 0.25 h\,$Mpc$^{1}$, the fit performs slightly worse on
equilateral configurations, in agreement with expectations that for a given
maximum $k$, equilateral triangles are the most nonlinear.

The usual fluid equations describing the largescale evolution of mass
density in the universe can be written as local in the density, velocity
divergence, and velocity potential fields. As a result, the perturbative
expansion in small density fluctuations, usually written in terms of
convolutions in Fourier space, can be written as a series of products of these
fields evaluated at the same location in configuration space. Based on this, we
establish a new method to numerically evaluate the 1loop power spectrum (i.e.,
Fourier transform of the 2point correlation function) with onedimensional
Fast Fourier Transforms. This is exact and a few orders of magnitude faster
than previously used numerical approaches. Numerical results of the new method
are in excellent agreement with the standard quadrature integration method.
This fast model evaluation can in principle be extended to higher loop order
where existing codes become painfully slow. Our approach follows by writing
higher order corrections to the 2point correlation function as, e.g., the
correlation between two secondorder fields or the correlation between a linear
and a thirdorder field. These are then decomposed into products of
correlations of linear fields and derivatives of linear fields. The method can
also be viewed as evaluating threedimensional Fourier space convolutions using
products in configuration space, which may also be useful in other contexts
where similar integrals appear.

The rapidly improving precision of measurements of gravitational lensing of
the Cosmic Microwave Background (CMB) also requires a corresponding increase in
the precision of theoretical modeling. A commonly made approximation is to
model the CMB deflection angle or lensing potential as a Gaussian random field.
In this paper, however, we analytically quantify the influence of the
nonGaussianity of largescale structure lenses, arising from nonlinear
structure formation, on CMB lensing measurements. In particular, evaluating the
impact of the nonzero bispectrum of largescale structure on the relevant CMB
fourpoint correlation functions, we find that there is a bias to estimates of
the CMB lensing power spectrum. For temperaturebased lensing reconstruction
with CMB StageIII and StageIV experiments, we find that this lensing power
spectrum bias is negative and is of order one percent of the signal. This
corresponds to a shift of multiple standard deviations for these upcoming
experiments. We caution, however, that our numerical calculation only evaluates
two of the largest bias terms and thus only provides an approximate estimate of
the full bias. We conclude that further investigation into lensing biases from
nonlinear structure formation is required and that these biases should be
accounted for in future lensing analyses.

We study the matter bispectrum of the largescale structure by comparing
different perturbative and phenomenological models with measurements from
$N$body simulations obtained with a modal bispectrum estimator. Using shape
and amplitude correlators, we directly compare simulated data with theoretical
models over the full threedimensional domain of the bispectrum, for different
redshifts and scales. We review and investigate the main perturbative methods
in the literature that predict the oneloop bispectrum: standard perturbation
theory, effective field theory, resummed Lagrangian and renormalised
perturbation theory, calculating the latter also at two loops for some triangle
configurations. We find that effective field theory (EFT) succeeds in extending
the range of validity furthest into the mildly nonlinear regime, albeit at the
price of free extra parameters requiring calibration on simulations. For the
more phenomenological halo model, we confirm that despite its validity in the
deeply nonlinear regime it has a deficit of power on intermediate scales, which
worsens at higher redshifts; this issue is ameliorated, but not solved, by
combined haloperturbative models. We show from simulations that in this
transition region there is a strong squeezed bispectrum component that is
significantly underestimated in the halo model at earlier redshifts. We thus
propose a phenomenological method for alleviating this deficit, which we
develop into a simple phenomenological "threeshape" benchmark model based on
the three fundamental shapes we have obtained from studying the halo model.
When calibrated on the simulations, this threeshape benchmark model accurately
describes the bispectrum on all scales and redshifts considered, providing a
prototype bispectrum HALOFITlike methodology that could be used to describe
and test parameter dependencies.

As galaxy surveys begin to measure the imprint of baryonic acoustic
oscillations (BAO) on largescale structure at the subpercent level,
reconstruction techniques that reduce the contamination from nonlinear
clustering become increasingly important. Inverting the nonlinear continuity
equation, we propose an Eulerian growthshift reconstruction algorithm that
does not require the displacement of any objects, which is needed for the
standard Lagrangian BAO reconstruction algorithm. In realspace DMonly
simulations the algorithm yields 95% of the BAO signaltonoise obtained from
standard reconstruction. The reconstructed power spectrum is obtained by adding
specific simple 3 and 4point statistics to the prereconstruction power
spectrum, making it very transparent how additional BAO information from
higherpoint statistics is included in the power spectrum through the
reconstruction process. Analytical models of the reconstructed density for the
two algorithms agree at second order. Based on similar modeling efforts, we
introduce four additional reconstruction algorithms and discuss their
performance.

Clustering of largescale structure provides significant cosmological
information through the power spectrum of density perturbations. Additional
information can be gained from higherorder statistics like the bispectrum,
especially to break the degeneracy between the linear halo bias $b_1$ and the
amplitude of fluctuations $\sigma_8$. We propose new simple, computationally
inexpensive bispectrum statistics that are near optimal for the specific
applications like bias determination. Corresponding to the Legendre
decomposition of nonlinear halo bias and gravitational coupling at second
order, these statistics are given by the crossspectra of the density with
three quadratic fields: the squared density, a tidal term, and a shift term.
For halos and galaxies the first two have associated nonlinear bias terms $b_2$
and $b_{s^2}$, respectively, while the shift term has none in the absence of
velocity bias (valid in the $k \rightarrow 0$ limit). Thus the linear bias
$b_1$ is best determined by the shift crossspectrum, while the squared density
and tidal crossspectra mostly tighten constraints on $b_2$ and $b_{s^2}$ once
$b_1$ is known. Since the form of the crossspectra is derived from optimal
maximumlikelihood estimation, they contain the full bispectrum information on
bias parameters. Perturbative analytical predictions for their expectation
values and covariances agree with simulations on large scales, $k\lesssim
0.09h/\mathrm{Mpc}$ at $z=0.55$ with Gaussian $R=20h^{1}\mathrm{Mpc}$
smoothing, for mattermattermatter, and mattermatterhalo combinations. For
halohalohalo crossspectra the model also needs to include corrections to the
Poisson stochasticity.

As confusion with lensing Bmodes begins to limit experiments that search for
primordial Bmode polarization, robust methods for delensing the CMB
polarization sky are becoming increasingly important. We investigate in detail
the possibility of delensing the CMB with the cosmic infrared background (CIB),
emission from dusty starforming galaxies that is an excellent tracer of the
CMB lensing signal, in order to improve constraints on the tensortoscalar
ratio $r$. We find that the maps of the CIB, such as current Planck satellite
maps at 545 GHz, can be used to remove more than half of the lensing Bmode
power. Calculating optimal combinations of different largescalestructure
tracers for delensing, we find that coadding CIB data and external
arcminuteresolution CMB lensing reconstruction can lead to significant
additional improvements in delensing performance. We investigate whether
measurement uncertainty in the CIB spectra will degrade the delensing
performance if no model of the CIB spectra is assumed, and instead the CIB
spectra are marginalized over, when constraining $r$. We find that such
uncertainty does not significantly affect Bmode surveys smaller than a few
thousand degrees. Even for larger surveys it causes only a moderate reduction
in CIB delensing performance, especially if the surveys have high (arcminute)
resolution, which allows selfcalibration of the delensing procedure. Though
further work on the impact of foreground residuals is required, our overall
conclusions for delensing with current CIB data are optimistic: this delensing
method can tighten constraints on $r$ by a factor up to $\approx2.2$, and by a
factor up to $\approx4$ when combined with external $\approx 3 \mu$Karcmin
lensing reconstruction, without requiring the modeling of CIB properties. CIB
delensing is thus a promising method for the upcoming generation of CMB
polarization surveys.