Likelihood analysis of small field polynomial models of inflation yielding a high Tensor-to-Scalar ratio.

Inflationary potentials, with Planckian field excursions, described by a 6th degree polynomial are studied. We solve the Mukhanov-Sasaki equations exactly and employ a probabilistic approach as well as multinomial fitting to analyse the results. We identify the most likely models which yield a tensor-to-scalar ratio r = 0.01 in addition to currently allowed Cosmic Microwave Background (CMB) spectrum and observables. Additionally, we find a significant inter-dependence of CMB observables in these models. This might be an important effect for future analyses, since the different moments of the primordial power spectrum are taken to be independent in the usual Markov chain Monte Carlo methods.

1 Introduction

In a previous article [1], a class of small field inflationary models which are able to reproduce the currently measured Cosmologic Microwave Background (CMB) observables, while also generating an appreciable primordial Gravitational Wave (GW) signal were studied. The existence of such small field models provides a viable alternative to the large field models that generate a high Tensor-to-Scalar ratio. Our exact analysis was shown to give accurate results [1]. Models which yield Tensor to Scalar ratio (r), of less than r ≲ 0.003 were previously studied in [1]. The initial study demonstrated a significant difference between analytical Stewart-Lyth [2, 3] estimates and the exact results. This result should be confronted with analyses such as in [4] where the Stewart-Lyth expression is relied upon, and [5] in which the authors use a Green’s function approach and perturbation theory, but assume the log of the input is well behaved. Our method extends and improves the method of the model building technique employed in [6, 7]. Previous analytical work [8-10] has shown that a fourth-order polynomial potential is sufficient to generate a high tensor-to-scalar ratio, even up to r ≳ 0.1. However, it was hard to realize this numerically. It was discovered in [1], that a fifth-order polynomial potential was required for generating 0.001 ≲ r ≲ 0.003. Furthermore a sixth-order polynomial seems to be required for a tensor-to-scalar ratio greater than r ≳ 0.003. A simple explanation is offered by observing (see Fig 1) that increasing r by factor ∼ 10, causes the e-folds per field excursion generated at the CMB window to decrease by a factor of ∼ 3. This means widening the CMB window, and losing the decoupling between the CMB window and the e-fold generating peak. Adding the 6th coefficient pushes the peak from ϕ ∈ [0.4, 0.5], to higher values of ϕ and decouples these regions.

Fig 1

A graph depicting as a function of the inflaton ϕ for two fifth-order polynomial models, and a sixth-order polynomial model.

For a fifth-order polynomial model with r = 0.001 (Blue line) the CMB window width is ∼8 e-folds, while the field changes by about Δϕ ∼ 0.1. Most of the e-folds are generated when ϕ reaches ∼0.4. When r is increased the CMB window widens, and approaches the e-fold generating peak (Red dots). While marginally affecting the CMB window width, the introduction of an additional coefficient, a6, allows to shift the peak to higher values of ϕ, thereby decoupling the CMB window and the e-fold generating peak (Green dash).

1.1 Conventions

In this article we follow the conventions of [11]. The Primordial Power Spectrum (PPS) is given by: Our conventions are: The PPS is expanded about the pivot scale k0. The pivot scale is usually set a-posteriori, at the scale in which the parameters {n − 1, α, β} are minimally dependent [12-14], in Planck + BICEP2 data analyses it is usually set at k0 = 0.05 hMpc−1.

2 Inflationary models

The small field models previously studied in [1] yielded results that are consistent with observable data up to values of r ≃ 0.003. While these values agree with the current limits on r set by Planck [15, 16], we are interested in studying models with higher r. The study of small field models is motivated by their appearance in many fundamental physics frameworks, effective field theory, supergravity [17] and string theory [18] in successive order of complexity. For models with r ≳ 0.003, significant running of running is found. This means that while three free parameters (corresponding to n, α, N) were previously needed, we now need an additional free parameter. Therefore we turn to a model of a degree six polynomial potential. Obviously considering higher degree models complicates the analysis by adding other tunable parameters. The potential is given by the following polynomial: It has been shown [6] that the potential can be written as: However, for simplicity, we express the potential as follows: with the subscript 0 denoting the value at the CMB point. By setting ϕ0 = 0; ϕ = 1 we limit ourselves to small field models in which Δϕ = 1 in Planck units, with little effect on CMB observables. According to the Lyth bound [19, 20], given a Tensor-to-Scalar ratio of r ≃ 0.01, the lower bound on the field excursion is approximately given by Δϕ4 ≳ 0.03 m. Here Δϕ4 is the field excursion while the first ∼4 efolds are generated. Our models satisfy this strict bound, as the first 4 efolds or so typically result in Δϕ4 ∼ 0.15 which is well above 0.03. The Lyth bound was further extrapolated [21] to cover the entire inflationary period. Applying this approach to models with r ∼ 0.01 yields Δϕ ≃ 2 m. However, in [7], it was shown that in models such as the ones we study, the value of Δϕ can be smaller because ϵ is non-monotonic. In this case, Δϕ = 1 m from the CMB point to the end of inflation is consistent with the Lyth bound.

