LiteBIRD science goals and forecasts. A case study of the origin of primordial gravitational waves using large-scale CMB polarization

We study the possibility of using the LiteBIRD satellite B-mode survey to constrain models of inflation producing specific features in CMB angular power spectra. We explore a particular model example, i.e. spectator axion-SU(2) gauge field inflation. This model can source parity-violating gravitational waves from the amplification of gauge field fluctuations driven by a pseudoscalar “axionlike” field, rolling for a few e-folds during inflation. The sourced gravitational waves can exceed the vacuum contribution at reionization bump scales by about an order of magnitude and can be comparable to the vacuum contribution at recombination bump scales. We argue that a satellite mission with full sky coverage and access to the reionization bump scales is necessary to understand the origin of the primordial gravitational wave signal and distinguish among two production mechanisms: quantum vacuum fluctuations of spacetime and matter sources during inflation. We present the expected constraints on model parameters from LiteBIRD satellite simulations, which complement and expand previous studies in the literature. We find that LiteBIRD will be able to exclude with high significance standard single-field slow-roll models, such as the Starobinsky model, if the true model is the axion-SU(2) model with a feature at CMB scales. We further investigate the possibility of using the parity-violating signature of the model, such as the TB and EB angular power spectra, to disentangle it from the standard single-field slow-roll scenario. We find that most of the discriminating power of LiteBIRD will reside in BB angular power spectra rather than in TB and EB correlations.


