Analyzing rotational bands in odd-mass nuclei
using Effective Field Theory and Bayesian methods

Our goal in this analysis is to use the information on the expected size of LECs to stablize the extraction of their values as we add more levels to the analysis, or as we use energy-level formulae computed at different EFT orders. At the same time, we want to estimate the expansion parameter, $Q$, of the theory, as well as the characteristic sizes for even and odd coefficients, $\bar{c}_{even}$ and $\bar{c}_{odd}$.

We want to obtain the posterior distribution for all the LECs that appear at order $k$, a set we collectively denote by ${\bf a}_{k}$. Here we will obtain the joint posterior pdf of ${\bf a}_{k}$, the expansion parameter $Q$, and the characteristic sizes. To do this we follow the successful endeavor by the BUQEYE collaboration in Refs. [33, 34, 35], and write the posterior, given experimental data, $\vec{y}_{exp}$, and prior information on the model, $P*$, as

Marginalization of this posterior distribution over $Q$, $\bar{c}_{even}$ and $\bar{c}_{odd}$ yields the posterior distribution for ${\bf a}_{k}$. Other marginalizations can be carried out to obtain posteriors for $Q$, $\bar{c}_{even}$ and $\bar{c}_{odd}$.

This joint posterior distribution tells us the probability of the LECs and the error model parameters given experimental data. We could use this posterior distribution to get other quantities or observables, such as the energy of a particular rotational level, which depend on the LECs or the error model parameters. These are now represented by distributions and not single numbers. Their distributions are called posterior predictive distributions (PPD). We write the PPD of an observable $\mathcal{O}$ as

where $\vec{\theta}$ represents the LECs and the error model parameters. Calculating the observable at each point in the parameter space $\vec{\theta}$ and then integrating over the parameters $\vec{\theta}$ allows one to carefully account for correlations between the parameters.

Using Bayes’ theorem, we can express the posterior (8) as

The terms on the right-hand side of Eq. (10) have the following interpretations:

1. 1.

   $\operatorname{pr}(\vec{y}_{exp}|{\bf a}_{k},Q,\bar{c}_{even},\bar{c}_{odd},P*)$ is the likelihood of the experimental data given specific values of both the LECs that appear in the energy formula at order $k$ and the parameters in our error model.

2. 2.

   $\operatorname{pr}({\bf a}_{k}|Q,\bar{c}_{even},\bar{c}_{odd},P*)$ is the prior distribution of the LECs given the parameters encoding the systematic expansion of the EFT.

3. 3.

   $\operatorname{pr}(Q|\bar{c}_{even},\bar{c}_{odd},P*)$ is the prior distribution of the expansion parameter given the characteristic sizes of even and odd natural coefficients.

4. 4.

   $\operatorname{pr}(\bar{c}_{even}|P*)$ and $\operatorname{pr}(\bar{c}_{odd}|P*)$ are the prior distributions of the even and odd characteristic sizes. (In Eq. (10) we assume an uncorrelated prior on $\bar{c}_{even}$ and $\bar{c}_{odd}$.)

5. 5.

   $\operatorname{pr}(\vec{y}_{exp}|P*)$ is the evidence, which we drop in what follows as it does not depend on the parameters we are interested in extracting and functions only as a normalization constant.

We now build the likelihood function accounting for the expected error between the experimental and theoretical values, for data on $K=1/2$ rotational bands. The corresponding likelihood for $K=3/2$ bands is built analogously. Following [34] we start by writing our observable (the energy of a particular rotation level) at order $k$ as

We choose the leading-order energy for each level, $A_{\rm rot}I(I+1)$, to be the reference scale $E_{\rm ref}$ for the observable. The dimensionless coefficients $c_{n}$ (see Eq. (7)) are assumed to be $\mathcal{O}(1)$. The theory error $\vec{\sigma}_{th}$ at any order is due to terms omitted from the summations in Eq. (11). Its most significant contribution comes from the first omitted term in the EFT expansion. Accounting only for this term yields an estimate for the theory error that is fully correlated across levels if $k+1$ is even, and anticorrelated for adjacent levels if $k+1$ is odd. To account for this correlation or anticorrelation we write the theory covariance matrix as the outer product of a vector representing the theory error, $\Sigma_{th}\equiv\vec{\sigma}_{th}\otimes\vec{\sigma}_{th}$. The vector $\vec{\sigma}_{th}$ contains the value of the first omitted term for each of the $m$ energy levels that enter the likelihood. We also account for experimental errors by writing the covariance matrix as

