# Primordial feature constraints from BOSS $+$ eBOSS

Thiago Mergulhão  Florian Beutler  John A. Peacock

###### Abstract

Understanding the universe in its pristine epoch is crucial in order to obtain
a concise comprehension of the late-time universe. Although current data in
cosmology are compatible with Gaussian primordial perturbations whose power
spectrum follows a nearly scale-invariant power law, this need not be the case
when a fundamental theoretical construction is assumed. These extended models
lead to sharp features in the primordial power spectrum, breaking its scale
invariance. In this work, we obtain combined constraints on four primordial
feature models by using the final data release of the BOSS galaxies and eBOSS
quasars. By pushing towards the fundamental mode of these surveys and using
the larger eBOSS volume, we were able to extend the feature parameter space
(i.e. the feature frequency $\omega$) by a factor of four compared to previous
analyses using BOSS. While we did not detect any significant features,
previous work showed that next-generation galaxy surveys such as DESI will
improve the sensitivity to features by a factor of 7, and will also extend the
parameter space by a factor of 2.5.

##  1 Introduction

Over the last few decades, we have witnessed major steps in cosmological
physics toward a better understanding of the universe and its nature. Projects
such as the Planck collaboration [1, 2], which probed the cosmic microwave
background (CMB) in an unprecedented way, and the Baryon Oscillation
Spectroscopic Survey (BOSS), which constructed a three-dimensional map of the
late universe distribution of galaxies, imposed strong constraints on the
standard cosmological model. However, there are still questions that need to
be addressed, many of them regarding the early universe and its dynamics – the
same desideratum present in particle physics to look for a more fundamental
description of nature is also present in cosmology, and in some sense, they
are fairly connected [3].

The primordial power spectrum of curvature perturbations is the seed for
structure formation, and its form depends on the physical processes that
operate at early times. In the standard model of cosmology, these
perturbations are: adiabatic; their statistical distribution is highly
Gaussian; and their power spectrum is a nearly scale-invariant power-law. Two
of the six free parameters in the standard model are associated with this: the
amplitude of the scalar perturbations, $A_{\rm s}$, and the spectral index
$n_{\rm s}$. Although this simple model allows us to explain a vast number of
observables, there are many reasons why this model by itself is not natural,
and there are even more ways to construct other viable scenarios. Among these,
we give the following examples111See [4] for a recent review about these and
other primordial feature models.:

  * •

Embedding inflation in grand unified theories can result in a series of phase
transitions [5] that can induce steps in the potential [6, 7], a mechanism
that gives a sudden burst to the inflaton’s kinetic energy. The impact of
these steps has been extensively analysed in the literature and they usually
cause localized oscillatory signals in the power spectrum [8, 9, 10, 11].

  * •

The inflaton can be only one degree of freedom in a multi-field scenario.
Considering this framework, the inflaton can go through a region in field
space where some of the massive degrees of freedom to which it is coupled
suddenly become massless, causing the production of heavy field quanta [12,
13, 14, 15, 16, 17]. Another scenario is when any of these other hidden fields
perform a sudden turn in their configuration space, which produces oscillating
features on top of the usual power-law prediction [18, 19, 20, 21, 22].

  * •

Inflation is known to magnify quantum fluctuations generated in the indefinite
past by many orders of magnitude, making it relevant for later cosmic
evolution. The initial wavelength of these fluctuations can be as small as the
Planck scale, and transPlanckian physics can modify the primordial power
spectrum, including global features [23, 24, 25, 26, 27, 28].

  * •

Logarithmically oscillating features can be generated in both natural [29] and
axion monodromy inflation [30]. Given the importance of axions in cosmology
[31] and the simplicity of their associated corrections to the primordial
power spectrum, constraining their impact with current and forthcoming data
may be vital for studying their existence.

The physical mechanism behind these primordial features (PF) is quite broad
and intrinsically connected to our fundamental descriptions of nature. Hence,
any detection of an additional inflationary signal lurking in data could
drastically change our picture of the early phases of the universe.
Consequently, many groups have been trying to understand the impact of these
primordial features on data and to develop pipelines to detect them. There is
already a plethora of analyses using WMAP [32, 33, 34, 35, 36, 37, 38, 39] and
Planck [8, 40, 2, 41, 42, 43, 44] datasets. Forecasts for future CMB surveys
are present in recent works [45, 46, 47]. Gravitational waves are one of the
main byproducts of inflation; the presence of PF could also affect the
stochastic gravitational wave background [48, 49, 50, 51, 52, 53], and [54]
showed a forecast for LISA. It is important to mention that even if inflation
is not assumed to be true, PF is one of the main signals that could
potentially help us to distinguish inflation from other scenarios [26, 55, 56,
44].

Additional primordial physics can affect the clustering of tracers, which
makes large-scale structure (LSS) a potential source of information on PF to
complement other datasets. These signals, if present, are susceptible to non-
linear effects arising from late-time structure formation. A well-known
example of this effect is the degradation of the Baryon Acoustic Oscillation
(BAO) peak due to non-linear clustering [57]. Nevertheless, much effort has
been made to address this problem. In [58, 59, 60, 61, 62] the impact of non-
linearities on PF was studied, providing analytical predictions using
perturbation theory. These results were tested with N-body simulations,
proving to be precise enough for analysis using current and future datasets.
In addition, by using perturbation theory one can partially undo the non-
linear effects by performing a field reconstruction [63, 64, 65]. Altogether,
this allowed the development of a robust analysis of PF with LSS [60]; this
showed that the final BOSS data release improved on the Planck constraints on
both linear and logarithmic features over a considerable volume in parameter
space.

The importance of LSS constraints on primordial features and the fast
progression of stage-IV spectroscopic surveys like DESI and Euclid motivate us
to start developing an in-depth search for PF using large-scale structure
data. In this work, we push further the analysis performed in [60] and obtain
combined PF credible intervals from the final BOSS galaxies and eBOSS quasars.
We show that the high number density of the former along with the greater
volume of the latter are complementary, allowing us to obtain credible
intervals for primordial features over a region in parameter space never
explored before. This work is structured in the following way: In Section 2 we
explain how the feature signals can be probed by a broadband analysis of the
galaxy power spectrum. In Section 3 we show all the primordial feature signals
that will be considered in this work and discuss their impact on the
distribution of galaxies. In Section 4 we explain the data we used. In Section
5 we scrutinize our methodology along with the statistical tools used to check
the significance of our results, which is presented in Section 6.

##  2 Galaxy clustering template

The presence of primordial features changes the initial power spectrum of
curvature perturbations. In this section, we show how it is possible to
propagate the primordial feature signal into the evolved galaxy field at
linear order.

###  2.1 Propagating the primordial feature to the evolved matter field

The primordial feature signals cannot be mimicked by any specific combination
of the cosmological parameters in the standard model. For this reason, we
follow [60] and assume that as a first approximation, one can simply
marginalize over the broadband signal in the galaxy power spectrum and then
look for residual signals. To see how PF propagate to the evolved matter field
we consider the linear evolution of the initial conditions through the
transfer function, $T(k)$:

$P(z,k)=k^{4}T(k)^{2}D(z)^{2}P_{\mathcal{R}}(k)\;,$ (2.1)

where $D(z)$ is the linear growth factor and $P_{\mathcal{R}}$ is the
primordial power spectrum. In the simplest models of inflation, this power
spectrum is given by

$P_{\mathcal{R},0}(k)=\frac{2\pi^{2}}{k^{3}}\mathcal{P}_{\mathcal{R},0}(k)=\frac{2\pi^{2}A_{\mathrm{s}}}{k^{3}}\left(\frac{k}{k_{\star}}\right)^{n_{\mathrm{s}}-1}\;.$
(2.2)

The contributions of primordial features, $\delta P_{\mathcal{R}}$, show up as
corrections to Eq. (2.2) in the form

$P_{\mathcal{R}}(k)=P_{\mathcal{R},0}(k)\left[1+\delta
P_{\mathcal{R}}(k)\right]\;.$ (2.3)

By plugging Eq. (2.3) into Eq. (2.1) it is clear that

$\displaystyle P(z,k)$ $\displaystyle=$ $\displaystyle P_{\rm
lin}(z,k)\left[1+\delta P_{\mathcal{R}}(k)\right]$ (2.4) $\displaystyle=$
$\displaystyle P^{\rm nw}_{\rm lin}(z,k)+P^{\rm w}_{\rm linBAO}(z,k)+\delta
P_{\mathcal{R}}P^{\rm nw}_{\rm lin}+P^{\rm w}_{\rm lin,BAO}\delta
P_{\mathcal{R}}(k)\;,$

where we split the usual linear matter field into a smooth part and the BAO
contributions, $P_{\rm lin}=P^{\rm nw}_{\rm lin}+P^{\rm w}_{\rm lin,BAO}$. Eq.
(2.4) is our final result for the evolved matter power spectrum with PF.

###  2.2 Galaxy clustering with primordial features

The galaxy power spectrum is introduced in the following way. In equation Eq.
(2.4) we factor out a linear smooth power spectrum, $P^{\rm nw}_{\rm lin}$,
that contains both BAO and PF oscillations, $\delta P_{\mathcal{R}}$:

$\tilde{P}(z,k)=P^{\rm nw}_{\rm lin}\left(1+\mathcal{O}_{\rm lin}+\delta
P_{\mathcal{R}}+\mathcal{O}_{\rm lin}\delta P_{\mathcal{R}}\right)\;,$ (2.5)

where $\mathcal{O}_{\rm lin}(k)\equiv P^{\rm w}_{\rm lin,BAO}/P^{\rm nw}_{\rm
lin}$. We included a tilde in the expression above to emphasize that this is
not the final result for the evolved power spectrum – non-linear displacements
will damp the oscillations and need to be taken into account (see §2.3). Eq.
(2.5) is the result for the matter field, and we need to connect it with our
observable, the tracer field. Since we are interested only in modulations
around the smooth power spectrum, we employ a phenomenological model instead
of the usual bias expansion in perturbation theory: we assume that the
connection between the matter and the galaxy field is given by broadband
information encapsulated in a linear bias, $B$, a spectral distortion due to
the non-linear velocity field, $F\left(k,\Sigma_{s}\right)$, and a series of
broadband parameters, $A(k)$:

$P^{\rm nw}_{\rm
g}(k)=B^{2}P^{\mathrm{nw}}(k)F\left(k,\Sigma_{s}\right)+A(k)\;,$ (2.6)

with

$\displaystyle A(k)$ $\displaystyle=$
$\displaystyle\frac{a_{1}}{k^{3}}+\frac{a_{2}}{k^{2}}+\frac{a_{3}}{k}+a_{4}+a_{5}k^{2}\;,$
(2.7) $\displaystyle F\left(k,\Sigma_{s}\right)$ $\displaystyle=$
$\displaystyle\frac{1}{\left(1+k^{2}\Sigma_{s}^{2}/2\right)^{2}}\;.$ (2.8)

The forms for $A(k)$ and $F(k,\Sigma_{s})$ are well motivated and commonly
used in clustering analysis [66, 67, 68, 69, 70, 71, 72].

###  2.3 The non-linear damping of oscillatory features: infrared resummation

Although Eq. (2.5) considers the impact of both spectral distortion due to the
velocity field and a set of broadband parameters to account for the broadband
non-linear clustering that we marginalise over, there is still one effect that
needs to be considered: the displacement of the galaxies due to long-
wavelength perturbations. In Eulerian perturbation theory both the overdensity
and velocity fields are considered to be small. However, since the usual time
scale of the system is the Hubble time, $H(a)^{-1}$, the galaxies end up
having significant displacements. In Eulerian perturbation theory, this
displacement field is considered to be small, but this is not the case. Thus,
the perturbations need to be resummed in order to have a consistent result
[73]. In the case of a linear primordial feature, the effect of the damping
will be the same as in the BAO case, and it transforms Eq. (2.5) into [58,
59]:

$P_{m}(k)=P^{\rm nw}(k)\left[1+\left(\mathcal{O}_{\rm lin}+\delta
P_{\mathcal{R}}+\mathcal{O}_{\rm lin}\delta
P_{\mathcal{R}}\right)e^{-k^{2}\Sigma_{\rm lin}^{2}}\right]\;,$ (2.9)

where

$\Sigma_{\operatorname{lin}}^{2}\equiv\frac{1}{6\pi^{2}}\int
dq\,P^{\mathrm{nw}}(q)\left[1-j_{0}\left(q\omega_{\operatorname{lin}}\right)+2j_{2}\left(q\omega_{\operatorname{lin}}\right)\right]\;.$
(2.10)

Eq. (2.9) is the final result for the matter field after the large-scale bulk
flow is considered. It can easily be extended to galaxy clustering:

$\displaystyle P_{g}(k)=P_{g}^{\mathrm{nw}}(k)\left\\{1+[\mathcal{O}_{\rm
lin}(k)+\delta P_{\mathcal{R}}^{X}(k)+\mathcal{O}_{\rm lin}(k)\delta
P_{\mathcal{R}}^{X}(k)]e^{-k^{2}\Sigma_{\mathrm{nl}}^{2}/2}\right\\}\;,$
(2.11)

where the primordial feature corrections are represented by $\delta
P_{\mathcal{R}}$. Before continuing, we emphasize that although the damping
term included above, $e^{-k^{2}\Sigma_{\mathrm{nl}}^{2}/2}$, is well motivated
for linear features, it is expected to have a scale-dependence for e.g. the
logarithmic feature. In [60], the authors showed that as a first
approximation, it is possible to use a scale-independent damping model for the
logarithmic feature. Here we follow [60] and approximate the non-linear
damping by using the scale-independent parameter $\Sigma_{\rm
nl}\approx\Sigma_{\rm lin}$. In [61] the authors showed how the non-linear
damping can be approximated analytically in perturbation theory for general
sharp features, and we will employ their results in future work. The simple
approximation that we make here has been tested in [62] using data from an
N-body simulation: they concluded that the template in Eq. (2.11) described
their data reasonably well. This work considered dark matter only and it would
be good to repeat it with a variety of galaxy-formation recipes; but our
expression in Eq. (2.11) is essentially the one that is commonly used to
describe damping of BAO features, so we do not expect its applicability to
depend on such details.