Introduction
A stochastic background of primordial gravitational waves (hereafter GWs) predicted by the inflationary paradigm [1,2] represents one of the main targets of ongoing experimental efforts in cosmology.The simplest models of inflation (realized by a single scalar field minimally coupled to gravity and slowly rolling down its potential) predicts several properties of the distributions of the scalar (density) fluctuations [3][4][5][6][7], which are in remarkable agreement with cosmological observations [8][9][10][11].Additionally, if the GW background is detected, it would provide definitive evidence for cosmic inflation [12][13][14][15].
The spectrum of inflationary gravitational waves (tensor modes) extends over about 21 decades in frequency and is measurable through several different means.Among them, the primordial B-modes of the cosmic microwave background (CMB) polarization represent the most promising probe to detect the inflationary stochastic GW background [16,17], and also the closest in time.Other possible options include pulsar timing arrays (hereafter PTA) and laser interferometers (see, e.g.Ref. [18] for a review).Even though recent PTA measurements presented strong evidence for the existence of a stochastic GW background [19][20][21][22][23][24], this signal is still compatible with an astrophysical origin, i.e., due to inspirals of supermassive black hole binaries [25][26][27].Nonetheless, the deviation of the observed signal from the expected astrophysical GW background could still have a primordial origin [28][29][30].
Whereas there has been no detection of the B modes from primordial tensor modes yet, CMB experiments currently hold the tightest constraints on their amplitude, customarily parametrized through the tensor-to-scalar ratio r parameter, i.e. the ratio of the amplitudes of the tensor and scalar primordial power spectra.An upper limit r < 0.036 at 95 % C.L., established by the BICEP/Keck collaboration at the pivot scale of k 0 = 0.05 Mpc −1 assuming a fixed cosmology [31], was shown to increase to r < 0.042 at 95 % C.L. when fitting also for ΛCDM parameters [32].With the addition of Planck PR3 data, the current upper limit reduces to r < 0.035 at 95 % C.L. [31,33].An even tighter limit is obtained using Planck PR4 data and likelihoods, together with BAO data, i.e. r < 0.032 at 95 % C.L. [32].Using a conditioned covariance matrix as advocated in Ref. [34], the upper limit from the same datasets was shown to increase to r < 0.037 at 95 % C.L. for a profile likelihood approach [35] and to r < 0.038 at 95 % C.L. for a Monte Carlo Markov chain approach [34].Finally, fitting also the slope of the tensor power spectrum (and adding the LIGO-Virgo KAGRA dataset to the previous ones), leads to an improved upper bound r < 0.028 at 95 % C.L. at the pivot scale k 0 = 0.01 Mpc −1 [36].Given the importance of this measurement, several B-mode experiments, such as the ground-based Simons Array [37], Simons Observatory (SO) [38], the South Pole Observatory [39], the Cosmology Large Angular Scale Surveyor (CLASS) [40] and CMB-S4 [41], the balloon-borne SPIDER [42] and the LiteBIRD (Lite satellite for the study of B-mode polarization and Inflation from cosmic background Radiation Detection) satellite [43], are currently targeting this very faint primordial signal or will do so towards the end of this decade.More specifically, LiteBIRD is a strategic large-class mission selected by the Japan Aerospace Exploration Agency (JAXA) to be launched in the early 2030s.It will orbit the Sun-Earth Lagrangian point L2 and will map the CMB polarization over the entire sky for three years, using three telescopes in 15 frequency bands between 34 and 448 GHz.LiteBIRD, being a full-sky satellite mission, will also have access to the very largest scale B modes produced during cosmic reionization, in addition to the smaller scale B modes produced during cosmic recombination.Instead, planned ground-based CMB experiments will only be able to access the recombination signature at smaller scales.With a detection limit at the level of r ≲ 10 −3 , LiteBIRD and CMB-S4 are the most sensitive among the experiments that have already passed the proposal stage.
However, simply detecting r on its own would not allow one to understand the origin of the primordial GW background.In the simplest scenario (i.e. standard single-field slow-roll inflation), the primordial scalar and tensor perturbations are produced by quantum vacuum fluctuations of the metric.The resulting GW background has three distinctive properties: (i) a nearly scale-invariant spectrum (with a slight red tilt given by the inflationary consistency relation); (ii) an almost Gaussian probability density function (pdf); and (iii) no net circular polarization (i.e., non-chiral or parity conserving polarization).In this simple framework, r can be directly related to the energy scale of inflation [44].However, all previous properties characterizing the vacuum-produced GW background do not necessarily hold if matter sources (i.e.excited extra particle content) are active during inflation.Checking for the scale dependence of the tensor spectrum, the presence of non-Gaussianities in the tensor modes and the existence of parity-violating correlations, both at the CMB and interferometer scales, therefore becomes a necessary step before any claim on the origin of these tensor modes can be made [18,33,[45][46][47][48].Other models that can violate the above properties of vacuum-produced GW background are those introducing non-minimal coupling of the inflaton [49][50][51], additional scalar fields [52][53][54][55][56], and modified gravity models [57].Non-zero spatial curvature and kinetically-dominated initial conditions for inflation could also play a role in amplification or suppression of the GW background [58,59].
Model-building with additional matter sources typically poses a challenge: sources are always at least gravitationally coupled to the inflaton sector, which results also in excitation of non-Gaussian scalar modes [60][61][62][63], potentially clashing with the bounds from CMB data [46,64,65].In this context, two of the most successful categories of models are based on Abelian [52,60,64,[66][67][68][69][70][71][72][73][74][75] and non-Abelian [76][77][78][79][80][81][82][83][84][85][86][87][88][89] gauge fields.In this paper, we will focus on a model belonging to the latter category, which sources tensor modes from an SU(2) gauge field coupled to a pseudoscalar "axionlike" field through a Chern-Simons term [81].Since the axion and gauge fields are spectators (i.e., their energy density is subdominant compared to that of the inflaton), inflation is achieved through a standard inflaton sector [81].However, the Chern-Simons coupling breaks the conformal invariance of massless free gauge fields coupled to gravity, allowing for tachyonic amplification of gauge field perturbations, a process controlled by the speed of the axion rolling along its potential.Gauge field fluctuations, in turn, lead to a peak in the primordial tensor power spectrum, with a characteristic Gaussian bump shape [90].This model is unique because it can source tensor modes from matter fields at linear order in the perturbed Einstein equations without breaking the statistical isotropy of the Universe.Therefore, unlike similar models such as axion-U(1) inflation [46,64], which can source tensor modes only at second order in gauge-field perturbations, the axion-SU(2) model typically produces a negligible amount of sourced scalars and scalar non-Gaussianity compared to the vacuum-produced ones, allowing for sizeable amplitudes of a sourced GW background without spoiling agreement with current CMB bounds.On the other hand, the axion-SU(2) setup produces a strongly non-Gaussian tensor signal due to self-coupling of the gauge field, providing another distinctive (and potentially crucial) feature of this model compared to the inflationary paradigm [64,69,83,84].Generalizations to SU(N) gauge theories have also been considered and lead to the same phenomenology [91].
From the perspective of future B-mode experiments, the axion-SU(2) setup appears to be a particularly interesting candidate to probe: it can source tensor modes exceeding the vacuum contribution by a factor of ∼ 5 on the reionization bump scales in the CMB Bmodes, while both contributions can still be comparable on the recombination bump scales [92].As we argue in section 2.5 of [43], this feature of the model highlights the benefits of a full-sky survey with access to the reionization bump, such as LiteBIRD, when trying to distinguish between sourced and vacuum origins of the primordial GW background.Our motivation in this work is to investigate how well LiteBIRD can test the properties (i) and (iii) of the vacuum GW background, i.e., the approximate scale-invariance of the spectrum and the parity symmetry.Specifically, we discuss whether LiteBIRD observations can exclude standard single-field slow-roll inflation if tensors fluctuations arise from the production of matter in a way that breaks the approximate scale-invariance of the spectrum and produces non-zero parity-violating correlations.If observed data were found to be inconsistent with predictions of standard single-field slow-roll models, it would provide a strong motivation to test also the property (ii), i.e., the Gaussianity of the stochastic GW background.The latter topic will be the subject of detailed investigation in a follow-up LiteBIRD collaboration paper.
This work is part of a series of papers that present the science achievable by the LiteBIRD space mission, expanding on the overview published in Ref. [43].In particular, we expand on the discussion presented in section 2.5 of Ref. [43] with a more quantitative approach, using LiteBIRD map-based simulations (including the component separation and angular power spectra estimation steps), to build a robust power spectrum covariance matrix and a likelihood based on the Hamimeche & Lewis [93] approximation.We then perform a parameter inference on a representative selection of parameter choices for the model, using a frequentist Monte Carlo approach based on the Feldman-Cousins prescription [94,95] (see also [35,42,46,96,97] for recent applications in cosmology) to account for the presence of physical boundaries on the parameters.
We also investigate another unique signature of the axion-SU(2) model: T B and EB parity-violating correlations in the CMB due to the Chern-Simons coupling.We find that the T B and EB spectra produced by this model cannot be detected by LiteBIRD, and that almost all the constraining power on the axion-SU(2) model resides in the BB spectrum.We finally assess the power of LiteBIRD in discriminating between the two possible mechanisms for the production of gravitational waves, including information from T T , EE, BB, T E, T B and EB CMB spectra in the covariance matrix.
The rest of the paper is structured as follows.In section 2 we review the axion-SU(2) model of [81] and the latest bounds available in the literature.In section 3 we describe the simulations and the method used.In section 4 we emphasize the benefits of a full-sky Bmode mission to test the origin of the stochastic GW background.In section 5 we quantify the constraining power of LiteBIRD on the parameters of the axion-SU(2) model.In section 6 we investigate the additional discriminating power of the T B and EB angular power spectra.We present our conclusions in section 7.

Spectator axion-SU(2) gauge field inflation model
We will consider the spectator axion-SU(2) gauge field inflation model of Ref. [81], based on the "chromo-natural" inflation model [98] (see also Refs.[45,78] for reviews).The Lagrangian for this model, contains a standard inflaton sector L inf , which always dominates the energy density and is responsible for inflation, and an axion field χ with a cosine-type potential where µ and f are dimensioned parameters and λ is the dimensionless coupling constant for an SU(2) gauge field coupled to the axion via a Chern-Simons term.The SU(2) gauge field A µ = a A a µ σ a (where σ a are the Pauli matrices and a = {1, 2, 3}) has a field strength tensor given by is the dual field strength tensor and ϵ µνρσ is the totally antisymmetric symbol with ϵ 0123 = 1.Here, √ −g M is the determinant of the metric tensor.This model can source tensor modes at linear order in the perturbed Einstein equations without breaking the observed statistical isotropy of the Universe, because the SU(2) gauge field establishes a homogeneous and isotropic background solution, Āa i = a(t)Q(t)δ a i [76,77].This solution is approached even if the Universe was initially highly anisotropic (i.e. is an attractor solution) [99][100][101][102].The perturbation around this solution gives scalar, vector, and tensor modes [76,77].In particular, tensor fluctuations of the gauge field, which are subject to tachyonic amplification near horizon crossing, source gravitons at the linear level in the stress-energy tensor.Only one of the helicities of the gauge field is amplified due to the parity-violating Chern-Simons coupling, leading to a chiral GW background with left-or right-handed circular polarization [78,80,85,86].The primordial power spectrum of the sourced tensor modes, assuming that only left-handed gravitational waves are amplified, has a log-normal shape 1 [90] controlled by the wavenumber k p where the spectrum peaks, the effective tensor-to-scalar ratio at the peak scale r * , the width of the Gaussian-shaped bump σ and the power spectrum of scalar curvature perturbations P ζ .The three parameters {r * , k p , σ} can be related to the parameters in the model Lagrangian2 {g, λ, µ, f } (Eq.2.1) [90,103].The peak wavenumber k p corresponds to the time t * at which χ is at the inflection point of the potential, χ(t * ) = πf /2, reaching its maximum velocity.We also define the dimensionless time-dependent mass parameter of the gauge field fluctuations as m Q (t) ≡ gQ(t)/H (with H the Hubble expansion rate during inflation), which is m * ≡ m Q (t * ) = (g 2 µ 4 /3λH 4 ) 1/3 at the inflection point.
We can also define a dimensionless effective coupling . The effective tensor-to-scalar ratio at the peak scale can be defined as where ] and M Pl is the Planck mass, and can take any positive value, in principle 3 .On the other hand, since χ can remain at the top of its cosine-type potential hill only for a limited amount of time due to its quantum fluctuations, σ and k p must satisfy the relation [90]: where ∆N is the number of e-folds during which the axion is rolling down its potential and k CMB is roughly equal to the largest observable CMB scale.Whereas Eq. 2.2 and the following discussion hold only for axion potentials of the cosine-type and those with an inflection point, other power spectrum shapes are possible for different V (χ) [103].In all previous equations, P ζ receives negligible sourced contributions for m Q ≥ √ 2 [80,81] and can therefore be assumed to be equal to the vacuum scalar power spectrum , with the amplitude A s , spectral index n s and pivot scale k 0 = 0.05 Mpc −1 .The tensor power spectrum from the inflaton sector (indicated by "vac") receives instead a sourced contribution: where P vac t (k) = A t (k/k 0 ) nt is the tensor power spectrum from quantum vacuum fluctuations, with the amplitude A t and spectral index n t .We also define the tensor-to-scalar ratio of purely vacuum fluctuations as r vac = A t /A s .
The theoretical self-consistency of the axion-SU(2) setup has been studied in a number of papers, focusing mainly on the backreaction effect of particle production from the background axion and gauge fields on the background evolution, which could possibly ruin the phenomenological success of this model [81,92,[104][105][106][107][108][109][110].Recently, the study on spin-2 particle production has been significantly improved by solving equations of motion with backreaction for a wide range of model parameters [92].According to this study, the amplitude of the sourced tensor modes can exceed by more than a factor O(10) the vacuum contribution  at the CMB scales, and by several orders of magnitude at smaller scales.Including the effect of non-Gaussian scalar perturbations produced in second order by sourced tensor modes [111,112] reduces the allowed ratio of sourced-to-vacuum tensors to a factor O(1) at ℓ ≳ 80, where both the scalar power spectrum and scalar non-Gaussianity are tightly constrained by CMB temperature data [9,11,113], and to a factor ∼ 5 at low multipoles ℓ ≲ 10, where the temperature constraints are weaker.
In light of these theoretical bounds, we choose two sets of parameters of the axion-SU(2) model with the purpose of obtaining B-mode spectra with reionization bump scales (ℓ ≲ 10) as different as possible from the Starobinsky model of inflation [114] (r vac = 0.00461, n t = −r vac /8), while having similar recombination bumps (ℓ ∼ 80 − 100) (Fig. 1).These two sets of parameters are r * , σ, k p , r vac = [0.023,1.1, 3.44 × 10 −4 Mpc −1 , 0.00461], which gives the "high reionization bump" model (dot-dashed orange curve in Fig. 1) and r * , σ, k p , r vac = [0.002,1.9, 0.03 Mpc −1 , 0.002], the "low reionization bump" model (dashed purple).These two choices update those in figure 4 of Ref. [43] with viable models according to the study in Ref. [92].The CMB angular power spectrum for the "high reionization bump" remains approximately the same as in Ref. [43], although we use r vac = 0.00461 instead of the almost negligible value of r vac = 10 −4 previously assumed.This is because the ratio of sourced-to-  vacuum tensors can be at most ∼ 5 at low multipoles [92].In this new choice of parameters, the vacuum fluctuations provide a similar recombination bump as the Starobinsky model, while the sourced modes enhance only the reionization bump.However, now a value of r * = 0.023 (about half of the previous r * = 0.041) is sufficient to obtain the amplitude of the power spectrum similar to the one in Ref. [43] at scales ℓ ≲ 10.On the other hand, the new "low reionization bump" model, due to the fact that the ratio of sourced-to-vacuum tensors can be at most O(1) at ℓ ≳ 80 [92], is significantly different from the one in [43].In this case, both vacuum and sourced tensors contribute equally (i.e., r vac = r * ) to produce a recombination bump similar to that of the Starobinsky model.Figure 1 also shows for reference the cosmic variance-only (including primordial and lensing B-mode variance) and total LiteBIRD ±1 σ binned error bars (including foreground residuals) as gray and blue regions, respectively.
In this work, we also check that all the parameter choices are consistent with observational constraints on the axion-SU(2) model from the analysis of current CMB datasets [115].However, as also noted in [115], the shape of the tensor power spectrum is very weakly constrained by the Planck and BICEP/Keck data: the upper limits on the model parameters are prior dominated and strongly affected by degeneracies.

Method, simulations and likelihood
In this section, we will briefly describe the map-based simulations with LiteBIRD specifications, the component separation method, the angular power spectrum estimation procedure, and the likelihood used in this work.A sketch of the pipeline from input maps to inference on inflationary model parameters is shown in Fig. 2.
Simulations.The simulation setup used in this paper is the same as described in section 5.2 of [43] to which we refer the reader for full details.Here, we summarize the main points.The input CMB and Galactic foreground maps in each LiteBIRD frequency channel are generated at a HEALPix [116] resolution of N side = 512.The input maps of the CMB are generated as random Gaussian realizations from the lensed CMB angular power spectra corresponding to the values of the Planck best-fit cosmological parameters [9] and r vac = 0. Note that simulations for r vac ̸ = 0 and/or other models of inflation are built from the component separation outputs of the maps with r vac = 0 as explained in the "Component separation" section below.The input maps of the Galactic foreground are instead generated using the PySM package [117,118], and more specifically the default model identified as d1s1, which includes templates of thermal dust and synchrotron emission [119][120][121].The frequency dependence of the thermal dust signal is modeled as a modified black-body with temperature and spectral index varying across the sky, following [119].Similarly, polarized synchrotron emission follows a power law with a spatially varying spectral index [122].The sky templates from PySM are then convolved with top-hat bandpasses, coadded, convolved with Gaussian circular beams, and subsequently added to instrumental (white) noise realizations (see Table 3 in Ref. [43] for bandpass width, beam size, and noise levels considered).This procedure results in 1000 realizations of CMB, noise and foregrounds, which constitute the input maps for the component separation procedure described below.
Component separation.Component separation is performed using the FGBuster4 code [123], implementing a parametric fitting approach.The foreground model assumed in the code includes three spectral parameters: the spectral index and temperature of a modified blackbody spectrum for the dust; and the spectral index of the power-law spectrum of synchrotron radiation.The three spectral parameters are fitted independently over patches defined by HEALPix pixels at a given N side resolution [124].To find a suitable balance between the statistical uncertainty on spectral parameters, which increases with the number of pixels (i.e., of free parameters in the fit), and the need to capture the actual spatial variability of the foreground sky, we adopted different values of N side for each parameter and according to the Galactic latitude [43].
The FGBuster code returns a map including post-component separation noise and foreground residuals (hereafter collectively named "noisy foreground residuals") at N side = 64 for each of the 1000 input CMB and noise map realizations.Since these output maps are independent of the CMB input signal [125], we can sum them to Gaussian map realizations of a lensed CMB signal computed from an arbitrary cosmological model.In this way, we can get 1000 component-separated maps for each choice of our inflationary model parameters.This feature will be used in sections 4, 5 and 6 to obtain constraints on the axion-SU(2) model without rerunning the computationally expensive component-separation code for each set of model parameters.To avoid correlations, we used 500 simulations to estimate the C ℓ covariance matrix with a fiducial model and the other 500 simulations as simulated data (see also below).
Angular power spectrum estimation.We compute auto-spectra of maps using two different estimators of angular power spectra at different scales, as in Ref. [126]: a pure pseudo-C ℓ s method (implemented in the NaMaster5 code [127]) for intermediate and small scales (ℓ ≥ 35) and a quadratic maximum likelihood (QML) estimator [128,129] at large scales (ℓ < 35).The pseudo-C ℓ s estimator is nearly optimal for smaller scales, but suboptimal at large scales [130,131].The QML method, despite being more computationally expensive compared to the pseudo-C ℓ s, produces nearly optimal variance estimates of the power spectrum at large scales.We use the QML implementation publicly available in the xQML6 code [129].The two estimators are computed separately on the maps and then joined to form a single auto power spectrum.As specified above, the noisy foreground residuals returned by FGBuster have a resolution of N side = 64; however, since the xQML implementation is memory intensive and therefore not feasible at this resolution, we downgrade the maps to N side = 16, after convolution with an apodizing kernel to avoid aliasing [9].
Sky masks.We use two different masks for the two power spectrum estimators described above: a binary mask obtained from the LiteBIRD polarization mask with sky fraction f sky = 49.5 % presented in section 5.2.4 of [43] for the QML estimator and the Planck Galactic plane mask appropriately apodized for the pseudo-C ℓ s (Fig. 3).We optimize the mask properties at both the spectrum and the parameter likelihood levels.First, we minimize both the variance of the BB spectrum and the residuals of the simulations with respect to the input spectrum divided by the error on the mean at each multipole, while varying the threshold of the binary mask in the QML part and the apodization scale and sky fraction in the pseudo-C ℓ s one.
We then verify that the mask configurations explored in the first step actually minimize the uncertainty and bias on the tensor-to-scalar ratio r vac for a fiducial Starobinsky model (see also discussion below).We find the optimal configuration by downgrading the LiteBIRD mask to N side = 16 and thresholding it at a value of 0.75, resulting in f sky = 47 % for the QML estimator, and by apodizing the Planck f sky = 60 % mask at a 10 • scale with the C 2 method [130], resulting in an effective f sky = 51 % for pseudo-C ℓ s.
Likelihood.We perform statistical inference on the inflationary model parameters using a likelihood based on the Hamimeche & Lewis (hereafter H&L) [93] approximation.Specifically, we add an offset term, as prescribed in Ref. [132], to deal with negative values of the estimated power spectra on large scales.We also use the Sellentin & Heavens (hereafter S&H) [133] correction to the H&L likelihood, accounting for the increased uncertainty in the parameters due to the finite number of simulations used in the estimation of the covariance matrix.The block BBBB of the C ℓ covariance matrix estimated from simulations is shown in Fig. 4  we checked that the Monte Carlo noise due to the limited number of simulations used to estimate the covariance does not affect parameter inference: adopting the conditioning strategy suggested in Ref. [34] leaves parameter estimates unchanged.All fits in this paper are performed by minimizing the negative logarithm of likelihood −2 log L by the multidimensional minimizer iMinuit7 .

Disentangling vacuum and sourced origins: benefits of a full-sky mission
In this section, we highlight the benefits of a full-sky B-mode mission with access to the reionization bump scales in helping to understand the origin of primordial GW background.In the event of a detection of primordial B modes, this feature could allow one to exclude quantum vacuum fluctuations of spacetime within the standard single-field slow-roll paradigm as the only source of the observed GW background, pointing instead towards production by matter sources, that is, the additional excited particle content active during inflation.From this perspective, the axion-SU(2) model can provide a useful benchmark: it can source tensor modes that exceed by a factor O(5) the vacuum contribution at the largest scales, producing a pronounced feature in the reionization bump, while producing at the same time a recombination bump essentially indistinguishable from the standard single-field slow-roll prediction [43] (Fig. 1).
Here we expand on the discussion in Ref. [43] in a more quantitative way.Suppose that the true sky was generated from the axion-SU(2) inflation.Then what if the observational data were fitted with a "wrong" model, in this case the standard single-field slow-roll one?To this end, we generate LiteBIRD simulations for the "high reionization bump" parameter choice of the axion-SU(2) model (corresponding to the orange dot-dashed curve in Fig. 1) using the procedure described in section 3, and fit the BB spectrum of each simulation for the tensorto-scalar ratio parameter r vac (as defined in section 2), assuming for simplicity that all other 0.000 cosmological parameters are fixed by T T , T E and EE.This is a good approximation, as the axion-SU(2) model provides negligible contributions to scalar perturbations, allowing us to break the degeneracy between the inflationary parameters and the ΛCDM parameters.The fit is performed three times for each simulation, taking three different ranges of multipoles of the BB spectrum: only the reionization bump range (ℓ = 2 − 30), only the recombination bump range (ℓ = 31 − 150) and finally the full range of multipoles (ℓ = 2 − 150).The resulting histograms are shown in Fig. 5, while we report in Table 1 the mean value and standard deviation associated with the fit for each case.As is evident, the fit recovers very different amplitudes for r vac in the three cases, potentially leading to incorrect interpretations within the standard single-field slow-roll scenario.As was expected from our choice of the axion-SU(2) parameters, the data point towards the Starobinsky model if the survey does not have access to the reionization bump, even though the true sky features a significant contribution from matter sources.We note that if we also fit, in addition to r vac , the underlying "high reionization bump" SU(2) simulations for the tensor spectral index n t , we expect n t to be partially degenerate with r * , σ and k p .In this case, the recovered n t would be very red, violating the standard single-field slow-roll consistency relation and pointing towards the presence of more than one field during inflation.Checking the Gaussianity of the signal would then be essential in discriminating between axion-SU(2) and other models (see section 1 for more details).

Constraining axion-SU(2) parameters with LiteBIRD
In this section, we explore the possibility of constraining the parameters of the axion-SU(2) model using the LiteBIRD satellite.We infer parameters using a frequentist approach based on the Monte Carlo simulations constructed in section 3, and account for the presence of physical boundaries on the parameters following the Feldman-Cousins prescription [94].We start in section 5.1 by describing the Feldman-Cousins approach and applying it to our case study, and in section 5.2 we explore the robustness of this construction.

Feldman-Cousins approach
Frequentist confidence intervals for the parameters of interest are built using the classical Neyman's construction [134].However, if the parameter has a physical boundary and its estimate is close to this boundary, Neyman's construction must be corrected as indicated by Feldman & Cousins [94] (hereafter FC).This approach is a standard staple in particle physics data analysis [135].Recently there has been renewed interest in using it for cosmology, being applied, for example, to Planck data [97], SPIDER B-mode data [42], early dark energy models [96], the tensor-to-scalar ratio from the Planck and BICEP/Keck data [35], and the search for signatures of axion-U(1) inflation in the current data [46].This technique is not plagued by prior volume effects that can arise during the marginalization procedure in the Monte Carlo Markov Chain (MCMC) approach, when performing inference on degenerate/unconstrained parameters [35,96,136].The Feldman-Cousins approach is therefore an ideal choice for our case study, since we expect the axion-SU(2) model parameters to be degenerate and hard to constrain simultaneously.
The FC recipe for our parameter of interest, r * , is the following: 1.For each of the input values in the physically-allowed range of interest, r in * , we generate a sufficient number of simulations.
2. We fit a model to each simulation, obtaining values of r mle * that maximize the likelihood given the input simulation.These can be used to build a histogram of r mle  = 0.015 as the observed value, obtained as the best fit to the simulation generated from the "high reionization bump" parameters (hereafter "ground truth" values) with r in * = 0.023.In Fig. 6, we show the FC construction (solid black lines) in which we fit the simulations only for r * (hereafter 1 parameter fit), fixing σ and k p to their ground truth in the template (Eq.2.2).In this case, we observe that LiteBIRD will be able to obtain a two-sided 95 % C.L. confidence interval of r * = 0.015 +0.04 −0.008 , if the observed sky has been generated from the "high reionization bump" model.In the next section, we will discuss the robustness of these constraints against a mismatch between the assumed values of σ and k p in the fit and the ground truth, as well as the degradation of the constraint in r * when we also treat σ and k p as free parameters in the fit.

Robustness of the constraints
In this section, we will explore how constraints on r * from the FC construction can change if we adopt "wrong" assumptions about σ and k p in the 1-parameter fit.We also assess the  degradation of the constraints on r * when σ and k p are additional free parameters in the fit, exploiting the full power of the FC construction: the parameter inference performed with this technique is free from prior volume effects due to degeneracies/unconstrained parameters in the Bayesian MCMC [35,96,136].
We start by examining the FC confidence belt built from the same simulations as in section 5.1 and fitting each one for r * ; however, in this case, we consider three different values of the peak scale k p = [10 −4 , 5 × 10 −4 , 5 × 10 −3 ] Mpc −1 as shown in Fig. 7, instead of fixing it to its ground truth, k p = 3.4 × 10 −4 Mpc −1 .We still fix σ to its ground truth; similar conclusions can be drawn by varying σ instead of k p .The behavior of the confidence belts in Fig. 7 can be easily interpreted.The larger r in * is, the larger the power in the underlying simulations at the reionization bump becomes.However, if we fix k p = 5 × 10 −3 Mpc −1 in the fit (green dashed), the template model can add power only at recombination bump scales (k ∼ 5 × 10 −3 Mpc −1 ): this results in a preference for small values r mle * around zero, since the power in simulations at recombination bump scales is already saturated by r vac .Similarly, for the k p = 5 × 10 −4 Mpc −1 case (dot-dashed blue), the confidence belt is slightly tighter compared to the initial one (solid black) as less power is allowed at smaller scales compared to reionization bump scales in the simulations, while the opposite is true for the k p = 10 −4 Mpc −1 case (dotted orange).The best choice of k p can be found by comparing the value of the likelihood for each model considered.Now, we check how much the constraints presented in section 5.1 are degraded when also σ and k p values re allowed to vary in the fit.In Fig. 8, we show the confidence belt for the 3-parameter fit (r * , σ, k p ).In this case, LiteBIRD will be able to put a 95 % C.L. upper  Table 2. 95 % C.L. confidence intervals or uppper limits on r * obtained from the FC construction in the 1-, 2-and 3-parameter fits.We assume as the underlying model the "high reionization bump" model (section 2).See section 5.2 for more details.
limit on r * ≤ 0.16, if the observed sky has been generated from the "high reionization bump" model.This significant degrading of the constraining power is mainly due to the degeneracy among the three parameters of the model given in Eq. 2.2.We also compare in Table 2 the 95 % C.L. confidence interval or upper limit on r * obtained from the FC construction in the 1-, 2-and 3-parameter fits.We report the observed value r obs * and the corresponding error bars or upper limits.For the 1-parameter fit, we fix σ and k p to their ground truth values, while in the 2-parameter case we additionally fit for σ while k p is fixed to the ground truth value.

Parity-violating T B and EB correlations
In this section, we investigate another unique signature of the axion-SU(2) inflation model when trying to disentangle it from standard single-field slow-roll models: T B and EB parity-violating correlations in the CMB produced by chiral gravitational waves [49,69,[138][139][140].We finally assess the full power of LiteBIRD in discriminating between two mechanisms for enhanced gravitational wave production, using all CMB spectra (T T , EE, BB, T E, T B and EB) in the covariance matrix.
Concerning T B and EB spectra, we improve and update the LiteBIRD forecast presented in [90] in a number of aspects.First, we relax the assumption used in [90] of the same CMB upper bound on r vac applying to all scales.Instead, we take into account that at reionization bump scales this bound becomes considerably relaxed compared to the one at at k = 0.05 Mpc −1 [32,126,141].This allows, in principle, for a larger signal from the axion-SU(2) model at these scales.Second, we update the forecast with new bounds on the SU(2) parameter space from backreaction [92] and the current upper limits from the Planck and BICEP data [31,32,35,115].Third, we show results for the full covariance matrix from realistic LiteBIRD simulations, instead of the simplified Fisher approach of [90].Fourth, since we found that the signal-to-noise ratio defined in [90] is strongly dependent on the fiducial model adopted in the covariance matrix, we use as a figure of merit ∆χ 2 between the standard single-field slow-roll fiducial (specifically, the Starobinsky model) and axion-SU(2) models, that is where M is the C ℓ covariance matrix, C fid ℓ is a vector containing the fiducial angular power spectra and similarly C SU(2) ℓ contains the theoretical spectra for the axion-SU(2) model.We quantify in Figs. 9 and 10 the discriminating power of LiteBIRD between the axion-SU( 2) and standard single-field slow-roll Starobinsky models (r vac = 0.00461) via ∆χ 2 values as a function of r * and σ given a fixed scale k p .For each point in this parameter space, we compute the quantity ∆χ 2 in Eq. 6.1: larger positive values of ∆χ 2 indicate an increased ability of LiteBIRD to discriminate axion-SU(2) from the fiducial model.The associated pvalue, i.e. the probability that the χ 2 statistic should exceed a particular value ∆χ 2 by chance, given that the null hypothesis (i.e. the fiducial model) is correct, can then be computed as [95] p(∆χ where χ 2 n (x) is the distribution for n degrees of freedom (hereafter d.o.f.).In our analysis, the number of d.o.f. will be equal to the total number of multipole bins in all spectra considered 9 .The p-value can be converted into an equivalent significance Z, defined such that a Gaussian distributed random variable, fluctuating Z standard deviations above its mean, has an uppertail probability equal to p, that is, Z = Φ −1 (1 − p), where Φ is the cumulative distribution of the standard Gaussian and Φ −1 its inverse quantile [95].For instance, Z = 5 (i.e. a 5 σ detection) corresponds to p = 2.87 × 10 −7 .We use this relation to draw reference contours of significance of 1 σ, 3 σ, 5 σ and 8 σ in Figs. 9 and 10.Note that this σ, although named the same, is not the parameter σ of the axion-SU(2) model, which is instead plotted on the x-axis of these figures. 9For this analysis we bin angular power spectra with equal width ∆ℓ = 10, resulting in 14 bins for each spectrum.We additionally checked that the distribution of ∆χ        In Fig. 9, we compare the contribution to ∆χ 2 from different blocks of the covariance matrix (the full covariance, BBBB, T BT B and EBEB blocks) in each panel.In Fig. 10, instead, we compare ∆χ 2 for different choices of k p that can be probed by CMB observations, using the full covariance matrix.The white area in both figures represents the portion of parameter space excluded by current upper limits on the tensor-to-scalar ratio by Planck and BICEP/Keck data at each scale [31,32,35,115,126].Note that we assume r * = 5r vac for all k p choices, except for k p = 5 × 10 −3 Mpc −1 for which we assume r * = r vac , following Ref.[92] (section 2).In Figs.11, 13 and 12, we also show the BB, |T B| and |EB| power spectra, respectively, used in Fig. 10.
We find that the contribution to the total ∆χ 2 is entirely dominated by BB: LiteBIRD will be able to exclude with high significance the Starobinsky model if the tensor fluctuations arise from axion-SU(2) inflation and if the typical features of this latter model occur at scales    observable by the CMB.Specifically, the significance is higher than 8 σ for parameter close to the current upper limits on axion-SU(2) parameters, while T B and EB remain under 1 σ significance in all cases considered, making it impossible for LiteBIRD to detect parityviolation generated by the axion-SU(2) model.Furthermore, as noted in Ref. [90], the T B contribution dominates over the EB one, making ∆χ 2 larger by a factor ∼ O(10 2 ).
Comparing the upper left and right panels of Fig. 9, we find that, for a given parameter choice, using only the BBBB block results in higher significance compared to the full covariance matrix.This is because differences between the fiducial and axion-SU(2) model exist only in the BB spectrum and marginally in T B and EB.Therefore, using the full covariance just increases the number of degrees of freedom, diluting the information content.Although the inference of parameters by extracting a single block from the full covariance matrix is not strictly correct, it is still useful in this context to understand in an approximate fashion the contribution of different angular power spectra to the final LiteBIRD sensitivity.
Our analysis can be improved.Since we are using the same sky mask, which was optimized for B modes (see section 3), for both temperature and polarization maps to compute the full covariance, the importance of the T B spectrum in the analysis could further increase   if we optimized masks to temperature field with a larger sky fraction.
In summary, we conclude that LiteBIRD will be able to exclude with high significance tensor fluctuations produced by vacuum fluctuations (specifically, in the Starobinsky model) if the GW background has been produced by the axion-SU(2) mechanism and the feature is sourced at the CMB scales.We also find that the BB spectrum will dominate the discriminating power of LiteBIRD, while T B and EB correlations would be of secondary importance.In the axion-SU(2) model, in fact, the production of detectable parity-violating correlations would imply an overproduction of B modes, which would conflict with current observational limits on tensor modes.

Conclusions
Enhanced primordial gravitational waves from gauge fields during inflation represent a new paradigm of the primordial universe that has been extensively studied in the literature [45,78].In this paper, we used realistic simulations to show that, thanks to the access to the reionization bump scales provided by a full-sky CMB mission, LiteBIRD can provide significant help in our effort to distinguish between inflation models, more specifically in excluding the production of the stochastic GW background within the standard vacuum fluctuations paradigm in favor of production by matter sources and vice versa.We also presented expected constraints on the model parameters of an SU(2) model with a characteristic "bump-like" feature in the reionization bump.In this case, LiteBIRD will be able to obtain a two-sided 95 % C.L. confidence interval on the bump feature amplitude r * .
The SU(2) model, in contrast to the standard inflationary scenario, also predicts parityviolating correlations in the CMB.We include the contribution of T B and EB in the full covariance matrix and assess the ability of LiteBIRD to disentangle standard single-field slowroll (specifically the Starobinsky model, a reference target for LiteBIRD) from the axion-SU(2) model, finding that this experiment will be able to distinguish them with high significance for a wide range of model parameters.We also find that the discriminating power of LiteBIRD will be determined mainly by BB, with T B and EB giving negligible contributions in almost   all the allowed parameter space.Detecting the parity-violating signal in T B and EB from the axion-SU(2) model remains out of sight for LiteBIRD.
We enforced bounds on sourced-to-vacuum tensor perturbations from the backreaction effect of spin-2 particle production on the background fields and from the non-Gaussian contribution to the scalar curvature perturbations produced at second order by sourced tensor modes, in order for the analytical template (Eq.2.2) to remain valid.However, we stress that lattice simulations [142,143] can be used to relax these bounds and study the system also in the strong backreaction regime.
In this paper, we used the power spectrum as our main observable, taking advantage of the strong scale dependence of the tensor spectrum and the presence of parity-violating correlations predicted in axion-SU(2) models.However, additional information is available in the bispectrum and higher order correlation functions, potentially allowing to distinguish vacuum produced from sourced tensors when the power spectrum does not present significant features at CMB scales [64,69,83,84].We plan to explore this possibility in a dedicated future LiteBIRD paper.
Finally, we mention the connection between this model and recent cosmic birefringence measurements.It has been shown in Ref. [144] that the large EB signal recently observed [145][146][147][148] cannot be explained by a parity-violating GW background from the axion-SU(2) model, due to the simultaneous amplification of the BB power spectrum that violates current upper bounds on the tensor-to-scalar ratio.If any signal is present in EB from primordial parity-violation due to matter sources during inflation, it is expected to be subdominant with respect to the observed cosmic birefringence signal.

Figure 2 .
Figure 2. Pipeline for inflationary model parameter constraints used in the paper.See sections 3, 5 and 6 for details.

Figure 3 .
Figure 3. Left panel: Binary mask (f sky = 47 %) used for the large-scale power spectrum estimation (QML) at N side = 16 resolution.Right panel: Planck Galactic-plane apodized mask (f sky = 51 %) at N side = 64 used for intermediate/small-scale pseudo-C ℓ s estimation (see section 3 for details).

Figure 6 .
Figure 6.FC construction for the 1-parameter fit (solid black lines), with σ and k p fixed to their ground truth.LiteBIRD will be able to obtain a two-sided 95 % C.L. confidence interval r * = 0.015 +0.04 −0.008 , if the observed sky has been generated from the "high reionization bump" model.The observed value r obs * = 0.015 is indicated as a vertical solid red line, and its intersection with the confidence belt with two red dots.The color bar shows the distribution of r mle * obtained by fitting the simulations as a function of the input r in * .See section 5.1 for details.

Figure 7 .
Figure 7. FC construction for the 1-parameter fit with k p fixed to values different from its ground truth, i.e. k p = [10 −4 , 5 × 10 −4 , 5 × 10 −3 ] Mpc −1 (dotted orange, dot-dashed blue and dashed green, respecitvely), compared to the ground truth confidence belt (solid black, same as in Fig. 6).We assume the ground truth for σ.The color bar shows the distribution of r mle * obtained by fitting the simulations (assuming the ground truth for k p and σ) as a function of the input r in * .See section 5.2 for details.

Figure 8 .
Figure 8. FC construction for the 3-parameter fit (solid black lines) with free parameters r * , σ and k p .LiteBIRD will be able to obtain a 95 % C.L. upper limit r * ≤ 0.16, if the observed sky has been generated from the "high reionization bump" model.The observed value r obs * = 0.05 is indicated as a vertical solid red line.The color bar shows the distribution of r mle * obtained by fitting the simulations as a function of the input r in * .See section 5.2 for details.
model simulations (i.e., the null hypothesis), where C sim, i ℓ is the spectrum estimated from the ith simulation, follows a χ 2 distribution with the expected number of d.o.f.

2 Figure 9 .
Figure 9.Comparison of ∆χ 2 (Eq.6.1) between the axion-SU(2) and standard single-field slowroll Starobinsky models from the full covariance matrix (top left panel), BBBB (top right), T BT B (bottom left) and EBEB (bottom right) covariance blocks for LiteBIRD, as a function of r * and σ given k p = 3.44 × Mpc −1 .The contours show 1 σ, 3 σ, 5 σ and 8 σ significance (note that this σ, although named the same, is not the parameter of the axion-SU(2) model, which is instead plotted on the x-axis).The white area is excluded by current upper limits on the tensor-to-scalar ratio.The red ⋆ symbol in the upper right panel indicates the "high reionization bump" model.Note that the color scale changes in every panel.

2 Figure 10 .
Figure 10.Same as Fig.9but each panel has a different k p and we use the full covariance matrix for all panels.Here we assume r * = 5r vac in the top two and bottom left panels, while the bottom right panel assumes r * = r vac , following Ref.[92].Note that the color scale changes in every panel and y-axis range changes in the bottom right panel.

Figure 11 .
Figure 11.BB power spectra for each of the parameter sets used in each panel of Fig.10.They correspond to spectra computed on a grid of 100 linearly spaced values of r * in the range 0.001 − 0.36 and 50 values of σ in the range 1 − 10 at each fixed k p value, excluding spectra not compatible with theoretical consistency arguments and observational bounds (see section 2 for details).Darker color lines correspond to larger r * and σ in this range.We also show for reference the BB spectrum for standard single-field slow-roll inflation obtained for r vac = 0.037 (solid red), saturating the current upper at scales 0.05 Mpc −1 , the "high reionization bump" axion-SU(2) model (dot-dashed orange, see section2) and the Starobinsky model (dotted black).

Figure 12 .
Figure 12.Same as Fig.11but for |T B| power spectra.We also show for reference the |T B| spectrum for the "high reionization bump" axion-SU(2) model (dot-dashed orange, see section 2).

Table 1 .
Mean Figure 5. Histograms of the inferred vacuum tensor-to-scalar ratio r vac values obtained from fitting an observed sky generated from the "high reionization bump" axion-SU(2) model (with parameters r * , σ, k p , r vac = [0.023,1.1, 3.44 × 10 −4 Mpc −1 , 0.00461]) in three different ranges of multipoles: only the reionization bump range (ℓ = 2 − 30, in red), only the recombination bump range (ℓ = 31 − 150, in yellow) and the full range (ℓ = 2 − 150, light blue).We also show for reference the r vac value for the Starobinsky model (r vac = 0.00461).values (indicated by rvac ) and standard deviation of the distributions of the r vac parameter in Fig. 5, obtained by fitting a true sky generated from the "high reionization bump" SU(2) model in three different ranges of multipoles.