When the coefficients {r0, a2, a3, a4} are fixed, the remaining coefficients are related by: The procedure of finding f1 and f2 was explained in detail for the degree 5 polynomial models in [1], and here we follow a similar procedure for the degree 6 models. So, ultimately, the model is parametrized by 5 parameters: the two physical parameters r0 and N and the three other parameters (a2, a3, a4) that are used to parametrize the n, α, β parameter space. It should be pointed out that N is not an observable, rather N ∼ 50 − 60 is a ‘soft’ constraint. Strictly speaking, N depends on the reheating temperature and only its maximum value can be determined. However for simplicity we treat N as an observable, in order to facilitate the study of a large sample of models.

3 Coefficient extraction methods

In this section we explain the two methods for calculating the most likely coefficients {a2, a3, a4}, given a large number of simulated models and the likelihood data for the CMB observables. This data is available through MCMC analysis of CMB data, such as the Planck data.

3.1 Likelihood assignment method—Gaussian extraction

To each potential, after calculating the observables n, α, β, we assign a likelihood. For each observable we calculate the likelihood according to the CosmoMC [22] likelihood analysis of the data sets used. We then assign the product of the likelihoods to the potential. A concrete example is the following: suppose we extract the trio (n, α, β) = (0.96, 0.011, 0.024), we look up the likelihoods: . We now multiply them, and so the likelihood attached to that specific model which yielded these observables is given by . We proceed to extract likelihoods for the different coefficients by a process of marginalization. The expectation is that this method will yield a (roughly) Gaussian distribution for each of the values of a2, a3, a4. The advantage of this method is in yielding not only the most likely value, but also the width of the Gaussian. This width can then be used as an indication for the level of tuning that is needed in these models.

3.1.1 Possible pitfalls

This method of likelihood assignment is vulnerable in two ways:

To be valid, this method requires a uniform cover of the relevant parameter space, by the potential parameters. If the cover significantly deviates from uniform, the results might be skewed by overweighting areas of negligible weight, or underweighting areas of significant weight. Fig 2 shows a mostly uniform cover.

Fig 2

Small field inflationary potentials which yield r = 0.01, as well as PPS observables within 68% and 99% confidence levels.

Every (n, α) pair is accessible using these models. The likelihood curves are results of a CosmoMC run with the latest BICEP2+Planck data.

Since only if the paired covariance is small, we must make sure that this is the case. In our underlying CosmoMC [22] analysis this is indeed the case. The covariance terms are, in general, one to two orders of magnitude smaller than the likelihoods at the tails of the Gaussian.

We also run the risk of false results if the fit we apply to the data points produced by the numerical analysis yields a large fitting error. However the fitting error of the polynomial function to the log PPS − log k data, is usually of the order of 10−6. This fitting is done over 30 data points generated by the MS equation numerical evaluation, for each potential. The error is calculated as , thus roughly speaking the error per data point is of the order of 10−6∼7. We conclude that the log PPS function is well fitted.

3.2 Multinomial fit

Another method for calculating the most likely coefficients is by fitting the simulated data with a multinomial function of the CMB observables. We aim to find a set of functions F such that, for example, a2 = F2(n, α, β). We assume that this function is smooth and thus can be expanded in the vicinity of the most likely CMB observables. Hence, we can find a set of multinomials (F2, F3, F4), such that: We have found that a quadratic multinomial is sufficiently accurate and that using a higher degree multinomial does not improve the accuracy significantly. Thus we may represent these by a symmetric bilinear form plus a linear term, as follows: where O = (n, α, β), B is the bilinear matrix, and the linear coefficient vector is A.

3.3 Pivot scale

So far, we discussed matching potentials and their resulting PPS around the CMB point. However, in order to correctly compare results of the PPS to observables, one has to take into account the pivot scale at which the CMB observables are defined. Since, in this case, the pivot scale is given by k0 = 0.05 hMpc−1, and the CMB point is at k ∼ 10−4 hMpc−1, the observables in the CMB point and k0 should be related in a simple way only if the spectrum varies slowly with k. This is not true for the case at hand. Two potentials can yield very different power spectra near the CMB point, and nevertheless yield the same observables at the pivot scale. These degeneracies, stem from our limited knowledge of the power spectra on small scales, and at the CMB point. For concreteness take two PPS functions, one that is well approximated by a cubic fit near the pivot scale, and the other that is well approximated only when we consider a quartic fit. Suppose, additionally, that these two PPS functions have the exact same first three coefficients, it follows that they yield the exact same observables {n, α, β}. However if we go to sufficiently small scales, or sufficiently high k values, these functions will diverge. This is also true at the large scale end, where the CMB point is set. Hence the degeneracy.