where we take $(\Sigma_{exp})_{ij}\equiv(\vec{\sigma}_{exp})_{i}^{2}\delta_{ij}$. The likelihood function is then

where $\vec{r}\equiv\vec{y}_{exp}-\vec{y}_{th}$ is the residual between the central experimental energy for a level and the theory result (11) and $m$ is the number of levels included in the likelihood estimation.

We note that since the theory error is the outer product of the theory error with itself, the theory covariance $\Sigma_{th}$ is singular. Including the experimental error solves this singularity problem for the covariance $\Sigma$. However, $\Sigma$ can still become ill-conditioned for higher values of $Q$ if the experimental errors are too small; numerical issues then arise when we try to invert the covariance matrix.

Including more terms in the estimate for the theoretical error produces a steeper peak in the likelihood function, see Fig. 1, which, in turn, restricts the values sampled for $Q$ to a narrower region. Because it precludes the sampler exploring large values of $Q$, this inclusion of more omitted terms in the model of the theoretical error solves the numerical problem of ill-conditioned matrices and gives a more accurate extraction of the LECs and the error-model parameters.

In what follows we estimate the theory error including omitted terms up to a certain cutoff order $k_{max}$. Our theory error estimate for the level with spin $I$ is then

where the $\bar{c}$ that is used here is $\bar{c}_{even}$ for even values of $l$ and $\bar{c}_{odd}$ otherwise. The $I$-dependence of the $l$th term is chosen to match that in Eq. (11), and is denoted here by $P_{l}(I)$, a polynomial of power $l$. We arrange the contributions to the theory error, (14) as the $p$ columns of a $m\times p$ matrix $\sigma_{th}$, where $p=k_{max}-k$ is the number of omitted terms. Each column in this matrix then corresponds to the theory-error structure, while each row corresponds to a different energy level. To obtain $\Sigma_{th}$ we then again take the outer product of $\sigma_{th}$ with itself, i.e., we construct an outer product in our $m$-dimensional data space, while also taking an inner product in order space. This results in the theory error associated with different orders being added in quadrature, while maintaining the correlation structure of the theory error across the data space.

The prior distributions for an order-$n$ LEC is taken to be a Gaussian with mean zero and standard deviation

encoding the EFT expectations for the sizes of the LECs arising from the power counting described in Sec. II. The standard deviation in Eq. (15) allows the possibility for even and odd LECs to have different typical sizes. Combining the Gaussian priors for the LECs yields

The LEC $E_{K}$ is just an energy shift and its size is not determined by the EFT power counting. We set the prior on it to be Gaussian with mean zero and a standard deviation, $\bar{E}$, that is wide enough to capture its value. The value for $\bar{E}$ is determined from the energy of the bandhead and $A_{\rm rot}$ by means of Eq. (3).

We choose not to impose any expectations regarding the size of expansion parameter in the prior for $Q$ and so take it to be flat between two limits:

Limiting $Q$ from above restricts the sampler from going to high values of $Q$, as they make the covariance matrix ill-conditioned and harder to invert. For all cases we check that the posterior for $Q$ is confined to values well below $Q_{cut}$.

The priors on the characteristic sizes $\bar{c}_{even}$ and $\bar{c}_{odd}$, are taken to be identical scaled-inverse-$\chi^{2}$ distributions

where the cutoff $\bar{c}_{cut}$ prevents numerical issues inverting the covariance matrix. The scaled-inverse-$\chi^{2}$ distribution, given by

is shown for different values of $\nu$ and $\tau$ in Fig. 2.

We stress that we chose identical priors for $\bar{c}_{even}$ and $\bar{c}_{odd}$ even though we expect the former to be larger than the latter based on previous analyses of data on rotational bands [19]. We did not want to bias our analysis by imposing this hierarchy on the prior, instead anticipating that it will emerge naturally in the posteriors for those parameters.