Another point that needs to be addressed is the inclusion of the isotropic BAO
scaling parameter, $\alpha$. In the standard analysis in Fourier space, this
parameter specifies the frequency or wavelength of the BAO wiggles and shows
up in Eq. (2.11) as:

$P_{g}(k)=P_{g}^{\mathrm{nw}}(k)\\{1+[\mathcal{O}_{\rm lin}(k/\alpha)+\delta
P_{\zeta}^{X}(k)+\mathcal{O}_{\rm lin}(k/\alpha)\delta
P_{\zeta}^{X}(k)]e^{-k^{2}\Sigma_{\mathrm{nl}}^{2}/2}\\}\;.$ (2.12)

To summarize, in addition to the PF-free parameters, we also need to fit the
following ten parameters per dataset:

$\\{B_{\rm NGC},B_{\rm
SGC},a_{1},a_{2},a_{3},a_{4},a_{5},\Sigma_{s},\Sigma_{\rm lin},\alpha\\}\;.$
(2.13)

Notice that we included two linear amplitudes, $B_{\rm NGC}$ and $B_{\rm
SGC}$. The reason is that for a given dataset, we have two parts, the north
galactic cap (NGC) and the south galactic cap (SGC). Therefore we allow them
to have different linear amplitudes, but the other parameters are shared (see
Section 4). This model follows the standard BOSS BAO analysis [68].

##  3 Primordial feature signals

We now present more detail on the specific PF models that we will seek to
constrain using LSS data.

###  3.1 Linear features

If some background inflationary quantities have a sudden change, the inflaton
field can have a feature that is spread over its entire power spectrum. A
general template for this signal is

$\delta\mathcal{P}_{\mathcal{R}}=A_{\rm lin}\sin(\omega_{\rm lin}k+\phi)\;,$
(3.1)

where $A_{\rm lin}$ is the feature amplitude, $\phi$ is some arbitrary phase
and $\omega_{\rm lin}$ is the frequency. In configuration space, a linear
feature appears as a second peak on the correlation function at a separation
of $s=\omega_{\rm lin}$. This is shown in Figure 1. This kind of signal was
already scanned for in both CMB [2] and LSS data [60], where the authors
showed that existing LSS datasets like BOSS can surpass Planck in parts of the
parameter space. However, the previous LSS analysis did not explore the entire
parameter space that can be constrained with the available data. The maximum
frequency range one can probe depends on the fundamental mode, $k_{f}$, which
is the minimum wave number available in a survey: for a survey with arbitrary
geometry and volume $V$, it can be approximated by $k_{f}=2\pi/V^{1/3}$. The
importance of this scale arises because the measured power spectrum is
convolved with the survey window function, which automatically damps any
possible structure in $P(k)$ that is much finer than $k_{f}$. Here we analyze
both the BOSS and eBOSS datasets using a binning in $k$ close to their
fundamental mode, thus allowing us to scan for PF with the highest frequencies
that can be probed given the limit set by the survey window. We will further
elucidate this point in Section 4.3. For linear features, we decided to use
the parameter range specified in Table 1. We emphasize that Eq. (3.1) is a
phenomenological template: although being simple, its global oscillation in
$k$ is theoretically unrealistic. Nevertheless, any relevant detection of Eq.
(3.1) will ultimately indicate the presence of a sinusoidal running in the
data. This conclusion will motivate the scan of other more realistic models
with similar running, so Eq. (3.1) is definitely a robust starting point for
primordial features analysis. We consider two of those realistic models in
§3.3 and §3.4.

Figure 1: left: The primordial curvature power spectrum in the vanilla model
of inflation (black dashed line) and in the linear feature case (solid blue
line). For this plot, we used $A_{\rm lin}=0.05$, $\omega_{\rm
lin}=250\,\mathrm{Mpc}$, $\phi=0$ and the BOSS high-$z$ window function;
right: An illustration of how the linear feature shows up in the evolved
galaxy field.

###  3.2 Logarithmic features

Figure 2: left: The primordial curvature power spectrum in the vanilla model
of inflation (black dashed line) and in the logarithmic feature case (solid
orange line). For this plot, we used $A_{\rm log}=0.05$, $\omega_{\rm
log}=35$, $\phi=0$ and the BOSS high-$z$ window function; right: An
illustration of how the logarithmic feature shows up in the evolved galaxy
field.

Axions are promising candidates for the inflaton. Contrary to the usual scalar
field, axions have a continuous shift symmetry that protects their potential
against radiative corrections, so they are safe from the eta-problem [74], and
the slow-roll conditions can be satisfied naturally. A known example of
inflation being driven by axions is axion monodromy inflation [4, 37, 75]. In
that model, the axion potential is:

$V(\phi)=\mu^{3}\phi-\Lambda^{4}\cos\left(\theta-\frac{\phi}{f}\right)\;,$
(3.2)

where $\mu$ is a mass scale. The presence of a sinusoidal modulation in Eq.
(3.2) includes a new mechanism in the theory in which primordial non-
Gaussianity can be generated: the oscillating perturbations can be in
resonance with the oscillating inflaton background, which can lead to what is
called Resonant non-Gaussianity [76, 77, 78, 79]. This parametric resonance
affects the primordial power spectrum, generating features [4, 37, 77]:

$\delta\mathcal{P}_{\mathcal{R}}(k)=\kappa\cos\left[\frac{\phi_{*}}{f}-\frac{\ln\left(k/k_{*}\right)}{\widetilde{f}}\right]\;,$
(3.3)

where $\phi_{*}$ is the background value when $k_{*}$ exits the horizon,
$\kappa=3\left(\Lambda^{4}/\mu^{3}f\right)(2\pi\tilde{f})^{1/2}$ and
$\smash{\tilde{f}=f\phi_{*}/M_{\rm Pl}^{2}}$. Given the theoretical importance
of axions in both cosmology and particle physics, constraining these
logarithmic signals is crucial. We scan for these features in galaxy
clustering using the template

$\delta\mathcal{P}_{\mathcal{R}}(k)=A_{\log}\sin\left(\omega_{\log}\ln\left[\mathcal{K}\right]+\phi\right)\;,$
(3.4)

where $\mathcal{K}=k/k_{\star}$ with the pivot scale
$k_{*}=0.05\,\mathrm{Mpc}^{-1}$, and $A_{\rm log}$, $\omega_{\rm log}$ are the
feature amplitude and frequency, respectively. In Figure 2, we show how this
primordial signal is present in both the primordial power spectrum, on the
left panel, and in the late-time galaxy power spectrum, on the right panel.
Notice that $\omega_{\rm log}$ has no units, so it can not be associated with
a specific value in configuration space.

Global logarithmic features have been widely tested with the CMB using
different datasets. See [35, 36, 37, 34, 38, 39] for an analysis using WMAP
and [2, 40] using Planck data. They have also been tested with LSS data [80,
60, 81]. In [60] it was shown that LSS can have stronger constraints than
Planck in parts of the parameter space, which motivates an in-depth analysis
of current and forthcoming data. Here we push the parameter space beyond what
has been explored in [60] and to the maximum frequency that is accessible with
current LSS datasets. A summary of the parameter space we scan is presented in
Table 1. When compared to previous LSS analysis [60] we have two differences.
Firstly, as briefly scrutinized at the end of §1, our power-spectrum binning
allows us to scan for features over a larger frequency range. Secondly, we
decided to adopt the phase parametrisation instead of two amplitudes because
it is more natural to impose a non-informative prior (see Appendix A for more
details).

###  3.3 Inflation with a step

Oscillatory features can also be generated when the inflaton potential passes
through a step in the potential [82]. Here we use the results of [8] to search
for these features in large-scale data and refer the reader to this paper for
more details. Hereafter this model will be referred to as step features.

Figure 3: left: The primordial curvature power spectrum in the vanilla model
of inflation (black dashed line) and the predictions for a step in the
potential case (solid purple line). For the step feature, we used
$\mathcal{A}_{s}=0.05$, $\omega_{s}=200\;\mathrm{Mpc}$ and $x_{s}=100$; right:
An illustration of how this feature shows up in the evolved galaxy field.

This feature signal was tested in previous work with CMB [2] data. Although
predictions for this model were presented in [83], it has never previously
been tested with LSS datasets. We perform a full MCMC analysis using the BOSS
and eBOSS datasets and test the pipeline with mock catalogues. Moreover, we
also consider the correct effect of the survey window function. As we will
discuss in Section 4.3, the window function convolution imposes a strong
limitation on our constraints when scanning for highly oscillatory features.
For this model, the impact is even more pronounced since the feature amplitude
also decreases for higher frequencies.

An example modification of the galaxy correlation function by this model and
its associated primordial power spectrum is shown in Figure 3. It is
analytically described by

$\displaystyle\ln\delta\mathcal{P}_{\mathcal{R}}^{\mathrm{s}}(k)=\ln\mathcal{P}_{\mathcal{R}}^{0}(k)+\mathcal{I}_{0}(k)+\ln\left(1+\mathcal{I}_{1}^{2}(k)\right)\;,$
(3.5)

where the kernels are given by:

$\displaystyle\mathcal{I}_{0}$ $\displaystyle=$
$\displaystyle\mathcal{A}_{\mathrm{s}}\mathcal{W}_{0}\left(k\omega_{s}\right)\mathcal{D}\left(\frac{k\omega_{s}}{x_{\mathrm{s}}}\right),$
(3.6) $\displaystyle\mathcal{I}_{1}$ $\displaystyle=$
$\displaystyle\frac{1}{\sqrt{2}}\left[\frac{\pi}{2}\left(1-n_{\mathrm{s}}\right)+\mathcal{A}_{\mathrm{s}}\mathcal{W}_{1}\left(k\omega_{s}\right)\mathcal{D}\left(\frac{k\omega_{s}}{x_{\mathrm{s}}}\right)\right]\;.$
(3.7)

The oscillatory terms $\mathcal{W}_{0}$ and $\mathcal{W}_{1}$ are

$\displaystyle\mathcal{W}_{0}(x)$ $\displaystyle=$
$\displaystyle\frac{1}{2x^{4}}\left[\left(18x-6x^{3}\right)\cos
2x+\left(15x^{2}-9\right)\sin 2x\right]\;;$ (3.8)
$\displaystyle\mathcal{W}_{1}(x)$ $\displaystyle=$
$\displaystyle-\frac{3}{x^{4}}(x\cos x-\sin x)\left[3x\cos
x+\left(2x^{2}-3\right)\sin x\right]\;;$ (3.9)

and $\mathcal{D}(x)$ is a damping term that controls the localization of the
feature:

$\mathcal{D}(x)=\frac{x}{\sinh x}\;.$ (3.10)

For more details about this model we refer the reader to [8]. To put these
corrections in the correct format for our pipeline we will rewrite them as

$\mathcal{P}_{\mathcal{R}}^{\mathrm{s}}(k)=\mathcal{P}_{\mathcal{R}}^{0}(k)\left[1+\mathcal{I}_{1}^{2}(k)\right]e^{\mathcal{I}_{0}(k)}\;.$
(3.11)

Since $\mathcal{I}_{0}$ is an oscillation around zero its exponential will be
an oscillation around unity, so we rewrite the model as

$e^{\mathcal{I}_{0}(k)}\equiv 1+\delta\mathcal{I}_{0}\;,$ (3.12)

so the equation to be tested becomes

$\mathcal{P}_{\mathcal{R}}^{\mathrm{s}}(k)=\mathcal{P}_{\mathcal{R}}^{0}(k)\left[1+\mathcal{I}_{1}^{2}(k)+\delta\mathcal{I}_{0}+\mathcal{I}_{1}^{2}(k)\delta\mathcal{I}_{0}\right]\equiv\mathcal{P}_{\mathcal{R}}^{0}(k)\left[1+\delta\mathcal{P}\right]\;.$
(3.13)

The priors we use for this model are presented in Table 1.

Figure 4: left: The primordial curvature power spectrum in the vanilla model
of inflation (black dashed line) and in the sound feature case (solid pink
line). For this plot we used $\mathcal{B}=-0.5$,
$|\tau_{f}|=175\,\mathrm{Mpc}$ and $\ln\beta=6$; right: An illustration of how
the sound feature shows up in the evolved galaxy field.

###  3.4 Variation in the speed of sound with a Gaussian profile

In an effective field theory framework, the inflaton field is necessarily
coupled to heavier fields. These fields can be integrated out, leading to an
effective single-field model [84, 85, 86]. In this effective field theory,
derivative couplings can lead to a sharp reduction in the speed of sound.
These variations generate features in both the primordial power spectrum and
bispectrum of the adiabatic modes [87, 88, 89]. In [90, 91, 92] the authors
computed the consequences for the primordial power spectrum of a time
variation in the sound speed. For a departure of the speed of sound from unity
they concluded that the primordial power spectrum is

$\begin{gathered}\frac{\Delta\mathcal{P}_{\mathcal{R}}}{\mathcal{P}_{\mathcal{R},0}}=\frac{1}{36}\left[\sin\left(2k\tau_{f}\right)+\frac{1}{k\tau_{f}}\cos\left(2k\tau_{f}\right)\right]\mathcal{D}_{S}-\frac{1}{72}\left[\frac{1}{\left(k\tau_{f}\right)^{2}}\sin\left(2k\tau_{f}\right)\right]k\frac{d}{dk}\mathcal{D}_{S}\end{gathered}\;,$
(3.14)

with the envelope $\mathcal{D}_{S}$ being

$\mathcal{D}_{S}=\frac{4\sqrt{\pi}\mathcal{B}k\tau_{f}}{\sqrt{\beta}}\exp\left(-\frac{k^{2}\tau_{f}^{2}}{\beta}\right)\;.$
(3.15)