A possible solution to this problem, is classifying the resulting power spectra by the level of minimal good fit. We define a good fit as one in which the cumulative relative error , is less than 10−7. Given a single power spectrum, we fit our result with a polynomial fit, increasing in order until the accumulated relative error is sufficiently small. The minimal degree polynomial fit that approximates the log(PPS) function to the aforementioned accuracy is called the minimal good fit.

We then study separately power spectra that are well fitted by cubic polynomials, quartic polynomials etc. In this way we make sure that we compare non-degenerate cases.

4 Monte Carlo analysis of Cosmic Microwave Background with running of running

In [11] it was shown that the inclusion of additional parameters, i.e., the running of the spectral index (α), and the running of the running (β) resolves much of the tension between different data sets. In this section we briefly discuss the effect of considering non-vanishing α and β on the most likely shape of the PPS. First we find n when it is the only free parameter. We then use n, and α as the free parameters, and finally we conduct an analysis with n, α, and β as the free parameters. The shape of the power spectrum changes significantly when running of running is considered.

The data sets that were used are the latest BICEP2+Planck baseline [15], along with the low l’s [23], low TEB and lensing likelihoods. The results of these analyses are given in Table 1, as well as in Fig 3. As expected the resulting power spectra converge at the pivot scale k0 = 0.05 hMpc−1. However for lower k’s, the resulting spectra diverge considerably, consistent with cosmic variance. Notably, the spectra also diverge at higher k’s. This indicates the inability of current observational data to constrain the models in this range of k’s. This inability is also demonstrated in Fig 4 where, for l > 1500, the most restrictive data cannot rule out models with significant running, or running of running. Fig 4 also shows that the three models are virtually indistinguishable in terms of the observed C’s.

Table 1

Results from 3 analyses of the latest BICEP2+Planck dataset, each adding a free parameter in the power spectrum.

The results shown are best fits, within the 68% confidence level for each analysis.

Parameter (68%)	free ns	free αs	free βs
log(1010 As)	3.1047 ± 0.0057	3.1073 ± 0.006	3.1061 ± 0.0065
ns	0.9751 ± 0.0045	0.973 ± 0.0057	0.9687-0.006+0.0051
αs	N/A	−0.009 ± 0.0067	0.008±0.013
βs	N/A	N/A	0.020 ± 0.013

Fig 3

Power spectra as recoverd using CosmoMC analysis with latest BICEP2+Planck data.

Allowed area (68% CL) for a fixed n analysis is shown (blue). Similarly a fixed α analysis (red), and a fixed β (green) are shown. The other colores are intersection areas. The pivot scale in this graph is at , where k0 = 0.05 hMpc−1. The apparent divergence in high k’s is due to the inability of Planck to constrain these k’s. With more data, it will be possible to differentiate between the three possibilities.

Fig 4

Power spectra in the C’s decomposition (upper panel), with a free n (thick blue line), free α (thin cyan dots), and free β (medium red dash).

The lower panel shows the realtive difference between the different cases. The relative difference (lower panel) is bound from above by ∼1%. Additionally the Planck observation error bars are shown.

The conclusion is that we will need additional accurate data from smaller cosmic scales to be able to differentiate between the three scenarios. These extra e-folds might come from future missions such as Euclid [24], or from μ-type distortion data [25, 26].

5 Results

We apply the methods discussed in Section 3 to the degree six polynomial inflationary potentials that yield r = 0.01. We calculate the most likely coefficients and extract the resulting most likely polynomial inflationary potential, the results are given in Tables 2 and 3. The PPS resulting from this inflationary potential is then calculated in order to confirm that the most likely coefficients reconstruct the most likely observables.

Table 2

A comparison of recalculated power spectra observables from results of the two extraction methods.

Observable	Recalculated
	Probabilistic Method	Multi-fit
ns	0.9687	0.9684
αs	0.0089	0.0076
βs	0.0176	0.020

Table 3

The most likely coefficients extracted by the process of likelihood assignment and marginalization, as well as by using the multinomial method.

	Gaussian extraction		Multinomial fit
r = 0.01	μ (average)	σ (standard deviation)	value
a2	0.0402	0.0156	0.01866
a3	−0.152	0.0414	−0.0235
a4	0.0215	0.1123	−0.3452

5.1 Results for degree six polynomials that yield r = 0.01

In Fig 2 we showed a cover for the joint likelihood map of n − α, of about 2000 potentials with r = 0.01. The cover is approximately uniform, thus we were able to assign likelihoods to every potential we study, as previously discussed.