The scaled-inverse-$\chi^{2}$ favors small values of $\bar{c}^{2}$ and has long tails. This allows the sampler to explore higher values of $\bar{c}^{2}$. The sharp decrease in this distribution for very small values of $\bar{c}^{2}$ could be a problem for cases where $\bar{c}_{odd}$ is much smaller than one. This is a concern in some $K=3/2$ bands where we expect smaller odd-order corrections to the leading-order energy than in $K=1/2$ bands.

In this subsection we show that lower-order LECs extracted for the selected rotational bands are stable across EFT orders and with the addition of high-energy data, provided that we account for enough omitted terms when treating the theory error. For most nuclei including omitted terms up to $k_{max}=10$, i.e., accounting for six omitted terms at ${\rm N}{}^{4}{\rm LO}$, was enough to stabilize the extraction of the LECs.

As an example, we show the stability of the extracted LEC, $A_{1}$, across number of levels included at different EFT orders in Fig. 4. In this figure, $I_{max}$ is the spin of the highest-energy level included in a particular analysis. The central values of the resulting posteriors are consistent with each other within 68% credible intervals, shown as error bars in the figure. Adding more levels to the analysis narrows the posteriors for the LECs up to a certain $I_{max}$, after which the widths of these distributions saturate. Fig. 4 also demonstrates striking agreement between the distributions obtained at low and high EFT orders: they are almost identical as long as omitted terms up to the same $k_{max}$ are accounted for in both analyses.

The importance of including more than one omitted term in the theory error estimate is evident in Fig. 5. The top and bottom panels of the figure show the way that posteriors for $B_{1}$ and $B$ evolve as $I_{max}$ increases. This is done using three error models that include different numbers of omitted terms. These results show that including more omitted terms in the model of the theory error removes the drifting and staggering of the central values.

For both cases the distributions at $k_{max}=10$ agree within errors as we go higher in $I_{max}$. The narrowing of the distribution as we go higher in $I_{max}$ is clearly seen in those two figures. In addition to having less data, the broadening of the error bands at low $I_{max}$ comes from the fact that including less levels in the analysis leads to highly correlated LECs. This allows the numerically larger errors on the lower-order LECs to contribute to the errors on the higher-order LECs, thereby enhancing them.

In Fig. 6 we show the decrease in the correlations between the LECs as $I_{max}$ increases. The high correlation between the LECs at low $I_{max}$ occurs because these analyses do not include enough data to constrain all LECs independently. Furthermore, the high correlation between the LECs at low $I_{max}$ also results in an unreliable extraction of the expansion parameter $Q$. This comes from the fact that at low energies the theory truncation error is very small compared to the experimental error. Indeed, adding more terms to our EFT error model (i. e., increasing $k_{max}$) leads to higher correlation between the LECs at low $I_{max}$. Thus, the number of levels required to reliably extract $Q$ increases with increasing $k_{max}$.

Starting instead at the low-$I$ end of the data: when we progressively remove the lowest-energy levels from the data set $D$ used to construct the likelihood we rapidly lose the ability to reliably extract the LECs. Figure 7 shows that the distribution for $E_{K}$ starts narrow and broadens as we remove levels from below. When we remove the six lowest energy levels the distribution of $E_{K}$ is exactly the same as the prior distribution: the likelihood is making no contribution to the $E_{K}$ posterior.

The previous results were nearly the same for all cases considered in this work. However, even for $k_{max}=10$, staggering and shifting of the LECs remains sizable for the $K=1/2$ bands in ¹⁸³W, ¹⁶⁷Tm and ²³⁵U. In ¹⁸³W and ¹⁶⁷Tm, these effects could be attributed to large expansion parameters, as they translate to large omitted contributions to the energies of the rotational levels. In ²³⁵U, the fermionic matrix elements could be larger than naively expected, causing the systematic expansion of the EFT to be questionable as discussed in Ref. [19].

For $K=3/2$ bands, we were able to extract stable LECs from the ¹⁵⁹Dy analysis by setting $I_{max}=15.5$. This extraction required us to consider omitted terms up to $k_{max}=12$. This is because the spin at which the EFT breaks in this nucleus is $I_{br}\approx 15.5$. (This, then, is the third case in which we do not use all the NNDC energy-level data available on a particular band.) ¹⁵⁷Gd is stable across orders and $I_{max}$, while ¹⁵⁵Gd exhibits shifting and staggering due to a larger expansion parameter, $Q\approx 0.07$.