The parameters $\mathcal{B}<0$ represents the feature amplitude, $\tau_{f}$ is
the characteristic time of the feature, which sets its frequency and the
parameter $\beta>0$ controls the sharpness of the feature. This PF model has
been tested by [92] and [93, 80] using CMB and LSS data, respectively. In
total, it has three free parameters: $\mathcal{B}$, $\tau_{f}$, and $\beta$,
which makes this model difficult to analyze. For this reason, we follow the
discussion in [91] to motivate a preferred region in parameter space to look
for this feature. In [94] it is argued that for the perturbative treatment of
the effective single-field dynamics to remain valid, the sharpness parameter
should satisfy $\ln\beta<14$ [91]. We then decided to scan over the region
$2.2<\ln\beta<10$. We summarize the priors in Table 1.

##  4 Datasets

In this work, we scan for signals of primordial features in large-scale
structure. Our results are obtained from a joint analysis of both the BOSS
low-z and high-z datasets as well as the eBOSS QSO sample. In this section, we
briefly describe the aspects of each dataset, which are summarized in Table 2.
222All the products used in this work can be downloaded at
https://fbeutler.github.io/hub/hub.html

Figure 5: We show how the survey window function affects primordial oscillatory features in large-scale structure data. We generated fictitious linear and logarithmic oscillations (fixed amplitude but different frequencies) and convolved them with the window functions we use in this work (blue and green lines for eBOSS and BOSS, respectively). What we plot above is the wiggly part of the resulting power spectrum. The red line is the spectrum with no primordial features, so only the usual BAO oscillations are present. Note that overall the window function smearing is more important for oscillations with large frequencies. For linear features, the BOSS window function completely erases the primordial feature signal for $\omega_{\rm lin}=2500\;\mathrm{Mpc}$. For logarithmic features, the window function effect is not straightforward, since the frequency is a function of $k$. Nevertheless, the smearing always occurs from small values of $k$. Regardless of the models considered in this work, all PF signals are erased in the limit $\omega\rightarrow\infty$. | Lin. feature | Log. feature | Step feature | Sound feature
---|---|---|---|---
$A_{\rm X}$ | [$-1$, 1] | [$-1$, 1] | … | …
$\omega_{\rm X}$ | [100, 4000] | [10, 360] | … | …
$\phi_{\mathrm{X}}/\pi$ | [$-0.5$, $0.5$] | [$-0.5$, $0.5$] | … | …
$A_{\rm step}$ | … | … | [$-1$, 1] | …
$\omega_{\rm step}$ | … | … | [100, 1000] | …
$x_{s}$ | … | … | [10, 650] | …
$\mathcal{B}$ | … | … | … | [$-1$, 0]
$\ln\beta$ | … | … | … | [2.2, 10]
$|\tau_{f}|$ | … | … | … | [100, 900]
Table 1: Ranges for the uniform priors we use for the four primordial feature
models we investigate in this paper. The parameters $\omega_{\rm lin}$,
$\omega_{\rm step}$, $\tau_{f}$ and $x_{s}$ are in units of Mpc.

###  4.1 BOSS data

The Baryon Oscillation Spectroscopic Survey (BOSS) was part of the Sloan
Digital Sky Survey III (SDSS-III) and, at the present moment, is the dataset
with the largest effective volume [95, 96, 97, 98]. Here, we use data release
12 (DR12) [99], which consists of 1 198 006 galaxies divided into two patches
in the sky, the North Galactic Cap (NGC) and the South Galactic Cap (SGC). We
follow the standard approach and split the data in two (independent) redshift
bins: the low-$z$ sample, with $0.2<z<0.5$, and the high-$z$ sample, with
$0.5<z<0.75$. As we will explain in §4.3, the $k$-binning choice is relevant
for the analysis of the primordial features, since smaller bins usually allow
us to probe PF with higher frequencies. We can not decrease the bin size
indefinitely, since it is limited by convolution with the window function:
this has a width of order the survey fundamental mode, $k_{f}$, and so modes
separated by less than this amount will be strongly correlated. Moreover, the
estimation of the covariance matrix will require more mock catalogues in order
to be reliable. With that in mind, we decided to use a binning that is close
to the survey fundamental mode defined in [100], which is $\Delta
k=0.002h\,\mathrm{Mpc}^{-1}$ for both NGC and SGC 333In principle one could
use a different binning for NGC and SGC since they have different volumes,
which allows us to probe primordial features with higher frequencies. But the
signal at such high frequencies is already dominated by the eBOSS sample so
this would not greatly affect the final result.. Since the inflationary signal
is isotropic, we only consider the power spectrum monopole in our analysis. In
recent work [81], the authors included the quadrupole in the analysis.
Although the feature signal is isotropic, the window function can leak some
information from the monopole into the higher-order multipoles. However, this
leakage is small and the quadrupole naturally has a smaller signal to noise
than the monopole. Moreover, the feature amplitude is expected to be small, so
the inclusion of the quadrupole will not bring significant additional
information. We also emphasize that an anisotropic analysis requires only a
trivial extension of the pipeline, which can easily be implemented in the
future.

Another product we will use is the MultiDark-Patchy mock catalogues that
reproduce the BOSS DR12 dataset [101]. These catalogues are obtained using
approximate solvers for the complicated gravitational collapse models,
reproducing a similar galaxy number density, survey geometry, and selection
function. We refer the reader to the reference above for more details. Their
use here is twofold. Firstly, we use the $997$ NGC and $1000$ SGC mock
catalogues to derive the covariance matrix required for the MCMC analysis.
Secondly, since these mocks were generated without primordial features, they
can be used as a validation test of our analysis pipeline for the case without
features.

To further improve the sensitivity to oscillatory signals in the BOSS data, we
perform a field reconstruction as described in [102]. This helps reduce the
effects of non-linear displacements, which as described in §2.3 can smear out
oscillations present in the power spectrum. In addition, the displacements can
also cause a change in the feature scale, shifting the oscillation frequency.
In [103], the authors tested the impact of a similar field reconstruction
scheme and concluded that it alleviates these two problems.

Figure 6: The oscillation half period for logarithmic features. Notice that
the frequency of the oscillation varies with $k$. The solid and dashed grey
lines represent the half period of linear oscillations with frequencies equal
to the Nyquist frequency for $\Delta k=0.002h\,\mathrm{Mpc}^{-1}$ and $\Delta
k=0.001h\,\mathrm{Mpc}^{-1}$, respectively.

###  4.2 eBOSS data

Sample | $[z_{\rm min},z_{\rm max}]$ | $z_{\rm eff}$ | $\Delta k\;[h\mathrm{Mpc}^{-1}]$ |  | $[k_{\rm min},k_{\rm max}]$
---
$\;[h\mathrm{Mpc}^{-1}]$
| Nyquist Freq.
---
$\;[\mathrm{Mpc}]$
BOSS low-$z$ | [0.2, 0.5] | 0.38 | 0.002 | [0.01, 0.3] |  2323
BOSS high-$z$ | [0.5, 0.75] | 0.61 | 0.002 | [0.01, 0.3] |  2323
eBOSS Quasars | [0.8, 2.2] | 1.52 | 0.001 | [0.02, 0.23] |  4647
Table 2: Summary of the data used in this work, presented in Section 4. For a
binned data vector there is a maximum frequency one can probe. This
fundamental limit due to the discrete nature of the data is called the Nyquist
frequency, defined as $\omega_{\rm Ny}=\pi/\Delta k$. We adopt the same
fiducial cosmology used in the BOSS and eBOSS official analysis. To convert
between units with and without the Hubble factor, we use $h_{\rm fid}=0.676$,
which is common to both fiducial cosmologies. In Figure 8 we also show BOSS
constraints beyond the Nyquist frequency defined above. In order to do this,
we switch the BOSS binning to be the same as eBOSS for frequencies above
$\omega_{\rm lin}\approx 2323\,\mathrm{Mpc}$.

The extended Baryon Oscillation Spectroscopic Survey (eBOSS) [104, 105] is
part of SDSS-IV [106]. Here we use the QSO sample, which is composed of $\sim
350\,000$ quasars. Although the number density is considerably lower when
compared to the BOSS dataset described previously, it covers a larger cosmic
volume. The survey volume is crucial when looking for primordial features, and
affects the analysis through the window function. A larger cosmic volume
allows us to probe features with higher frequencies/smaller wavelength. The
eBOSS QSO dataset is also divided between two galactic caps (NGC and SGC) and
has no redshift overlap with the BOSS galaxy sample. We also use QSO mock
catalogues to validate our pipeline for that dataset [107]. These mocks were
obtained by implementing the Extended Zel’dovich (EZ) approximate N-body
simulation technique [108]. The eBOSS QSO sample has a cosmic volume that is
considerably larger than the BOSS sample, which allows a smaller $k$-binning:
$\Delta k=0.001h\,\mathrm{Mpc}^{-1}$. We use 1000 mock catalogues to compute
the eBOSS QSO covariance matrix.

###  4.3 Window function

All data collected from surveys are formed as a combination of the true
density fluctuations multiplied by weights that account for the survey
geometry and systematic effects. In Fourier space this weight appears as a
convolution, modifying the shape of the observed power spectrum. We estimate
the window function from the random catalogue and convolve the theory with it,
following the steps described in [109]. In summary, the steps are:

  * •

Compute the theory in Fourier space and Hankel transform it to obtain the
theoretical correlation function.

  * •

Calculate the convolved theoretical correlation function by multiplying it
with the window function. For the monopole, this multiplication is trivial
because we are ignoring the quadrupole.

  * •

Perform an inverse Hankel transform on the outcome of the last step to obtain
the convolved power spectrum. This is the quantity we use to compare with the
data.

In Figure 5 we show a series of examples of the influence of the window
function for both linear (first four panels on the top) and logarithmic
features (last four panels). The greyish blue line is the usual BAO wiggle,
and the other lines represent the results with additional primordial features
with fixed amplitude $A_{\rm lin,log}=0.075$. In general, the window function
smears out any oscillatory signal, and this smearing is stronger for higher-
frequency oscillations. This trend can be seen for both linear and logarithmic
oscillations. For instance, for linear frequencies $\omega_{\rm lin}\approx
800\,\mathrm{Mpc}$ the resulting oscillatory pattern is almost the same for
all survey geometries, with the damping being stronger for BOSS than eBOSS. As
the frequency is increased, the BOSS window function tends to degrade the
oscillations with respect to eBOSS even further, and the signal is completely
removed for $\omega_{\rm lin}\approx 2500\,\mathrm{Mpc}$. This limiting
frequency is higher for the eBOSS QSO dataset because it has a greater cosmic
volume. For small frequency values, the BOSS constraints will overcome the
eBOSS ones, simply because BOSS has a greater signal-to-noise; but because
BOSS loses sensitivity to high frequencies first, we expect eBOSS constraints
to dominate at high frequencies. For the logarithmic features, a similar trend
is present. In Figure 6 we plot the separation between zeroes for a
logarithmic oscillation with different frequencies. The horizontal lines
represent the Nyquist frequencies, $\omega_{\rm Ny}=\pi/\Delta k$, for $\Delta
k=0.002h\,\mathrm{Mpc}^{-1}$ and $0.001h\,\mathrm{Mpc}^{-1}$. A survey window
function with a fundamental mode close to these values will erase all signal
whose separation of zeros is below $\omega_{\rm Ny}$. Therefore, according to
Figure 6, the smearing due to the window function for logarithmic oscillations
will occur from lower values of $k$. This is confirmed in Figure 5, and as the
frequency increases the signal starts to disappear first for small values of
$k$.

##  5 Statistical methods and analysis pipeline

In this section, we explain how we test the primordial feature models and
compare them with the reference single-field slow-roll model.

| BOSS high-$z$ | BOSS low-$z$ | eBOSS QSO
---|---|---|---
$\alpha_{\rm iso}$ | 0.9887 $\pm$ 0.0089 | 0.997 $\pm$ 0.0105 | 1.023 $\pm$ 0.020
Official results | 0.9887 $\pm$ 0.0087 | 1.000 $\pm$ 0.010 | 1.025 $\pm$ 0.020
Table 3: Results for a preliminary BAO isotropic analysis using our pipeline
with no primordial features. The measurements of $\alpha_{\rm iso}$ are in
very good agreement with the ones quoted by previous works [68, 66]

###  5.1 Performance metrics

If the hypothesis of primordial features is correct, then there will be two
kinds of oscillations in the power spectrum: the first one is the usual baryon
acoustic oscillation (BAO) feature, whose physics is well-known and which is
one of the main science drivers for modern galaxy redshift survey experiments.
The second one is the primordial feature, which couples to the BAO and can be
written as in Eq. (2.12). We take as reference model the simple situation in
which only the BAO oscillation exists; we call this model $\mathcal{M}_{0}$.
The main goal is then to answer the following question: Given some primordial
feature model, $\mathcal{M}_{1}$, from the options set out in Section 3 and
the best existing clustering data, how likely is $\mathcal{M}_{1}$ to be true
relative to $\mathcal{M}_{0}$?. This is a model comparison problem, which can
be solved using standard Bayesian methods [110, 111, 112]. If
$\mathcal{M}_{0}$ and $\mathcal{M}_{1}$ are taken as equally likely a priori,
then the relative probability (odds ratio) between the two models is given by
the Bayes factor, $B_{01}$:

$B_{01}=\frac{p(D|\mathcal{M}_{0})}{p(D|\mathcal{M}_{1})}\;,$ (5.1)

where $p(D|\mathcal{M}_{i})$ is the evidence of model $i$, given some data
$D$. Usually, the computation of the evidence is numerically costly, since it
involves the evaluation of a multi-dimensional integral. Nevertheless, note
that the primordial feature models are extensions of $\mathcal{M}_{0}$, and
they coincide when the feature amplitude, $A$, vanishes. In other words
$\mathcal{M}_{0}$ is nested inside $\mathcal{M}_{1}$, and it has been shown
that Eq. (5.1) can then be evaluated using the Savage-Dickey ratio [113, 111]:

$B_{01}=\frac{p(A|D,\mathcal{M}_{1})}{p(A|\mathcal{M}_{1})}\Big{|}_{A=0}\;,$
(5.2)

