Bayesian approach to model-based extrapolation of nuclear observables

Background: The mass, or binding energy, is the basis property of the atomic nucleus. It determines its stability and reaction and decay rates. Quantifying the nuclear binding is important for understanding the origin of elements in the universe. The astrophysical processes responsible for the nucleosynthesis in stars often take place far from the valley of stability, where experimental masses are not known. In such cases, missing nuclear information must be provided by theoretical predictions using extreme extrapolations. To take full advantage of the information contained in mass model residuals, i.e., deviations between experimental and calculated masses, one can utilize Bayesian machine-learning techniques to improve predictions.

Purpose: To improve the quality of model-based predictions of nuclear properties of rare isotopes far from stability, we consider the information contained in the residuals in the regions where the experimental information exist. As a case in point, we discuss two-neutron separation energies $S_{2n}$ of even-even nuclei. Through this observable, we assess the predictive power of global mass models towards more unstable neutron-rich nuclei and provide uncertainty quantification of predictions.

Methods: We consider 10 global models based on nuclear density functional theory with realistic energy density functionals as well as two more phenomenological mass models. The emulators of $S_{2n}$ residuals and credibility intervals (Bayesian confidence intervals) defining theoretical error bars are constructed using Bayesian Gaussian processes and Bayesian neural networks. We consider a large training dataset pertaining to nuclei whose masses were measured before 2003. For the testing datasets, we considered those exotic nuclei whose masses have been determined after 2003. By establishing statistical methodology and parameters, we carried out extrapolations toward the $2n$ dripline.

Results: While both Gaussian processes and Bayesian neural networks reduce the root-mean-square (rms) deviation from experiment significantly, GP offers a better and much more stable performance. The increase in the predictive power of microscopic models aided by the statistical treatment is quite astonishing: The resulting rms deviations from experiment on the testing dataset are similar to those of more phenomenological models. We found that Bayesian neural networks results are prone to instabilities caused by the large number of parameters in this method. Moreover, since the classical sigmoid activation function used in this approach has linear tails that do not vanish, it is poorly suited for a bounded extrapolation. The empirical coverage probability curves we obtain match very well the reference values, in a slightly conservative way in most cases, which is highly desirable to ensure honesty of uncertainty quantification. The estimated credibility intervals on predictions make it possible to evaluate predictive power of individual models and also make quantified predictions using groups of models.

Conclusions: The proposed robust statistical approach to extrapolation of nuclear model results can be useful for assessing the impact of current and future experiments in the context of model developments. The new Bayesian capability to evaluate residuals is also expected to impact research in the domains where experiments are currently impossible, for instance, in simulations of the astrophysical $r$ process.

Figure 1

The experimental $S_{2n}(Z,N)$ datasets for even-even nuclei used in our study: AME2003 [48, 49] (light dots, 537 points), additional data that appeared in the AME2016 evaluation [50, 51] (dark dots, 55 points), and JYFLTRAP [52] (stars, 4 points).

Figure 2

Residuals of $S_{2n}(Z,N)$ for the six global mass models with respect to the testing dataset AME2003. The rms values of $\delta (Z,N)$ in MeV are marked: for AME2003 (upper number) and AME2003-H (lower number).

Figure 3

Similar as in Fig. 2 but for $S_{2n}(Z,N)$ residuals smoothed with the Gaussian folding function to emphasize long-range systematic trends.

Figure 4

Posterior distributions of the GP parameters, with the posterior mean and standard deviation listed.

Figure 5

Residuals of $S_{2n}(Z,N)$ for the six global mass models with respect to the testing dataset (AME2016-AME2003): $\delta (Z,N)$ (dots) and the GP emulator ${\delta }^{\mathrm{GP}}(Z,N)$ (circles).

Figure 6

Residuals of $S_{2n}(Z,N)$ for the six global mass models with respect to the testing dataset (AME2016-AME2003): $\delta (Z,N)$ (dots) and the BNN emulator ${\delta }^{\mathrm{BNN}}(Z,N)$ (circles).

Figure 7

Empirical coverage probability for the six models used in our study as functions of multiples of the standard deviation (s.d.). The reference curve corresponds to the Gaussian quantiles, i.e., the probability that new testing data fall into the corresponding CI according to the domain. That is, at 1 s.d., the reference curve is 0.68 (and 68% of the testing data should fall into the CI); at 2 s.d. it is 0.95, and so on.

Figure 8

Extrapolations of $S_{2n}$ for the even-even Sn chain calculated with the six global mass models with statistical correction ${\delta }^{\mathrm{em}}$ and one-sigma CIs obtained with GP($T+H$) and BNN($T+H$). Experimental (circles) and extrapolated (asterisks) values from AME2016 [51] are marked.

Figure 9

Extrapolations of $S_{2n}$ for the even-even Sn chain calculated with DD-PC1 with statistical GP($T+H$) and BNN($T+H$) approaches. One-sigma and 1.65-sigma CIs are marked.

Figure 10

Predictions and confidence intervals given by DD-PC1+BNN($T+H$) for the Sn chain for (a) two MCMC runs using 100 000 samples (after 10 000-sample tuning) and (b) two MCMC runs using 1 000 000 samples (after 100 000-sample tuning). The curve (1) in the top panel corresponds to the BNN($T+H$) curve of Fig. 8.

Figure 11

Posterior sample mean and standard deviation of $S_{2n}$ for several Sn isotopes predicted in DD-PC1. The statistical calculations were carried out with GP($T+H$) and BNN($T+H$) methods trained on the AME2016+JYFLTRAP dataset. The results are shown as a function of the number of MCMC samples. The number of samples used initially for tuning was 10 000 (panels a and b) and 100 000 (panel c).

Figure 12

Extrapolations of $S_{2n}(Z,N)$ for the Sn chain corrected with GP($T+H$) and one-sigma CIs, combined for four representative models. The different models are consistent overall once the statistical correction and uncertainty are taken into account.