The values of the LECs and the error-model parameters at ${\rm N}{}^{4}{\rm LO}$ for the nuclei considered in this work are given in Tabele 1.

In addition to using the scaled-inverse-$\chi^{2}$ distribution as a prior for $\bar{c}_{even}$ and $\bar{c}_{odd}$ we tried truncated Gaussians with mean zero and standard deviations $\sigma=7$ and $\sigma=3$ respectively for all cases. These truncated Gaussian priors allow for smaller values of the characteristic sizes. But the standard deviations were chosen to still allow values for $\bar{c}_{odd}$ and $\bar{c}_{even}$ larger than those resulting from scaled-inverse-$\chi^{2}$ priors with $\nu=1$ and $\tau^{2}=1$.

The change in prior for $\bar{c}_{odd}$ and $\bar{c}_{even}$ does not significantly change the posteriors for the LECs: the corresponding central values differ by less than 1%, and are consistent with each other within the 68% credible intervals. Central values of the posteriors for $Q$ differ by less than 15%, and were similarly consistent.

The strongest dependence on the prior is that exhibited by the posteriors for $\bar{c}_{odd}$ and $\bar{c}_{even}$: the central values differ in some cases by more than 50%. However, even these values are consistent with each other within 68% credible intervals, since the posteriors for the characteristic sizes are broad.

The changes in the posteriors of $Q$ on one hand and $\bar{c}_{odd}$ & $\bar{c}_{even}$ on the other are anticorrelated. We only care about combinations of them to set the size of the theory error and the expected size of the LECs. Thus, the dependence of the theory error and the expected size of the LECs on the prior for the characteristic sizes is less profound. The difference in sizes of the theory error resulting from the chosen priors is less than 20% for all cases except ¹⁵⁷Gd, where the difference is about 40%.

For all results that follow we used the scaled-inverse-$\chi^{2}$ distribution with $\nu=1$ and $\tau^{2}=1$ as the prior for both $\bar{c}_{odd}$ and $\bar{c}_{even}$, in keeping with the naturalness assumption.

Figure 8 shows the PPDs of the energy residuals as a function of the spin $I$ for two cases considered in this work. These distributions are calculated using Eq. (9). In each figure, translucent blue lines connect energy residuals resulting from different LECs sets sampled from the posterior distribution in Eq. (8). The solid black line represents the median of the PPD, and the dashed lines encompass the region between the 16th and 84th percentiles. The dark and light red bands show the truncation error and the experimental error added in quadrature at 68% and 95% credible levels respectively. To calculate the truncation error, we consider a theory error that accounts for six omitted terms. They are added in quadrature, just as they are in the likelihood defined in Sec. III. This calculation was done using the median values of $Q,\bar{c}_{even}$, and $\bar{c}_{odd}$. The dependence of the size of the theory error on the prior on $\bar{c}_{even}$ and $\bar{c}_{odd}$ is small in these cases: the theory error changes by about 10% when the prior is changed.

The correlation coefficient written in the legend in Fig. 8 is the largest between any LEC and any error-model parameter for the shown analysis. When this value is small, the truncation error and the propagated LEC error could in principle be added together in quadrature.

In viewing Fig. 8 it is important to remember that the truncation error on the energy residuals is highly correlated across levels. This comes from the high correlation between levels when building the correlation matrix that goes into the likelihood. This correlation also flows into a correlation between levels in the PPD of the energies. A correlation plot between two energy levels, like the ones in Fig. 9, gives a 2D cut of this multi-dimensional correlation.

In both panels we see the importance of accounting for more than one omitted term in the theory error. This is clearly shown in the reverse in the direction of the correlation from a negative to a positive correlation when going from the orange ellipse to the red ellipse. The orange ellipse is obtained when we account for only one omitted term, while the red ellipse includes the effect of two omitted terms. After accounting for six omitted terms the green ellipse is obtained and the 68% ellipse in principle expands. This is more clearly seen when we go to high-energy levels plotted in the lower panel in Fig. 9. Note also that for lower-energy levels the correlation is smaller since the experimental error dominates over the truncation error, and we assumed that the experimental errors are not correlated across energy levels.