where $p(A|D,\mathcal{M}_{1})$ is the marginal posterior for $A$ under model
$\mathcal{M}_{1}$ given the data and $p(A|\mathcal{M}_{1})$ is the prior for
$A$ in that model. Hence the only ingredient we need in order to compute the
odds ratio is the marginal posterior of the feature amplitude for each model.
We estimate the posterior using samples obtained via the emcee package [114].
For each analysis, we run four parallel chains and use the Gelman-Rubin
criterion for convergence. After the chains converge, we apply a 40% burn-in
and combine them into one final chain. Further analysis is performed using the
package getdist [115]. For step and sound features, the amplitude is positive-
definite, so the samples we obtain near $A=0$ might not be representative. For
these models, we therefore computed the Bayes factor using the pocoMC software
[116, 117]. When necessary, we will also quote the minimum $\chi^{2}$ for
different models. This value is found using the iminuit software [118], which
is a Python code based on the C$++$ library Minuit [119]. Our main results are
obtained by the joint constraints from all three data products we use: BOSS
low-$z$, high-$z$ and the eBOSS QSO sample, which we assume to be
uncorrelated, so the final combined posteriors are obtained by multiplying
each individual posterior.

In addition to the Bayes factor, we also compute the 95% credible intervals
for the primordial feature amplitude as a function of its frequency. To do
this, we first split the MCMC samples into frequency bins to obtain the binned
posterior. Hence, we analyze the marginal posterior on the feature amplitude.
We start integrating it around $|A|=0$ and stop when the cumulative
probability reaches $95\%$. Thus, for a given frequency bin $i$, we will have
a value $A_{i}$ such that $p(|A|\leq A_{i}|D)=95\%$. High values of $A_{i}$
indicate the marginal posterior on $A$ for the frequency bin $i$ is distant
from the origin, which can be associated with a possible detection.

Figure 7: The combined result of the analysis when looking for linear,
logarithmic, step and sound features using mock catalogues of BOSS high-$z$,
low-$z$, and eBOSS QSO sample. In the left panels, we show the median of the
posterior for each mock individually (light colours), as well as their sample
mean (dark colour). In the right panels, we show the 95% credible interval for
the PF amplitude, $|A|$, obtained by integrating the posterior.

###  5.2 Validating the pipeline

The mocks are constructed so that the only feature on top of the smooth
broadband power spectrum is the BAO signal. However, apparent primordial
features can be generated by noise, and it is convenient to use an ensemble of
mocks to estimate the probability of observing a given signal just by chance,
particularly given the need to allow for the look-elsewhere effect as
discussed below in §6.4. We perform such a test for all models that we
introduced in Section 3 and discuss the results here. The comparison of the
mock analysis with data is performed in Section 6.

Starting with linear features, we scan the range $\omega_{\rm
lin}\in[100,4000]\;\mathrm{Mpc}$ and estimate a posterior on $A_{\rm lin}$
using the MCMC samples with a bin size of $\Delta\omega_{\rm
lin}=10\;\mathrm{Mpc}$. Finally, we compute the statistical metrics we
introduced in 5.1. The result of this procedure applied to mocks is shown in
the upper panel of Figure 7. On the left, we show the median of the
marginalised posterior on $A_{\rm lin}$, whereas, on the right, we show its
95% credible interval. The light blue lines are the result for each mock, and
the dark blue is the mean of 100 mocks. Notice that the credible interval on
$A_{\rm lin}$ is approximately constant for $\omega_{\rm
lin}\approx[100,500]\;\mathrm{Mpc}$, but tends to increase for higher
frequencies. What is causing this degradation is the window function
convolution, as discussed in Section 4.3. Also, note that in both left and
right figures, there is an increase in noise at $\omega_{\rm
lin}=150\;\mathrm{Mpc}$, which arises when the primordial feature becomes
approximately degenerate with the BAO oscillation. Another point worth
mentioning is that the 95% credible interval seems to saturate for
$\omega_{\rm lin}\approx 4000\;\mathrm{Mpc}$. It happens because the window
function degradation cannot occur perpetually, since we are also imposing a
uniform prior on the amplitude (see Table 1). At some point, this degradation
saturates the prior, and, consequently, the scatter of $A_{\rm lin}$ around
the origin. We can note that the scatter in median $A$ values rises with
frequency in a similar way to the increase in 95% confidence limits, but that
the occasional mode has a preferred value above that limit: this is exactly
the look-elsewhere effect in operation, although the effect seems to be modest
(see §6.4 for a detailed discussion).

We now move to the results for logarithmic features. In this case, we scan for
features with $\omega_{\rm log}\in[10,360]$ and bin the resulting posterior
using a bin size of $\Delta\omega_{\rm log}=1$. The results are presented in
the middle panel of Figure 7. Each light orange line represents the outcome of
a single mock, whereas the dark orange line represents the mean of 25 mocks.
The window function degradation is present once again, but the way it acts on
the amplitude differs from the case of linear oscillations (see the discussion
at the end of Section 4.3). For the frequency range considered here, the error
is not expected to saturate the prior. The degeneracy between logarithmic
features and the BAO is not as prominent as in the case of the linear feature.

Next, we consider step features. The results for this model are in the second-
last panel of Figure 7. We scan the range $\omega_{\rm
step}\in[100,800]\;\mathrm{Mpc}$, and bin the posterior with a bin size
$\Delta\omega_{\rm step}=10\;\mathrm{Mpc}$. Although this model is made of
linear oscillations, its error increases more rapidly than for the linear
features. The reason for this is that according to Eqs. 3.8 and 3.9, the
overall oscillation amplitude is decreased as we move to higher values of
$\omega_{\rm step}$. This, when added to the natural window function effect,
limits the constraining power of this model. The Planck 2018 inflation paper
quotes a best-fit step model with a frequency $\omega_{\rm step}\sim 1000$
[120]. As we will see in the next section, the current LSS constraining power
in this range of parameter space is quite poor.

Finally, for sound features, we scanned the range $|\tau_{f}|\in[75,900]$ and
binned the posterior with a bin size of $\Delta\tau_{f}=10\;\mathrm{Mpc}$. The
prior on the other parameters are in Table 1. The results are presented in the
bottom panel of Figure 7. As for the step feature, discussed in the previous
paragraph, the degradation of this PF amplitude due to non-linear damping,
window-function convolution, and envelope function poses a strong obstacle
when searching for this feature in LSS datasets. Our forecast with the mocks
indicates that, even at small frequencies, the $95\%$ credible interval on
$\mathcal{B}$ almost recovers the prior for all values of $\tau_{f}$. We thus
already expect the data to be non-informative, and the model comparison should
be inconclusive, with $\ln B_{01}\approx 0$. This will be confirmed later on
with BOSS and eBOSS data in Section 6.

Before moving to constraints from the real data, we decided to carry on a last
sanity check. Since the primordial feature analysis is an extension of an
isotropic BAO pipeline, we studied whether we could recover BAO constraints
present in the literature. Our measurements are presented in Table 3. Starting
with the BOSS high-$z$, we find the same median for $\alpha_{\rm iso}$ as
quoted in the official analysis [68], but with an error approximately 2%
larger. We also found a good agreement for BOSS low-$z$, with our error bar
approximately $5\%$ larger. For eBOSS, we recover the same error as in the
official analysis [66] and the median of our result is 0.1$\sigma$ away from
the official result. We, therefore, conclude our BAO pipeline is consistent
with previous results.

##  6 Large-Scale Structure constraints

In this section, we present the large-scale structure constraints on
primordial features obtained by combining both BOSS and eBOSS data. For each
model we compute the general metrics we introduced in Section 5. We scrutinize
our results in Section 7.

Figure 8: top-left: Median of the marginalised posterior on $A_{\rm lin}$ for
mocks and data. We also show the 1 and 2$\sigma$ regions from the scatter of
the 100 mocks (light blue and grey backgrounds). top-right: The 95% credible
interval on $A_{\rm lin}$ for both data and mocks. The background represent
the same as before. bottom: How each data product contributes to the combined
result for the credible interval. The current LSS data discard the blue
coloured region for linear features at 95% confidence level.

###  6.1 Linear features

In Figure 8 we show the results for linear features. The upper-left panel
shows the median of the marginalised posterior for $A_{\rm lin}$ as a function
of the feature oscillation. Notice that we included BOSS constraints above the
Nyquist frequency given in Table 2, which is $\omega_{\rm
lin}=2323\,\mathrm{Mpc}$. These were obtained by switching the binning of the
BOSS power spectrum to be the same as for eBOSS. These results extend the
range of the combined constraint, and they emphasize the different scaling of
the error with the frequency. We use a bin size of $\Delta\omega_{\rm
lin}=10\,\mathrm{Mpc}$ starting with the upper left plot, where we see that
although the observed medians for $A_{\rm lin}$ lie inside the scatter from
the mocks, they are not as noisy, especially as we move to higher frequencies.
A consequence is that the credible interval for the data should be tighter
than almost all mocks, which is corroborated in the upper-right plot, where we
compare the 95% credible intervals. We conclude that the observed dataset is
an outlier of the mocks, and such a discrepancy was reported previously in the
literature [60]. We discuss the consequences of this discrepancy between data
and mocks in the next section. In the lower panel, we show how each dataset
contributes to the total constraint on $A_{\rm lin}$. At small frequencies,
where the window function convolution is not predominant, the BOSS dataset
dominates the constraint. At $\omega_{\rm lin}\approx 2500\;\mathrm{Mpc}$
there is a cross-over between BOSS and eBOSS, and thereafter eBOSS dominates
over BOSS. Finally, in the top panel of Figure 13 we show the Bayes factor for
linear PF as a function of their frequency.

Figure 9: top-left: Mean of the marginalised posterior on $A_{\rm log}$ for
mocks and data. We also show the mean result of the mocks. We also show the 1
and 2$\sigma$ regions from the scatter of the 100 mocks (light blue and grey
backgrounds) top-right: The 95% credible interval on $A_{\rm log}$ for both
data and mocks. The background represent the same as before. bottom: How each
data product contributes to the combined result for the credible interval. The
current LSS data discard the orange coloured region for logarithmic features
at 95% confidence level.

Figure 10: top-left: Mean of the marginalised posterior on $A_{\rm step}$ for
mocks and data. We also show the mean of the mocks as well as the 1 and
2$\sigma$ regions from the scatter of the 100 mocks (light blue and grey
backgrounds) top-right: The 95% credible interval on $A_{\rm step}$ for both
data and mocks. The background represents the same as before. bottom: How each
data product contributes to the combined result for the credible interval. The
current LSS data discard the lavender coloured region for step features at 95%
confidence level.

###  6.2 Logarithmic features

In Figure 9, we show the results for logarithmic features. In the upper-left
and upper-right panels, we show the median and the 95% credible interval for
$A_{\rm log}$, respectively, as a function of frequency. Here we used
$\Delta\omega_{\rm log}=1$. In the bottom panel, we show how each dataset
contributes to the combined credible interval. There is no clear cross-over
between BOSS and eBOSS, and both datasets are equally informative for
$\omega_{\rm log}\sim 300$. There is a bump in the credible interval for
$\omega_{\rm log}\sim 280$, which indicates that around this frequency the
null hypothesis may be excluded with high significance. In Figure 14, we
provide additional information about this peak. On the left panel, we show the
$\Delta\chi^{2}$ relative to the vanilla inflation model for all datasets, and
in the right panel, we show the marginalized posterior for that frequency bin.
The Bayes factor for logarithmic features is shown in the middle panel of
Figure 13.

###  6.3 Step and sound features

We show the large-scale structure constraints on step features in Figure 10,
following the same structure as linear and logarithmic features. For this
result, we used a bin size of $\Delta\omega_{\rm step}=10\,\mathrm{Mpc}$. For
this PF model, data and mocks seem to agree well. All datasets have similar
constraining power for this model, saturating the prior at the same frequency
$\omega_{*}\approx 1230\,\mathrm{Mpc}$. The Bayes factor for step features as
a function of the frequency is shown in the bottom panel of Figure 13.

In Figure 11, we display the two-dimensional marginalized posterior for
$\tau_{f}$ and $\ln\beta$ for each dataset and their combined result. Notice
that the current LSS datasets have favoured regions in parameter space around
some values of $\tau_{f}$, which are the darker horizontal contours in this
figure.

Figure 11: The $1\sigma$ and $2\sigma$ contours of the marginalized posterior
on $\tau_{f}$ and $\ln(\beta)$ for sound features using the most recent LSS
datasets. Although we have a weak constraint on the PF amplitude,
$\mathcal{B}$, this plot shows preferred oscillatory modes (or preferred PF
scale, $\tau_{f}$), which are the darker horizontal contours. The Bayes factor
for this model obtained from BOSS low-$z$, high-$z$ and eBOSS Quasars are $\ln
B_{01}=0.620,0.093,0.045$, respectively. It is therefore neither favoured nor
discarded by current LSS data.

###  6.4 The look-elsewhere effect

Figure 12: left: The CDF of the maximum values of $|\bar{A}|$ for individual
frequency channels for linear (upper) and logarithmic (bottom) feature. The
black-dashed line is the CDF measured from the data, whereas the coloured
lines are the best-fits of the template
$\tilde{P}(<x)=1-\exp\left[{-(x/b)^{a}}\right]$. The values of $(a,\;b)$ for
linear and logarithmic features are $(1.61,\;1.18)$ and $(1.44,\;1.06)$,
respectively; right: How statistically relevant the events in the data are, at
individual frequency bins, in light of the look-elsewhere effect. The red
point indicate the $3.99\sigma$ event we see in the data for logarithmic
features at $\omega_{\rm log}\approx 286$. When the look-elsewhere effect is
taken into account, the significance of this event is reduced to to
$1.22\sigma$. Similarly, the blue point indicate a $3.31\sigma$ event for
linear features at $\omega_{\rm lin}\approx 2605\,\mathrm{Mpc}$. The look-
elsewhere effect reduces its significance to $0.39\sigma$.