By a process of marginalization, as discussed in Section 3, we extract the most likely coefficients, which yield the likeliest observables. This process is represented graphically in Fig 3. The results are shown in Table 3 and the PPS reconstruction is shown in Fig 5. The advantage of this method is that it also determines the deviation from the average value. This can be used as an indicator for the level of tuning that is required to construct the most likely small field model. A discussion of tuning in field theoretic models can be found in [27], as well as in [28] and [29]. In most cases the tuning level can be viewed as simply , which in this case are given by (0.375, 0.27, 5.5) for a2, a3, a4 respectively. The width of the Gaussian fits for {a2, a3, a4} are {0.015, 0.041, 0.112} respectively This is shown in Fig 6. These widths represent the effective area in parameter space that yields observables within the 68% CL. Which is another measure of the tuning required in the sixth-order polynomial models.

Fig 5

Reconstruction of the PPS from the most likely potential with r = 0.01, as calculated by the multinomial (reverse fitting) method (X’s and red dash), as well as the probabilistic method (circles and blue line).

The CMB observables are well within the 68% confidence levels of the MCMC analysis for both. However, the probabilistic method seems to yield more precise results.

Fig 6

The calculated likelihoods for the coefficients {a2, a3, a4}, in models with r = 0.01.

The most likely coefficients are given by: a2 = 0.04, a3 = −0.15, a4 = 0.02. The tuning level for each coefficient is given by Barbieri-Giudice measure and is (0.375, 0.27, 5.5).

Recalculating the CMB observables that this most likely model yields, we find n = 0.9687, α = 0.0089, β = 0.0176. These values are very close to the most likely values found from the previously discussed MCMC analysis of the BICEP2+Planck data. The resulting scalar index fits the most likely value in Table 1 exactly, while α and β deviate from these values by no more than 12.5%. We found that this is a relic of the binning method. Adding more models to the simulated data and refining the binning process results in even better proximity to the desired values.

Using the method of multinomial evaluation (3.2), we found the multinomial coefficients for each of the model degrees of freedom. For instance for a2 we have: Since n ∼ 1, and α and β are of the order of 10−2, the above result suggests that a2 is primarily dominated by n. Similarly, we have found that a3 is dominated by a linear combination of α, and β, and a4 is primarily dominated by β. This method yields the most likely CMB observables with comparable accuracy to the previous method upon recalculation: n = 0.9684, α = 0.0077, β = 0.020. Table 2 contains a comparison between the results of the two methods.

5.2 Most likely potentials

Since n is better constrained, we opt for the analysis that yields a more precise value of n. The leading 6th degree polynomial which yields r = 0.01 at the proper pivot scale, is thus given by: Upon initial investigation, it seems these models produce a relatively flat tensor power spectrum. This might motivate future research of the tensor power spectrum predictions and constraints.

6 Observable dependence

An interesting finding is an inter-dependence of the three observables n, α, β. For models that yield r = 0.01, there is a quadratic relation between the observables, such that β = β(n, α) this is evident in Fig 7. It should be stressed that this is a phenomenon associated with the models and not with the observational data. This is supported by the small n, α, β paired-covariance found in [11], implying weak dependence among observables in the data itself.

Fig 7

Dependence of β on the other observables, exposes an approximate quadratic relations between α and β.

The width of the resulting band indicates the deviation from a quadratic relation, which is correlated to n.

One might think that the previous findings in [16, Figure 23] indicate that n and α are dependent. However, the graph shows a dependence between α and n which is the scalar index evaluated at k = 0.002 hMpc−1. Taking some initial n evaluated at k0, it follows that n evaluated at some other scale, depends on the index running α. Indeed, when one examines the color coding in [16, Figure 23], which represents n at the pivot scale, it is clear that n and α are independent.

7 Summary and outlook

A large sample of potentials that yield r = 0.01 and conform to the allowed observable values was successfully generated. The sample provided a uniform cover of the allowed region of parameters which enabled us to assign likelihoods to each of the potentials and extract each coefficient’s likelihood (3.1). Another approach was implemented, representing each coefficient as a multinomial function of the observables (3.2), which yielded similar results. A most likely small field potential giving rise to r = 0.01 and (n ≃ 0.9694, α ≃ 0.009, β ≃ 0.0175), was identified, and its power spectrum simulated. An interesting inter-dependence of (n, α, β) was found in these models, which may have some bearing on future MCMC analysis. We hope to perform such an MCMC analysis, either with priors that include this dependence or with a numerical scheme that reflects it.

The Planck collaboration may soon release additional analysis products, and BICEP3 is also expected to release results in the near future. Thus, it might be possible to check our prediction for the tensor-to-scalar ratio. However, ruling out models of the class discussed in this paper might be a more difficult task due to the lack of constraining power of current observations in the range l > 1500. We expect that either S4 cosmology or μ-distortion data will be able to resolve this in the foreseeable future by adding observational data on smaller scales.