In Figs. 10, 11 and 12 we compare the marginalized posterior distributions of the LECs, on the y-axis, with their expected sizes from the EFT power counting, on the x-axis. Since we also extract the theory error parameters from the sampler and they are highly correlated among themselves, we calculate the expected size from the distributions of the error model parameters using Eq. (9). We notice that the error on the distribution of the LECs is very small compared to the error on the expected sizes that comes from the distribution of the theory error parameters.

As these graphs are model-checking graphs, and since the estimates of LEC sizes plotted on the x-axis are meant as order-of-magnitude estimates, we do not expect perfect linear correlations. Nevertheless, Fig. 10 shows that, for all $K=1/2$ bands considered, the LEC $\Delta A$ agrees with its expected size within error bands. This result is surprisingly better than expected. In contrast, the size of $\Delta A$ for $K=3/2$ bands is larger than expected, especially for ¹⁵⁵Gd (see blue points in Fig. 10). There are two factors that could contribute to this. First, the $K=3/2$ bands have larger fermionic matrix elements. This could hinder the systematic expansion of the EFT. Second, the $K=3/2$ bands have relatively larger expansion parameters, see Fig. 13.

The same discussion applies to the results in Figs. 11 and 12, where we see good agreement between the LECs and their expected sizes for $K=1/2$ bands. The disagreement with power-counting estimates for $K=3/2$ bands at $\mathrm{N}^{3,4}\mathrm{LO}$ is less of a concern than the one at $\mathrm{N}^{2}\mathrm{LO}$ seen in Fig. 10, since these higher-order LECs are smaller than their expected sizes. This doesn’t undermine the convergence of the EFT expansion.

We also note here that the scale of the x-axis is prior dependent and could change by more than 50% in some nuclei, depending on the choice of prior on $\bar{c}_{even}$ and $\bar{c}_{odd}$. For ¹⁵⁷Gd changing the prior on $\bar{c}_{even}$ and $\bar{c}_{odd}$ to a truncated normal allowed for smaller values of $\bar{c}_{odd}$ and $A_{3}$ was then equal to the expected size (i. e.,the point for ¹⁵⁷Gd then falls exactly on the line in figure 11). This did not happen when the truncated normal is chosen as a prior for the analysis in ¹⁵⁵Gd and ¹⁵⁹Dy; this may occur because there is strong $\mathrm{N}^{5}\mathrm{LO}$ energy-level staggering present in the data for these nuclei.

The size of $\bar{c}_{even}$ and $\bar{c}_{odd}$ is constrained by both the sizes of the LECs and the size of the theory error. In a good systematic expansion the tension between those factors on setting the size of the $\bar{c}$’s would be small and one number apiece would suffice to represent the even and odd order corrections. However when the systematic expansion is hindered, as in the case for $K=3/2$ bands due to large fermionic matrix elements, this tension becomes clear. One example of this is seen in Figs. 10 and 12 for $K=3/2$ bands. There $\Delta A$ is large and favors large values of $\bar{c}_{even}$, however, $B$ is small and favors smaller values of $\bar{c}_{even}$. The eventual result is a compromise. This tension may be exacerbated by the truncation error also providing information on the size of the $\bar{c}$’s.

In Fig. 13 we see a clear correlation between the extracted values of $Q$ and those that are expected based on each nucleus’ single-particle and vibrational energy scales, $E_{\rm sp}$ and $E_{\rm vib}$. The expected $Q$ is the larger of $E_{\rm rot}/E_{\rm sp}$ and $E_{\rm rot}/E_{\rm vib}$, while the extracted $Q$ comes from sampling the posterior in Eq. 10. This extracted $Q$ is what actually determines the convergence of the EFT expansion. It is markedly smaller than would be naively expected. The break-down scale of the theory is thus higher than naively expected: our rotational EFT works to much higher $I$ than energy-scale arguments would suggest. This could occur because coupling between the higher rotational states explicitly included in the EFT and the high-energy states that are not explicitly included in our EFT is hindered by the large difference in angular momentum between them.