In all of this, it is important to note that considering a number of different
frequency bins makes our analysis susceptible to the ‘look-elsewhere’ effect
[121, 122]. In other words, an analysis of real data might indicate a
preferred value of $|A|$ for some frequency that is well above the 95% limit
derived from analyzing mocks – but this need not constitute a true detection
of primordial features if there was no reason to focus on that frequency in
advance. The more frequency channels we inspect, the more likely it is to have
one with an outlying value of $|A|$. In light of this effect, the statistical
significance of events at a specific frequency is reduced. The impact of the
look-elsewhere effect can be understood by analyzing the distribution of the
maximum amplitudes that can be measured under the null hypothesis,

$|\bar{A}|^{\max}=\max_{\phi,\omega}|\bar{A}|\left(\omega,\phi\right)\;,$
(6.1)

where $\bar{X}\equiv\hat{X}/\rm{Var}(X)^{\frac{1}{2}}$. We are thus interested
in computing

$P_{N_{\rm eff}}\left(|\bar{A}|^{\max}\geq x\right)=1-P_{N_{\rm
eff}}\left(|\bar{A}|^{\max}<x\right)=1-\left[P_{1}\left(|\bar{A}|^{\rm
max}<x\right)\right]^{N_{\rm eff}}\;,$ (6.2)

where $N_{\rm eff}$ is the effective number of independent frequency channels,
$P_{N}\left(|\bar{A}|^{\max}<x\right)$ is the cumulative distribution function
(CDF) of the distribution of maximum $|\bar{A}|$ over a total of $N$ frequency
channels, and $P_{1}\left(|\bar{A}|^{\max}<x\right)$ is the same but for a
single frequency bin. The safest way to evaluate Eq. (6.2) is to have an
ensemble of reliable mocks with no PF, from which we can construct the
probability distribution function (PDF) of $|\bar{A}|^{\rm max}$, and
consequently compute its CDF, $P(|\bar{A}|^{\rm max}<x)$, to be used in Eq.
(6.2). But unless the mocks are completely realistic, this approach may not be
reliable; we therefore make a more indirect argument, focusing on the
effective number of independent frequency bins, $N_{\rm eff}$.

To estimate $N_{\rm eff}$, we measure both $P_{N}$ and $P_{1}$ from an
ensemble of mocks, plug them into the last equality of Eq. (6.2) and solve for
$N_{\rm eff}$. For linear and logarithmic features, we found that $N_{\rm
eff}\approx 220$ and $215$, respectively. The associated effectively
independent frequency bins are $\Delta\tilde{\omega}_{\rm lin}\approx
17.70\;\mathrm{Mpc}$ and $\Delta\tilde{\omega}_{\rm log}\approx 1.72$, both a
little larger than the adopted values of $\Delta\omega_{\rm
lin}=10\;\mathrm{Mpc}$ ($N=390)$ and $\Delta\omega_{\rm log}=1$ ($N=370$),
respectively.

Figure 13: The combined Bayes factor (from BOSS and eBOSS datasets) for
linear, logarithmic, and step features. Peaks above the black dashed line
favour the PF model over the vanilla inflationary model in that frequency
band. The opposite holds for points below the line. According to Jeffreys’
scale [123], moderate and strong evidence is associated with $|\ln
B_{01}|\approx 2.5$ and $5.0$, respectively. Apart from a single frequency bin
for logarithmic features, the current LSS data strongly favour the vanilla
inflationary scenario over the PF signals we studied in this work. We used the
following bins for the frequencies: $\Delta\omega_{\rm lin}=10\,\mathrm{Mpc}$,
$\Delta\omega_{\rm step}=10\,\mathrm{Mpc}$ and $\Delta\omega_{\rm log}=1$. For
sound features the posterior is dominated by the prior,
$\ln\left(B_{01}\right)\approx 0$ for all $\tau_{f}$, hence we decided to omit
this model from this study of the Bayes factor.

Now, as shown in Appendix B, the eBOSS and BOSS mocks are not a perfect
representation of the data, and they have a noisier constraint on the
amplitudes of PF and thus a different single-channel amplitude distribution,
$P_{1}(|A|)$. We did however note that $P_{1}(|\bar{A}|<x)$ for both mocks and
data is well-fitted by the form
$\tilde{P}(<x)=1-\exp\left[{-(x/b)^{a}}\right]$. Our approach is therefore to
use the $P_{1}$ function measured from data but to allow for the look-
elsewhere effect by using the $N_{\rm eff}$ values quoted above, obtained from
the mocks. We believe this to be a defensible procedure, since the degree of
correlation between frequency channels is set by the survey window, and should
not depend on the detailed nature of the mock fluctuations. For convenience,
we convert the $p$-value computed using Eq. (6.2) into the corresponding
$S-\sigma$ point of a Gaussian (via $S=2^{1/2}\operatorname{Erf}^{-1}[1-p]$,
where $\operatorname{Erf}^{-1}$ is the inverse error function).

In Figure 12, in the left panel, we show how well $P_{1}(|A|<x)$ measured
using the data is fitted by $\tilde{P}$. We use this to map the real
significance of events in the data at individual frequency bins, and we
present the result in the right panel of Figure 12. The $x$-axis in that
figure represents the values of $|A|$ normalized to its standard deviation,
and in the $y$-axis we show its real significance when the look-elsewhere
effect is taken into account. We also included two points representing the
strongest event we detected in the data for linear and logarithmic features.
Note how the look-elsewhere effect reduces the significance of single-channel
events.

|  | $\Sigma_{\rm nl}\,[h^{-1}\,\mathrm{Mpc}]$ |
---|---|---|---
| BOSS low-$z$ | BOSS high-$z$ | eBOSS Quasars
Mean of mocks | $6.00\pm 0.26$ | $5.66\pm 0.27$ | $13.57\pm 1.46$
Data | $4.1^{+1.6}_{-1.3}$ | $2.6^{+1.3}_{-1.8}$ | $2.44^{+0.94}_{-2.2}$
Table 4: The value of $\Sigma_{\rm nl}$ for each dataset compared to their
respective mocks. For each dataset, we performed a total of 100 fits from the
mocks and a single one from the data. The figure we quote above for the mocks
is the average best fit and its standard deviation, and for the data, we quote
the $1\sigma$ range obtained from the posterior. For all datasets the value
fitted from the mocks is larger than the data. This inconsistency between the
mocks and the data poses problems with the interpretation of our results (see
§7 for more details).

##  7 Discussion

In this section, we comment on the results presented in Section 6. First, we
scrutinize the comparisons between data and mocks and explain why they do not
have similar constraints for global oscillatory features. Finally, we comment
on the credible intervals and the Bayes Factors we obtained for each model,
which constitutes our main results.

The results for linear and logarithmic features in Figures 8 and 9,
respectively, show that the observed dataset is an outlier of the mocks. In
[60], the authors found the same trend, and the reason is the same as in the
present work: mock catalogues prefer higher values for the non-linear damping
parameter, $\Sigma_{\rm nl}$. This leads to a higher scatter of the mean
amplitude and a more substantial error. In Table 4, we compare the best-fit
values of $\Sigma_{\rm nl}$ between data and mocks. The error we quote for
mocks is the sample standard deviation of the best fit of the 100 mocks we
analyzed. We see that the observed result is $\sim 5\sigma$ away from the
mocks for BOSS, and we find similar results for eBOSS. It reveals that the
non-linear damping in the mocks is not consistent with the data. These mocks
were generated using approximate methods, so we speculate that it may have
induced non-physical smearing of features, probably due to some resolution
issue. For BOSS Patchy mocks, the problem is less severe (see also the
discussion in §8.3 of [68]). Therefore, it explains why the disagreement we
see in Figure 8 is not that extreme for low PF frequencies, which is the
region in parameter space dominated by BOSS. In Appendix B, we show that if
$\Sigma_{\rm nl}$ is fixed to the best-fit extracted from the data or N-body
simulations, the PF amplitudes and errors fitted from the mocks are consistent
with the data for all frequencies.

Another trend in our results is the increase of the error on the PF amplitude
with its frequency. As we discussed in §4.3, this is caused by the window
function convolution, which will smear out any oscillatory signal present in
the data. The smearing will be stronger for higher frequencies. Moreover,
since the cosmic volume of the three datasets we used are different, the
window function impact will be different. This effect is perceptible at the
bottom panel of Figures 8, 9 and 10, where we show how each dataset
contributes to the combined constraint. Notice that this effect completely
deteriorates the BOSS results for linear features at higher frequencies,
causing the eBOSS QSO sample to dominate over BOSS.

In Figures 8, 9 and 10, we also show the median and credible interval on the
PF amplitude for the linear, logarithmic, and step models. These are the
combined results of BOSS and eBOSS datasets. The coloured regions in the
credible interval plots indicate the volume in parameter space excluded by
current data at 95% significance. In Figure 13, we show the Bayes factor,
$-\ln B_{01}$, for these three models as a function of the frequency, which is
the main result of this work. The Bayes factor is interpreted as the odds
ratio between two models 444For instance, $B_{01}=2$ means that the model
$\mathcal{M}_{0}$ is twice as likely to be true relative to $\mathcal{M}_{1}$,
given the data.. Peaks above the dashed line in Figure 13 indicate that
current data prefer the PF over the vanilla model, and the opposite holds for
points below that line. The significance of these points can be interpreted
according to the Jeffreys’ scale [123]: moderate and strong evidence is
associated with $|\ln B_{01}|\approx 2.5$ and $5.0$, respectively. Note that
$B_{01}$ asymptotically approaches unity for high frequencies, a trend that
will be present for all PF models. It happens because, for all oscillatory PF
models, there is an associated limiting frequency that can be probed with a
given cosmic volume, $\omega_{*}$. Therefore, the error on the PF amplitude,
$\sigma\left(|A|\right)$, will increase indefinitely as
$\omega\rightarrow\omega_{*}$ . In this regime, the data become less
informative and consequently the posterior starts to be dominated by the
prior, making the model comparison inconclusive [111].

Starting with linear features, in the upper panel of Figure 13, we see that
apart from some specific bins of $\omega_{\rm lin}$ the featureless
inflationary model is favoured over the linear feature model. According to
Jeffreys’ scale, all the significant events support the simple inflationary
model. The data have four 3$\sigma$ events and no 4$\sigma$ events. None of
the 3$\sigma$ events are in the range $\omega_{\rm lin}\leq
900\;\mathrm{Mpc}$, in agreement with [60]. Regarding the events at larger
frequencies, we find that approximately 20% and 7% of the mocks have at least
one 3$\sigma$ and 4$\sigma$ event, respectively. It indicates that the peaks
we see are not statistically significant, since we see similar events in the
mocks with modest frequency.

Figure 14: left panel: The $\chi^{2}$ profile for logarithmic features in the
range $\omega_{\rm log}\in[10,380]$. The black solid line is the value for the
featureless model. The more striking downward peaks at $\omega_{\rm
log}\approx 286$ for the BOSS data sets are associated with the combined (BOSS
+ eBOSS) $\ln B_{01}\approx-4$ feature in Figure 9. right panel: The
marginalised posterior for $A_{\rm log}$ in the frequency bin $286<\omega_{\rm
log}<287$, where the logarithmic feature is favoured the most over the
featureless model. The pairs of vertical lines represent the $1\sigma$,
$2\sigma$, $3\sigma$, and $4\sigma$ limits, respectively.

The Bayes factor for logarithmic features is shown in the middle panel of
Figure 13. For $\omega_{\rm log}\sim 286$, we find $\ln B_{01}\approx-4$,
indicating that for this frequency band, there is strong evidence of
primordial logarithmic features. To investigate this event more deeply, we
show in the left panel of Figure 14 the $\Delta\chi^{2}$ profile between the
featureless inflationary model and the logarithmic features models. The grey-
shaded area represents the frequency range where the event occurs. The best-
fit for $\omega_{\rm log}\in[286,287]$ produces
$\Delta\chi^{2}=-8.32,-10.01,-4.20$ for BOSS low-$z$, high-$z$ and eBOSS,
respectively. In the right panel of the same figure, we show the marginalized
posterior for $A_{\rm log}$ with $\omega_{\rm log}\in[286,287]$. The vertical
lines there represent the $1\sigma$, $2\sigma$, $3\sigma$, and $4\sigma$
intervals around the peak of the posterior. The null hypothesis is discarded
with a spurious significance of $\sim 3.9\sigma$. The real significance is
obtained after the look-elsewhere effect is taken into account. Following the
steps in §6.4, we conclude that this effect dilutes this significance to
approximately $1.22\sigma$, which is far from allowing us to claim any
detection, regardless of the substantial Bayes factor.

In the lower panel of Figure 13 we show our results for step features. The
vanilla model of inflation is favoured for all frequency bins. The Planck
paper on inflation [120] quotes a best-fit with $\omega_{\rm step}\approx
1230\;\mathrm{Mpc}$, but we cannot reach frequencies close to this value with
current data.

The last model we comment on is the sound feature. The Bayes factors for this
model are $\ln B_{01}=0.620,0.886$, and $0.045$ for BOSS low-$z$, high-$z$,
and eBOSS Quasars, respectively. According to Jeffreys’ scale, these values
indicate that the model comparison is inconclusive – sound features are
neither discarded nor preferred by current LSS datasets. The best fits for
these same datasets lead to an improvement of
$\Delta\chi^{2}=-6.81,-8.76,-11.15$. Finally, in Figure 11, we show the
marginalized posterior for the characteristic time scale of this PF,
$\tau_{f}$, and its sharpness, $\ln\beta$ (see Eq. 3.14). A takeaway from this
figure is that the current LSS data have a series of preferred modes in the
posterior (darker horizontal contours). This is an improvement over previous
results using the WiggleZ data [93]. Currently, the most stringent constraints
for sound features are from the CMB [91, 92, 124], and since LSS alone can
currently restrict regions in the $\tau_{f}$ parameter, regardless of not
having strong constraints on the amplitude, a combined analysis would be
ideal. We leave this task for future work.

##  8 Summary and conclusions

Understanding the nature of inflation is one of the major open questions in
modern cosmology. To address this problem, many theories and frameworks were
developed (e.g. see [4]) over recent decades. These extended models usually
predict additional corrections to the power spectrum of curvature
perturbations beyond the results of single-field slow-roll inflation models.
Therefore, detecting any secondary signal lurking in the data would be a
window into new physics happening during the early universe. In this paper, we
obtained combined large-scale structure constraints on these additional
signals, called primordial features (PF), using the BOSS DR12 and eBOSS QSO
datasets. In [60] it was shown that LSS can lead to tighter constraints on PF
than Planck, so we expect our results to be close to the optimal constraints
that can be obtained with current cosmological datasets. They also showed that
next-generation galaxy surveys such as DESI will improve the sensitivity to
features by a factor of 7. Furthermore, the fact that the DESI quasar sample
will have a larger cosmic volume will also extend the frequency range of
linear features by a factor of 2.5.

We implemented a pipeline to search for PF signals considering an approximate
method for the non-linearities present in the late-time cosmological probes,
as discussed in Section 2. This framework proved to be versatile and robust
since it can be used to search for any oscillatory PF and is a simple
extension of a typical BAO analysis. We used it to search for four different
signals of PF, introduced in Section 3: two phenomenological oscillatory
templates (linearly and logarithmically spaced oscillations); PF due to the
presence of a step in the inflaton’s potential; and PF due to a variation in
the inflaton’s speed of sound. In Section 5, we presented the statistical
metrics we used to compute the significance of any potential detection and
applied these to mocks. This exercise revealed that the mocks we used are not
a perfect reproduction of the data – they prefer higher values of the non-
linear damping, $\Sigma_{\rm nl}$ (see Table 4 and Appendix B). We were able
to circumvent this problem when estimating the impact of the look-elsewhere
effect, as discussed in §6.4, which allowed us to estimate the global
significance of the events we see in the data.

For linear features, we scanned the range $\omega_{\rm
lin}\in[100,4000]\;\mathrm{Mpc}$, and our results are presented in Figure 8.
We showed that although BOSS DR12 has a higher number density than the eBOSS
QSO sample, the latter can lead to stronger constraints at high frequencies
($\omega_{\rm lin}>2500\;\mathrm{Mpc}$) due to the smaller fundamental mode in
eBOSS. For logarithmic features, we scanned the range $\omega_{\rm
log}\in[10,380]$, and the results are given in Figure 9. We found one event
with a Bayes factor of $\ln B_{01}\approx-4$, which apparently favours the
feature model relative to a featureless model. The null hypothesis for this
frequency bin is eliminated with a significance of $3.9\sigma$. However, when
the look-elsewhere effect is taken into account, the significance is reduced
to approximately $1.22\sigma$. The results for step features are given in
Figure 10, and we found no significant indication that the data prefer this
model with respect to the featureless inflationary model. In Figure 13 we
summarize the Bayes factor for all these three models as a function of their
frequency. Finally, for sound features, we found global Bayes factors of $\ln
B_{01}=0.620,0.886$ and $0.045$ for BOSS low-$z$, high-$z$, and eBOSS Quasars.
According to the Jeffreys’ scale, these models are thus neither favoured nor
discarded by current LSS datasets alone. Nevertheless, in Figure 11 we showed
that LSS constraints have preferred values for the typical scale of the speed
of sound reduction, $\tau_{f}$. The CMB constraints for this model also
present a set of preferred values of $\tau_{f}$, so combining both LSS and CMB
constraints will lead to a significant reduction in the allowed parameter
space.

With this work, we are one step closer to extracting the maximum possible
information about the primordial universe using large-scale structure, leading
to constraints competitive with (or, in some cases, surpassing) the CMB. Our
approach provides a baseline, which we will seek to improve before applying it
to the forthcoming DESI Y1 data. From the theoretical side, for instance, we
plan to search for other primordial feature models and use more realistic
modelling for the non-linear damping of features. Nevertheless, in all of
this, the most thrilling improvement will be in the data itself, with a number
density of objects being one order of magnitude larger, populating a volume
approximately five times larger [125, 126, 127, 128]. This will ultimately
lead to an improvement in the precision of PF amplitudes by a factor of
$\sim\,$7 (see Figure 10 of [60]) and, perhaps more importantly, allow us to
scan a frequency range $\sim\,$2.5 times larger, up to $\omega_{\rm lin}\sim
10\,000\,\mathrm{Mpc}$.

## Acknowledgments

This project has received funding from the European Research Council (ERC)
under the European Union’s Horizon 2020 research and innovation program (grant
agreement 853291). FB is a Royal Society University Research Fellow.

##  Appendix A Comparing the parametrisation

In the literature, there are results considering two distinct modelling
approaches of the oscillatory features. In this appendix, we scrutinize their
differences and how to correctly impose a non-informative prior.

In [120], the search for linear and logarithmic features is performed assuming

$\delta P=A_{\rm X}\sin(\omega_{\rm X}k+\phi)\;,$ (A.1)

where $X=\mathrm{lin,log}$. In others words, the parametrisation is

$\delta P=A^{\rm X}_{\sin}\sin(\omega_{X}k)+A^{\rm
X}_{\cos}\cos(\omega_{X}k)\;.$ (A.2)

The equations relating $(A_{\rm X},\phi)$ with $(A^{\rm X}_{\sin},A^{\rm
X}_{\cos})$ are

$\displaystyle A_{X}^{2}$ $\displaystyle=$ $\displaystyle(A^{\rm
X}_{\sin})^{2}+(A^{\rm X}_{\cos})^{2}\;,$ (A.3) $\displaystyle\tan{\phi}$
$\displaystyle=$ $\displaystyle\frac{A^{\rm X}_{\cos}}{A^{\rm X}_{\sin}}\;.$
(A.4)

Figure 15: upper panel: Comparing the 95% credible intervals on $A_{\rm lin}$
using $(A^{\rm lin}_{\rm sin},A^{\rm lin}_{\rm cos})$ and the $(A_{\rm
lin},\phi)$ parametrisation. These two parametrisations are connected by a
non-linear transformation, so to impose a uniform prior in one of them leads
to a non-uniform prior in the other. bottom panel: How the prior
$\mathcal{U}(-0.5,0.5)$ for $A^{\rm sin}_{\rm lin}$ and $A^{\rm cos}_{\rm
lin}$ are changed when transformed into the $(A_{\rm lin},\phi)$
parametrisation.

Figure 16: The influence of $\Sigma_{\rm nl}$ on linear feature constraints
from the eBOSS Quasar sample and from the BOSS high-$z$ galaxy sample. The
value of $\Sigma_{\rm nl}$ is fixed to the same value used on the official BAO
analysis, [66, 68].

Since the equations above are non-linear, one should be careful when assigning
the priors. A non-informative prior on the oscillations means that we have no
knowledge about its amplitude and phase, so the natural choice for a non-
informative prior would be, for instance, $p(A_{\rm X})=\mathcal{U}(-0.5,0.5)$
and $p(\phi)=\mathcal{U}(-\pi,\pi)$. According to the equations above, it is
evident that these priors will be non-uniform in the parametrisation $(A^{\rm
X}_{\sin},A^{\rm X}_{\cos})$. In [60], the authors decided to impose a uniform
prior on $(A^{\rm X}_{\sin},A^{\rm X}_{\cos})$, but according to the
discussion above they are non-uniform on $(A_{\rm X},\phi)$, and therefore
informative. In the upper panel of Figure 15 we show how a uniform prior on
$A_{\rm sin}$ and $A_{\rm cos}$ is transformed when mapped into the phase
parametrisation. To check the consistency of using an uniform prior on
$(A^{\rm X}_{\sin},A^{\rm X}_{\cos})$, we did the following. We first
performed the analysis using the BOSS data and assuming $p(A^{\rm\rm
lin}_{\sin},A^{\rm lin}_{\cos})=\mathcal{U}(-0.5,0.5)$. The 95% credible
intervals for $A^{\rm lin}$ are represented by the solid lines in Figure 15.
Next, we performed another analysis considering the phase parametrisation and
assuming $p(A_{\rm lin})=\mathcal{U}(-0.5,0.5)$ and
$p(\rm\phi)=\mathcal{U}(-\pi,\pi)$ and the results are represented by the
dashed lines in the same figure. The credible intervals for the uniform
parametrisation on $(A^{\rm X}_{\sin},A^{\rm X}_{\cos})$ are more stringent,
indicating that it’s prior is slightly informative. Finally, we transformed
the priors $p(A^{\rm\rm lin}_{\sin},A^{\rm lin}_{\cos})=\mathcal{U}(-0.5,0.5)$
into priors on $(A_{\rm lin},\phi)$ using Eqs. A.3 and A.4, and the results
are represented by the dotted lines in the bottom panel of Figure 15. We
conclude that the most conservative choice is to use the parametrisation
$(A_{\rm lin},\phi)$.

##  Appendix B The impact of non-linear clustering

In Section 6 we showed that the mock constraints on global oscillatory
features are noisier than the data. Here we show some further tests we
performed to understand what may be causing this.

The main purpose of the mock catalogues produced for the BOSS and eBOSS
datasets is to compute covariance matrices. Since performing high-resolution
simulations is computationally costly and the number of catalogues needed to
have reliable covariance matrix estimations is high, of the order of
$\sim\mathcal{O}(10^{3})$ mocks, it is necessary to simplify the process of
catalogue generation. The mocks we used in this paper were generated by
assuming approximate gravity solvers to obtain the large-scale dark matter
field [129, 108]. To obtain the BOSS Patchy mocks the dark-matter particles
were displaced from Lagrangian to Eulerian positions using second order
Lagrangian perturbation theory (2LPT) [101]. To obtain the EZmocks for eBOSS
[107], on the other hand, the effective Zel’dovich approximation was assumed
(i.e. linear Lagrangian perturbation theory). The tracers are introduced later
on top of the matter field following a biasing scheme, along with the
systematic effects inherent in each dataset 555Some of these systematic
effects include survey geometry and veto masks, radial selection, fibre
collision and redshift assignment failures (e.g. see [130, 131]).. Both
methods were widely tested and generated catalogues precise enough for
estimating the covariance matrix for data analysis in cosmology [132, 133,
134]. Nevertheless, as we pointed out before, the whole pipeline for mock
generation requires a series of approximations in both its theoretical and
processing aspect. For instance, 2LPT does not completely capture the cosmic
web dynamics up to the scales used to generate the models. Moreover, the
resolutions of the density fields are low ($\sim 5h^{-1}\;\mathrm{Mpc}$ for
eBOSS EZmocks and $\sim 2.5h^{-1}\;\mathrm{Mpc}$ for BOSS Patchy mocks),
meaning that below these scales the tracers are positioned following a cloud-
in-cell distribution. These and other limitations entail catalogues that
although good enough for estimating the covariance matrix, do not truthfully
reproduce all ingredients present in the data (see [135, 101, 107]). One
consequence is that the BAO signal in these mocks is weaker than in data [68,
66]. Another consequence, and this is the most relevant for constraining
primordial features, is the fact these mocks also favour larger non-linear
damping, $\Sigma_{\rm nl}$, due to the smooth galaxy distribution inside grid
cells. This is emphasized in Figure 16. We show the results for linear
features, the isotropic BAO parameter, $\alpha_{\rm iso}$, and the non-linear
damping parameter for both eBOSS and BOSS high-$z$ products as a function of
the feature frequency. The solid green line represent the mean of the mocks
with varying $\Sigma_{\rm nl}$, whereas the solid red line is the same but
with the damping parameter fixed to the value used in the official BAO
analysis of each dataset. For the amplitude of linear feature, fixing
$\Sigma_{\rm nl}$ improves the mock constraint by a factor of $\sim 1.5$,
leading to the same error as predicted by the data. The impact of fixing
$\Sigma_{\rm nl}$ when analysing the data is not important because the value
we obtain by varying it is the same as the value we fix. For both datasets, on
the right, we show the scatter of $\Sigma_{\rm nl}$ for different mocks and
their mean. Notice that in both cases none of the mocks have a $\Sigma_{\rm
nl}$ compatible with the data.

We finish this section by emphasizing that having mock catalogues with
features in the primordial power spectrum is crucial for learning how these
signals can be measured in realistic scenarios and to understand how
systematic effects can interfere with this process. Our findings show that the
usual approximate methods to obtain numerically cheap mock catalogues may not
be ideal to accomplish this, which requires more precise gravity solvers or
full N-body simulations.

## References

  * [1] N. Aghanim, Y. Akrami, M. Ashdown, J. Aumont, C. Baccigalupi, M. Ballardini et al., _Planck 2018 results-vi. cosmological parameters_ , _Astronomy & Astrophysics_ 641 (2020) A6.
  * [2] Y. Akrami, F. Arroja, M. Ashdown, J. Aumont, C. Baccigalupi, M. Ballardini et al., _Planck 2018 results-x. constraints on inflation_ , _Astronomy & Astrophysics_ 641 (2020) A10.
  * [3] A. Achúcarro et al., _Inflation: Theory and Observations_ , 2203.08128.
  * [4] J. Chluba, J. Hamann and S. P. Patil, _Features and New Physical Scales in Primordial Observables: Theory and Observation_ , _Int. J. Mod. Phys. D_ 24 (2015) 1530023, [1505.01834].
  * [5] A. Vilenkin and L. H. Ford, _Gravitational Effects upon Cosmological Phase Transitions_ , _Phys. Rev. D_ 26 (1982) 1231.
  * [6] J. A. Adams, G. G. Ross and S. Sarkar, _Multiple inflation_ , _Nucl. Phys. B_ 503 (1997) 405–425, [hep-ph/9704286].
  * [7] J. A. Adams, B. Cresswell and R. Easther, _Inflationary perturbations from a potential with a step_ , _Phys. Rev. D_ 64 (2001) 123514, [astro-ph/0102236].
  * [8] V. Miranda and W. Hu, _Inflationary steps in the planck data_ , _Physical Review D_ 89 (2014) 083529.
  * [9] C. Dvorkin and W. Hu, _Generalized Slow Roll for Large Power Spectrum Features_ , _Phys. Rev. D_ 81 (2010) 023518, [0910.2237].
  * [10] W. Hu, _Generalized Slow Roll for Non-Canonical Kinetic Terms_ , _Phys. Rev. D_ 84 (2011) 027303, [1104.4500].
  * [11] P. Adshead, C. Dvorkin, W. Hu and E. A. Lim, _Non-Gaussianity from Step Features in the Inflationary Potential_ , _Phys. Rev. D_ 85 (2012) 023531, [1110.3050].
  * [12] N. Barnaby, Z. Huang, L. Kofman and D. Pogosyan, _Cosmological Fluctuations from Infra-Red Cascading During Inflation_ , _Phys. Rev. D_ 80 (2009) 043501, [0902.0615].
  * [13] D. J. H. Chung, E. W. Kolb, A. Riotto and I. I. Tkachev, _Probing Planckian physics: Resonant production of particles during inflation and features in the primordial power spectrum_ , _Phys. Rev. D_ 62 (2000) 043508, [hep-ph/9910437].
  * [14] G. Dvali, A. Gruzinov and M. Zaldarriaga, _Cosmological perturbations from inhomogeneous reheating, freezeout, and mass domination_ , _Phys. Rev. D_ 69 (2004) 083505, [astro-ph/0305548].
  * [15] L. Kofman, A. D. Linde, X. Liu, A. Maloney, L. McAllister and E. Silverstein, _Beauty is attractive: Moduli trapping at enhanced symmetry points_ , _JHEP_ 05 (2004) 030, [hep-th/0403001].
  * [16] A. E. Romano and M. Sasaki, _Effects of particle production during inflation_ , _Phys. Rev. D_ 78 (2008) 103522, [0809.5142].
  * [17] D. Green, B. Horn, L. Senatore and E. Silverstein, _Trapped Inflation_ , _Phys. Rev. D_ 80 (2009) 063533, [0902.1006].
  * [18] L. Kofman and D. Pogosyan, _Nonflat perturbations in inflationary cosmology_ , _Physics Letters B_ 214 (1988) 508–514.
  * [19] M. Braglia, D. K. Hazra, L. Sriramkumar and F. Finelli, _Generating primordial features at large scales in two field models of inflation_ , _JCAP_ 08 (2020) 025, [2004.00672].
  * [20] L. C. Price, J. Frazer, J. Xu, H. V. Peiris and R. Easther, _MultiModeCode: An efficient numerical solver for multifield inflation_ , _JCAP_ 03 (2015) 005, [1410.0685].
  * [21] X. Gao, D. Langlois and S. Mizuno, _Influence of heavy modes on perturbations in multiple field inflation_ , _JCAP_ 10 (2012) 040, [1205.5275].
  * [22] A. Achucarro, J.-O. Gong, S. Hardeman, G. A. Palma and S. P. Patil, _Features of heavy physics in the CMB power spectrum_ , _JCAP_ 01 (2011) 030, [1010.3693].
  * [23] U. H. Danielsson, _A Note on inflation and transPlanckian physics_ , _Phys. Rev. D_ 66 (2002) 023511, [hep-th/0203198].
  * [24] J. Martin and R. Brandenberger, _On the dependence of the spectra of fluctuations in inflationary cosmology on transPlanckian physics_ , _Phys. Rev. D_ 68 (2003) 063513, [hep-th/0305161].
  * [25] V. Bozza, M. Giovannini and G. Veneziano, _Cosmological perturbations from a new physics hypersurface_ , _JCAP_ 05 (2003) 001, [hep-th/0302184].
  * [26] X. Chen, _Primordial Features as Evidence for Inflation_ , _JCAP_ 01 (2012) 038, [1104.1323].
  * [27] J. R. Fergusson, H. F. Gruetjen, E. P. S. Shellard and B. Wallisch, _Polyspectra searches for sharp oscillatory features in cosmic microwave sky data_ , _Phys. Rev. D_ 91 (2015) 123506, [1412.6152].
  * [28] J. R. Fergusson, H. F. Gruetjen, E. P. S. Shellard and M. Liguori, _Combining power spectrum and bispectrum measurements to detect oscillatory features_ , _Phys. Rev. D_ 91 (2015) 023502, [1410.5114].
  * [29] K. Freese, J. A. Frieman and A. V. Olinto, _Natural inflation with pseudo \- Nambu-Goldstone bosons_ , _Phys. Rev. Lett._ 65 (1990) 3233–3236.
  * [30] L. McAllister, E. Silverstein and A. Westphal, _Gravity Waves and Linear Inflation from Axion Monodromy_ , _Phys. Rev. D_ 82 (2010) 046003, [0808.0706].
  * [31] M. Kawasaki and K. Nakayama, _Axions: Theory and Cosmological Role_ , _Ann. Rev. Nucl. Part. Sci._ 63 (2013) 69–95, [1301.1123].
  * [32] J. Hamann, A. Shafieloo and T. Souradeep, _Features in the primordial power spectrum? A frequentist analysis_ , _Journal of Cosmology and Astroparticle Physics_ 2010 (Apr., 2010) 010, [0912.2728].
  * [33] A. Shafieloo, T. Souradeep, P. Manimaran, P. K. Panigrahi and R. Rangarajan, _Features in the Primordial Spectrum from WMAP: A Wavelet Analysis_ , _Phys. Rev. D_ 75 (2007) 123502, [astro-ph/0611352].
  * [34] P. D. Meerburg, R. A. M. J. Wijers and J. P. van der Schaar, _WMAP7 constraints on oscillations in the primordial power spectrum_ , _Mon. Not. R. Astron. Soc._ 421 (Mar., 2012) 369–380, [1109.5264].
  * [35] J. Martin and C. Ringeval, _Superimposed oscillations in the WMAP data?_ , _Phys. Rev. D_ 69 (2004) 083515, [astro-ph/0310382].
  * [36] J. Martin and C. Ringeval, _Exploring the superimposed oscillations parameter space_ , _JCAP_ 01 (2005) 007, [hep-ph/0405249].
  * [37] R. Flauger, L. McAllister, E. Pajer, A. Westphal and G. Xu, _Oscillations in the CMB from Axion Monodromy Inflation_ , _JCAP_ 06 (2010) 009, [0907.2916].
  * [38] M. Aich, D. K. Hazra, L. Sriramkumar and T. Souradeep, _Oscillations in the inflaton potential: Complete numerical treatment and comparison with the recent and forthcoming CMB datasets_ , _Phys. Rev. D_ 87 (2013) 083526, [1106.2798].
  * [39] H. Peiris, R. Easther and R. Flauger, _Constraining Monodromy Inflation_ , _JCAP_ 09 (2013) 018, [1303.2616].
  * [40] R. Easther and R. Flauger, _Planck Constraints on Monodromy Inflation_ , _JCAP_ 02 (2014) 037, [1308.3736].
  * [41] V. Miranda and W. Hu, _Inflationary Steps in the Planck Data_ , _Phys. Rev. D_ 89 (2014) 083529, [1312.0946].
  * [42] M. Braglia, X. Chen and D. K. Hazra, _Primordial standard clock models and CMB residual anomalies_ , _Phys. Rev. D_ 105 (2022) 103523, [2108.10110].
  * [43] P. D. Meerburg and D. N. Spergel, _Searching for oscillations in the primordial power spectrum. II. Constraints from Planck data_ , _Phys. Rev. D_ 89 (2014) 063537, [1308.3705].
  * [44] X. Chen and M. H. Namjoo, _Standard Clock in Primordial Density Perturbations and Cosmic Microwave Background_ , _Phys. Lett. B_ 739 (2014) 285–292, [1404.1536].
  * [45] D. Chandra and S. Pal, _Investigating the constraints on primordial features with future cosmic microwave background and galaxy surveys_ , _Journal of Cosmology and Astroparticle Physics_ 2022 (Sept., 2022) 024, [2205.01164].
  * [46] J. Hamann and A. Malhotra, _Constraining primordial tensor features with the anisotropies of the Cosmic Microwave Background_ , 2209.00827.
  * [47] M. Braglia, X. Chen, D. K. Hazra and L. Pinol, _Back to the features: assessing the discriminating power of future CMB missions on inflationary models_ , 2210.07028.
  * [48] J. Fumagalli, G. A. Palma, S. Renaux-Petel, S. Sypsas, L. T. Witkowski and C. Zenteno, _Primordial gravitational waves from excited states_ , _Journal of High Energy Physics_ 2022 (Mar., 2022) 196, [2111.14664].
  * [49] J. Fumagalli, S. Renaux-Petel and L. T. Witkowski, _Oscillations in the stochastic gravitational wave background from sharp features and particle production during inflation_ , _JCAP_ 08 (2021) 030, [2012.02761].
  * [50] M. Braglia, X. Chen and D. K. Hazra, _Probing Primordial Features with the Stochastic Gravitational Wave Background_ , _JCAP_ 03 (2021) 005, [2012.05821].
  * [51] J. Fumagalli, S. e. Renaux-Petel and L. T. Witkowski, _Resonant features in the stochastic gravitational wave background_ , _JCAP_ 08 (2021) 059, [2105.06481].
  * [52] L. T. Witkowski, G. Domènech, J. Fumagalli and S. Renaux-Petel, _Expansion history-dependent oscillations in the scalar-induced gravitational wave background_ , _JCAP_ 05 (2022) 028, [2110.09480].
  * [53] J. Fumagalli, G. A. Palma, S. Renaux-Petel, S. Sypsas, L. T. Witkowski and C. Zenteno, _Primordial gravitational waves from excited states_ , _JHEP_ 03 (2022) 196, [2111.14664].
  * [54] J. Fumagalli, M. Pieroni, S. Renaux-Petel and L. T. Witkowski, _Detecting primordial features with LISA_ , _JCAP_ 07 (2022) 020, [2112.06903].
  * [55] X. Chen, _Fingerprints of primordial universe paradigms as features in density perturbations_ , _Physics Letters B_ 706 (Dec., 2011) 111–115, [1106.1635].
  * [56] X. Chen and C. Ringeval, _Searching for Standard Clocks in the Primordial Universe_ , _JCAP_ 08 (2012) 014, [1205.6085].
  * [57] M. Crocce and R. Scoccimarro, _Nonlinear Evolution of Baryon Acoustic Oscillations_ , _Phys. Rev. D_ 77 (2008) 023533, [0704.2783].
  * [58] Z. Vlah, U. Seljak, M. Y. Chu and Y. Feng, _Perturbation theory, effective field theory, and oscillations in the power spectrum_ , _JCAP_ 03 (2016) 057, [1509.02120].
  * [59] A. Vasudevan, M. M. Ivanov, S. Sibiryakov and J. Lesgourgues, _Time-sliced perturbation theory with primordial non-Gaussianity and effects of large bulk flows on inflationary oscillating features_ , _JCAP_ 09 (2019) 037, [1906.08697].
  * [60] F. Beutler, M. Biagetti, D. Green, A. Slosar and B. Wallisch, _Primordial features from linear to nonlinear scales_ , _Physical Review Research_ 1 (2019) 033209.
  * [61] S.-F. Chen, Z. Vlah and M. White, _Modeling features in the redshift-space halo power spectrum with perturbation theory_ , _JCAP_ 11 (2020) 035, [2007.00704].
  * [62] M. Ballardini, R. Murgia, M. Baldi, F. Finelli and M. Viel, _Non-linear damping of superimposed primordial oscillations on the matter power spectrum in galaxy surveys_ , _JCAP_ 04 (2020) 030, [1912.12499].
  * [63] Y. Shi, M. Cautun and B. Li, _New method for initial density reconstruction_ , _Phys. Rev. D_ 97 (2018) 023505, [1709.06350].
  * [64] A. Ota, H.-J. Seo, S. Saito and F. Beutler, _Modeling iterative reconstruction and displacement field in the large scale structure_ , _Phys. Rev. D_ 104 (2021) 123508, [2106.00146].
  * [65] H.-J. Seo, A. Ota, M. Schmittfull, S. Saito and F. Beutler, _Iterative reconstruction excursions for Baryon Acoustic Oscillations and beyond_ , _Mon. Not. Roy. Astron. Soc._ 511 (2022) 1557–1573, [2106.00530].
  * [66] R. Neveux et al., _The completed SDSS-IV extended Baryon Oscillation Spectroscopic Survey: BAO and RSD measurements from the anisotropic power spectrum of the quasar sample between redshift 0.8 and 2.2_ , _Mon. Not. Roy. Astron. Soc._ 499 (2020) 210–229, [2007.08999].
  * [67] H.-J. Seo, F. Beutler, A. J. Ross and S. Saito, _Modeling the reconstructed BAO in Fourier space_ , _Mon. Not. Roy. Astron. Soc._ 460 (2016) 2453–2471, [1511.00663].
  * [68] BOSS collaboration, F. Beutler et al., _The clustering of galaxies in the completed SDSS-III Baryon Oscillation Spectroscopic Survey: baryon acoustic oscillations in the Fourier space_ , _Mon. Not. Roy. Astron. Soc._ 464 (2017) 3409–3430, [1607.03149].
  * [69] W. E. Ballinger, J. A. Peacock and A. F. Heavens, _Measuring the cosmological constant with redshift surveys_ , _Mon. Not. Roy. Astron. Soc._ 282 (1996) 877–888, [astro-ph/9605017].
  * [70] H. Magira, Y. P. Jing and Y. Suto, _Cosmological redshift-space distortion on clustering of high-redshift objects: correction for nonlinear effects in power spectrum and tests with n-body simulations_ , _Astrophys. J._ 528 (2000) 30, [astro-ph/9907438].
  * [71] C. Park, M. S. Vogeley, M. J. Geller and J. P. Huchra, _Power spectrum, correlation function, and tests for luminosity bias in the CfA redshift survey_ , _The Astrophysical Journal_ 431 (Aug., 1994) 569.
  * [72] J. A. Peacock and S. J. Dodds, _Reconstructing the linear power spectrum of cosmological mass fluctuations_ , _Mon. Not. Roy. Astron. Soc._ 267 (1994) 1020–1034, [astro-ph/9311057].
  * [73] L. Senatore and M. Zaldarriaga, _The IR-resummed Effective Field Theory of Large Scale Structures_ , _JCAP_ 02 (2015) 013, [1404.5954].
  * [74] D. Baumann and L. McAllister, _Inflation and string theory_.  Cambridge University Press, 2015.
  * [75] R. Flauger, L. McAllister, E. Silverstein and A. Westphal, _Drifting oscillations in axion monodromy_ , _Journal of Cosmology and Astroparticle Physics_ 2017 (2017) 055.
  * [76] X. Chen, R. Easther and E. A. Lim, _Generation and Characterization of Large Non-Gaussianities in Single Field Inflation_ , _JCAP_ 04 (2008) 010, [0801.3295].
  * [77] R. Flauger and E. Pajer, _Resonant Non-Gaussianity_ , _JCAP_ 01 (2011) 017, [1002.0833].
  * [78] X. Chen, _Folded Resonant Non-Gaussianity in General Single Field Inflation_ , _JCAP_ 12 (2010) 003, [1008.2485].
  * [79] X. Chen, R. Easther and E. A. Lim, _Large Non-Gaussianities in Single Field Inflation_ , _JCAP_ 06 (2007) 023, [astro-ph/0611645].
  * [80] G. A. Palma, D. Sapone and S. Sypsas, _Constraints on inflation with LSS surveys: features in the primordial power spectrum_ , _JCAP_ 06 (2018) 004, [1710.02570].
  * [81] M. Ballardini, F. Finelli, F. Marulli, L. Moscardini and A. Veropalumbo, _New constraints on primordial features from the galaxy two-point correlation function_ , 2202.08819.
  * [82] J. Adams, B. Cresswell and R. Easther, _Inflationary perturbations from a potential with a step_ , _Physical Review D_ 64 (2001) 123514.
  * [83] G. A. Palma, D. Sapone and S. Sypsas, _Constraints on inflation with lss surveys: features in the primordial power spectrum_ , _Journal of Cosmology and Astroparticle Physics_ 2018 (2018) 004.
  * [84] C. P. Burgess, M. W. Horbatsch and S. P. Patil, _Inflating in a Trough: Single-Field Effective Theory from Multiple-Field Curved Valleys_ , _JHEP_ 01 (2013) 133, [1209.5701].
  * [85] A. Achucarro, J.-O. Gong, S. Hardeman, G. A. Palma and S. P. Patil, _Effective theories of single field inflation when heavy fields matter_ , _JHEP_ 05 (2012) 066, [1201.6342].
  * [86] A. Achucarro, V. Atal, S. Cespedes, J.-O. Gong, G. A. Palma and S. P. Patil, _Heavy fields, reduced speeds of sound and decoupling during inflation_ , _Phys. Rev. D_ 86 (2012) 121301, [1205.0710].
  * [87] J.-O. Gong, K. Schalm and G. Shiu, _Correlating correlation functions of primordial perturbations_ , _Phys. Rev. D_ 89 (2014) 063540, [1401.4402].
  * [88] Y.-F. Cai, W. Zhao and Y. Zhang, _Cosmic microwave background power asymmetry from primordial sound speed parameter_ , _Phys. Rev. D_ 89 (2014) 023005, [1307.4090].
  * [89] N. Bartolo, D. Cannone and S. Matarrese, _The Effective Field Theory of Inflation Models with Sharp Features_ , _JCAP_ 10 (2013) 038, [1307.3483].
  * [90] A. Achucarro, J.-O. Gong, S. Hardeman, G. A. Palma and S. P. Patil, _Mass hierarchies and non-decoupling in multi-scalar field dynamics_ , _Phys. Rev. D_ 84 (2011) 043502, [1005.3848].
  * [91] A. Achucarro, V. Atal, B. Hu, P. Ortiz and J. Torrado, _Inflation with moderately sharp features in the speed of sound: Generalized slow roll and in-in formalism for power spectrum and bispectrum_ , _Phys. Rev. D_ 90 (2014) 023511, [1404.7522].
  * [92] A. Achúcarro, V. Atal, P. Ortiz and J. Torrado, _Localized correlated features in the CMB power spectrum and primordial bispectrum from a transient reduction in the speed of sound_ , _Phys. Rev. D_ 89 (2014) 103006, [1311.2552].
  * [93] B. Hu and J. Torrado, _Searching for primordial localized features with CMB and LSS spectra_ , _Phys. Rev. D_ 91 (2015) 064039, [1410.4804].
  * [94] D. Cannone, N. Bartolo and S. Matarrese, _Perturbative Unitarity of Inflationary Models with Features_ , _Phys. Rev. D_ 89 (2014) 127301, [1402.2258].
  * [95] K. S. Dawson, D. J. Schlegel, C. P. Ahn, S. F. Anderson, É. Aubourg, S. Bailey et al., _The baryon oscillation spectroscopic survey of sdss-iii_ , _The Astronomical Journal_ 145 (2012) 10.
  * [96] D. J. Eisenstein, D. H. Weinberg, E. Agol, H. Aihara, C. A. Prieto, S. F. Anderson et al., _Sdss-iii: Massive spectroscopic surveys of the distant universe, the milky way, and extra-solar planetary systems_ , _The Astronomical Journal_ 142 (2011) 72.
  * [97] S. A. Smee, J. E. Gunn, A. Uomoto, N. Roe, D. Schlegel, C. M. Rockosi et al., _The multi-object, fiber-fed spectrographs for the sloan digital sky survey and the baryon oscillation spectroscopic survey_ , _The Astronomical Journal_ 146 (2013) 32.
  * [98] SDSS collaboration, J. E. Gunn et al., _The 2.5 m Telescope of the Sloan Digital Sky Survey_ , _Astron. J._ 131 (2006) 2332–2359, [astro-ph/0602326].
  * [99] B. Reid et al., _SDSS-III Baryon Oscillation Spectroscopic Survey Data Release 12: galaxy target selection and large scale structure catalogues_ , _Mon. Not. Roy. Astron. Soc._ 455 (2016) 1553–1573, [1509.06529].
  * [100] F. Beutler and P. McDonald, _Unified galaxy power spectrum measurements from 6dFGS, BOSS, and eBOSS_ , _JCAP_ 11 (2021) 031, [2106.06324].
  * [101] F.-S. Kitaura et al., _The clustering of galaxies in the SDSS-III Baryon Oscillation Spectroscopic Survey: mock galaxy catalogues for the BOSS Final Data Release_ , _Mon. Not. Roy. Astron. Soc._ 456 (2016) 4156–4173, [1509.06400].
  * [102] D. J. Eisenstein, H.-j. Seo, E. Sirko and D. Spergel, _Improving Cosmological Distance Measurements by Reconstruction of the Baryon Acoustic Peak_ , _Astrophys. J._ 664 (2007) 675–679, [astro-ph/0604362].
  * [103] Y. Li, H.-M. Zhu and B. Li, _Non-linear reconstruction of features in the primordial power spectrum from large-scale structure_ , _Mon. Not. Roy. Astron. Soc._ 514 (2022) 4363–4378, [2102.09007].
  * [104] K. S. Dawson et al., _The SDSS-IV extended Baryon Oscillation Spectroscopic Survey: Overview and Early Data_ , _Astron. J._ 151 (2016) 44, [1508.04473].
  * [105] A. J. Ross et al., _The Completed SDSS-IV extended Baryon Oscillation Spectroscopic Survey: Large-scale structure catalogues for cosmological analysis_ , _Mon. Not. Roy. Astron. Soc._ 498 (2020) 2354–2371, [2007.09000].
  * [106] M. R. Blanton, M. A. Bershady, B. Abolfathi, F. D. Albareti, C. A. Prieto, A. Almeida et al., _Sloan digital sky survey iv: Mapping the milky way, nearby galaxies, and the distant universe_ , _The Astronomical Journal_ 154 (2017) 28.
  * [107] C. Zhao et al., _The completed SDSS-IV extended Baryon Oscillation Spectroscopic Survey: 1000 multi-tracer mock catalogues with redshift evolution and systematics for galaxies and quasars of the final data release_ , _Mon. Not. Roy. Astron. Soc._ 503 (2021) 1149–1173, [2007.08997].
  * [108] C.-H. Chuang, F.-S. Kitaura, F. Prada, C. Zhao and G. Yepes, _EZmocks: extending the Zel’dovich approximation to generate mock galaxy catalogues with accurate clustering statistics_ , _Mon. Not. Roy. Astron. Soc._ 446 (2015) 2621–2628, [1409.1124].
  * [109] M. J. Wilson, J. A. Peacock, A. N. Taylor and S. de la Torre, _Rapid modelling of the redshift-space power spectrum multipoles for a masked density field_ , _Mon. Not. Roy. Astron. Soc._ 464 (2017) 3121–3130, [1511.07799].
  * [110] R. Trotta, _Bayesian Methods in Cosmology_ , 1, 2017.  1701.01467.
  * [111] R. Trotta, _Applications of Bayesian model selection to cosmological parameters_ , _Mon. Not. Roy. Astron. Soc._ 378 (2007) 72–82, [astro-ph/0504022].
  * [112] C. R. Jenkins and J. A. Peacock, _The power of Bayesian evidence in astronomy_ , _Mon. Not. Roy. Astron. Soc._ 413 (2011) 2895, [1101.4822].
  * [113] J. M. Dickey, _The weighted likelihood ratio, linear hypotheses on normal location parameters_ , _Annals of Mathematical Statistics_ 42 (1971) 204–223.
  * [114] D. Foreman-Mackey, D. W. Hogg, D. Lang and J. Goodman, _emcee: The mcmc hammer_ , _PASP_ 125 (2013) 306–312, [1202.3665].
  * [115] A. Lewis, _GetDist: a Python package for analysing Monte Carlo samples_ , 1910.13970.
  * [116] M. Karamanis, F. Beutler, J. A. Peacock, D. Nabergoj and U. Seljak, _Accelerating astronomical and cosmological inference with preconditioned Monte Carlo_ , _Mon. Not. Roy. Astron. Soc._ 516 (2022) 1644–1653, [2207.05652].
  * [117] M. Karamanis, D. Nabergoj, F. Beutler, J. A. Peacock and U. Seljak, _pocoMC: A Python package for accelerated Bayesian inference in astronomy and cosmology_ , 2207.05660.
  * [118] H. Dembinski and P. O. et al., _scikit-hep/iminuit_ , .
  * [119] F. James and M. Roos, _Minuit: A System for Function Minimization and Analysis of the Parameter Errors and Correlations_ , _Comput. Phys. Commun._ 10 (1975) 343–367.
  * [120] Planck Collaboration, Y. Akrami, F. Arroja, M. Ashdown, J. Aumont, C. Baccigalupi et al., _Planck 2018 results. X. Constraints on inflation_ , _Astron. Astrophys._ 641 (Sept., 2020) A10, [1807.06211].
  * [121] E. Gross and O. Vitells, _Trial factors for the look elsewhere effect in high energy physics_ , _Eur. Phys. J. C_ 70 (2010) 525–530, [1005.1891].
  * [122] A. E. Bayer and U. Seljak, _The look-elsewhere effect from a unified Bayesian and frequentist perspective_ , _JCAP_ 10 (2020) 009, [2007.13821].
  * [123] H. Jeffreys, _The theory of probability_.  OUP Oxford, 1998.
  * [124] G. Cañas Herrera, J. Torrado and A. Achúcarro, _Bayesian reconstruction of the inflaton’s speed of sound using CMB data_ , _Phys. Rev. D_ 103 (2021) 123531, [2012.04640].
  * [125] DESI collaboration, R. Zhou et al., _Target Selection and Validation of DESI Luminous Red Galaxies_ , _Astron. J._ 165 (2023) 58, [2208.08515].
  * [126] C. Hahn et al., _DESI Bright Galaxy Survey: Final Target Selection, Design, and Validation_ , 2208.08512.
  * [127] A. Raichoor et al., _Target Selection and Validation of DESI Emission Line Galaxies_ , _Astron. J._ 165 (2023) 126, [2208.08513].
  * [128] E. Chaussidon et al., _Target Selection and Validation of DESI Quasars_ , _Astrophys. J._ 944 (2023) 107, [2208.08511].
  * [129] F.-S. Kitaura, G. Yepes and F. Prada, _Modelling Baryon Acoustic Oscillations with Perturbation Theory and Stochastic Halo Biasing_ , _Mon. Not. Roy. Astron. Soc._ 439 (2014) 21, [1307.3285].
  * [130] BOSS collaboration, A. J. Ross et al., _The clustering of galaxies in the completed SDSS-III Baryon Oscillation Spectroscopic Survey: Observational systematics and baryon acoustic oscillations in the correlation function_ , _Mon. Not. Roy. Astron. Soc._ 464 (2017) 1168–1191, [1607.03145].
  * [131] P. Laurent et al., _Clustering of quasars in SDSS-IV eBOSS : study of potential systematics and bias determination_ , _JCAP_ 07 (2017) 017, [1705.04718].
  * [132] C.-H. Chuang et al., _nIFTy Cosmology: Galaxy/halo mock catalogue comparison project on clustering statistics_ , _Mon. Not. Roy. Astron. Soc._ 452 (2015) 686–700, [1412.7729].
  * [133] L. Blot et al., _Comparing approximate methods for mock catalogues and covariance matrices II: Power spectrum multipoles_ , _Mon. Not. Roy. Astron. Soc._ 485 (2019) 2806–2824, [1806.09497].
  * [134] M. Lippich et al., _Comparing approximate methods for mock catalogues and covariance matrices – I. Correlation function_ , _Mon. Not. Roy. Astron. Soc._ 482 (2019) 1786–1806, [1806.09477].
  * [135] M. Pellejero-Ibañez, A. Balaguera-Antolínez, F.-S. Kitaura, R. E. Angulo, G. Yepes, C.-H. Chuang et al., _The bias of dark matter tracers: assessing the accuracy of mapping techniques_ , _Mon. Not. Roy. Astron. Soc._ 493 (2020) 586–593, [1910.13164].