Bayesian Monte Carlo Evaluation of Imperfect ( n , 233 U) Data and Model

. Conventional nuclear data evaluation methods using generalized linear least squares make the following assumptions: prior and posterior probability distribution functions (PDFs) of all model parameters and data are normal (Gaussian); the linear approximation is su ﬃ ciently accurate to minimize the cost function (even for nonlinear models); the model (e.g., of neutron cross section) and experimental data (including covariance data) are without defect and prior PDFs of parameters and measured data are known perfectly . Neglect of covariance between model parameters and measured data in conventional evaluations contributes to imperfections. These assumptions are inherent to the generalized linear least squares minimization method commonly used for resolved resonance region neutron cross section evaluations but are often not justiﬁed due to the presence of non-normal PDFs, nonlinear models (e.g., R -matrix formalism), and inherent imperfections in data and models (e.g. imperfect covariance data). Here, these assumptions are removed in a mathematical framework of Bayes’ theorem, which is implemented using the Metropolis-Hastings Monte Carlo method. Most importantly, new parameters are introduced to parameterize discrepancies between the theoretical model and measured data to quantify judgement about discrepancies or imperfections in a reproducible manner. An evaluation of 233 U in the eV region using the ENDF-B / VIII.0 library and transmission data (Guber et al.) is presented, and posterior parameters are compared to those obtained by conventional evaluation methods. This example illustrates the e ﬀ ects of removing the most harmful assumption: that of model-data perfection.


Introduction
Conventional nuclear data (ND) evaluation tools using generalized linear least squares (GLS) such as SAMMY [1] typically incorporate the following common assumptions: 1. The theoretical model perfectly describes the experimental data, and the prior probability distribution function (PDF) of model parameters and data is known perfectly.
2. Model parameters and data obey a normal (Gaussian) joint-PDF.
3. The linear approximation is sufficiently accurate to determine the posterior PDF, even for nonlinear models. * e-mail: brownjm@ornl.gov * * e-mail: arbanasg@ornl. gov This manuscript has been authored by UT-Battelle, LLC, under contract DE-AC05-00OR22725 with the US Department of Energy (DOE). The US government retains and the publisher, by accepting the article for publication, acknowledges that the US government retains a nonexclusive, paid-up, irrevocable, worldwide license to publish or reproduce the published form of this manuscript, or allow others to do so, for US government purposes. DOE will provide public access to these results of federally sponsored research in accordance with the DOE Public Access Plan (http://energy.gov/downloads/doe-public-access-plan).
The first assumption may be the most significant contributor to unrealistically small uncertainties of evaluated cross sections, which are propagated to neutron transport applications. Evaluators must often adjust evaluated ND uncertainties to provide realistic estimates of uncertainty in reconstructed cross sections and applications.
The new Bayesian framework presented here eliminates all three conventional assumptions and provides evaluators with a tool to quantitatively incorporate their expert judgment during evaluations in a scientifically repeatable way. By eliminating the first assumption, this framework may effectively appear similar to methods used to quantify the model defect [Schnabel et al., e.g. [2]]. However, as stated above, the first assumption is more general because it recognizes another aspect: namely, the imperfections inherent to the (generalized 1 ) data and their covariances. These data contain imperfections which likely exist even when a model happens to be perfect. One example is the conventional tendency to neglect prior covariances between model parameters and the measured data, which could be seen as one of several other potential sources of imperfections in the generalized data. 1

Theory
To simplify the derivation, a generalized data vector [3], z, is used to denote the concatenation of a vector of model parameters, P, and a vector of data, D, that is: In this work we restrict Bayesian prior PDFs to be of the normal kind, 2 namely, as defined by the expectation value z , computed as and by its covariance matrix, C, computed analogously as as expected for normal PDFs. To derive an expression for the Bayesian posterior PDF without making Assumption 1, it is useful to define an auxiliary variable, δ, as as well as its covariance matrix, where T (P) is a theoretical model 3 . With these definitions, it can be demonstrated that making Assumption (1) amounts to constraining the posterior expectation values of δ and its posterior covariance matrix to be identically equal to zero. This suggests that the removal of Assumption (1) amounts to setting the posterior expectation values of δ and its covariance matrix to some arbitrary finite constant values, namely, δ and ∆ , which are hereby introduced as a new tool for evaluators to quantify their confidence in the quality of (generalized) data or the validity of the theoretical model. In this notation, a prime symbol upon an expectation value indicates that a posterior PDF (also indicated by a prime symbol) was used to compute that expectation value in contrast to the (non-primed) prior expectation values and the (non-primed) prior PDF. With these definitions, the Bayesian posterior PDF becomes where the likelihood function can be expressed as [4,5] where λ and Λ are the parameters determined by constraining the posterior expectation value of δ and its covariance matrix to be equal to the δ and ∆ , respectively. For linear models, expressions for λ and Λ in terms of δ and ∆ can be found analytically, whereas for nonlinear models, they could be found using an iterative algorithm presented in [4]. The posterior PDF could then be simplified by adding the exponent in the likelihood function to the one in the prior PDF, resulting in where is a generalized cost function. The conventional cost function, χ 2 (P), is recovered by setting δ = 0 and ∆ = 0, yielding λ = 0 and Λ = 0, so that the likelihood function effectively becomes a Dirac delta function, δ Dirac (T (P) − D). This allows data, D, to be integrated out of the posterior PDF. This integration over data, D, yields a (conventional) posterior PDF of parameters P alone, where is the conventional cost function, z| D=T (P) ≡ P T (P) (13) has been defined for convenience, and where z and C are defined by Eqs. (3) and (4), respectively. Often conventional evaluation methods then make Assumption (2), namely that the posterior PDF is a normal one, by approximating the posterior expectation values of parameters, P , by the values that minimize χ 2 (P), and by approximating the posterior covariance by the Hessian of the χ 2 (P) evaluated at the minimum. Finally, Assumption (3) is introduced with with the use of first-order derivatives during minimization of χ 2 (P) by the iterative Newton-Raphson method.

Evaluation of low energy RRR of 23U
The framework described above can be implemented several ways. Simply incorporating the cost function given in Eq. (10) into an existing GLS evaluation code would remove the assumption of "perfection" described above (note that this is not done in this work). To remove the assumptions of PDF normality and model linearity, simple Monte Carlo (MC) methods can be applied. This work used the well-known Metropolis-Hastings (MH) algorithm [6][7][8], a subset of Markov Chain MC methods. The combination of MC methods and the new generalized cost function effectively removes the three assumptions listed in the introduction.
Although this framework was designed to work with arbitrary models, the primary goal in this effort was to evaluate model parameters for a cross section model in the RRR, where the Reich-Moore [9] approximation to Rmatrix theory [10] is implemented in the code SAMMY. Because it is the most current, the ENDF/B-VIII.0 [11] evaluation of 233 U is used to obtain the starting mean values and uncertainties for resonance parameters P. D is set to the data produced by Guber et al. [12], and the covariance C of the resulting hypervector z is a combination of ENDF/B-VIII.0 variance and variance reported by Guber et al. For simplicity in this example, all offdiagonal elements of C are equal to zero, even though the posterior PDF will result in a posterior covariance C with non-zero off-diagonals. Because it is known a priori that the cross section model and experimental data cannot perfectly match (because of data or model defects), small values (0.01) have been assigned to the diagonal of Λ to incorporate that information into the likelihood function of δ. In this example, the posterior expectation value of δ, δ , is zero, therefore λ is set to zero in the likelihood function.

Uncertainty on Model-Data Residuals
To quantify the expectations of the difference between models and data, all elements of λ = 0, the diagonal elements of Λ = 0.01, and off-diagonal Λ were set equal to zero. Setting λ = 0 is essentially stating that the model should match the data perfectly, but by setting the diagonal of Λ equal to 0.01, the PDF of the likelihood becomes finite, allowing imperfections. This means that small variations of δ away from zero will be penalized a small amount (within the variance of 0.01), and large variations will be penalized by more significant cost. Six resonance parameters were varied for three resonances: the energy E i and neutron width γ n,i for each. Starting with these prior data, the MHMC Bayesian algorithm is iterated until convergence of a posterior parameter set. The posterior mean cross section model T (P) is shown in Fig. 1 compared to the mean of prior models and GLS posterior (SAMMY). The lower plot in the figure shows the relative standard deviation of the model of posterior parameters.
As seen in the lower plot ( Fig. 1) of relative uncertainty in transmission, the uncertainties propagated from the parameters to the model are proportional between the GLS and MHMC solutions (of course, this ignores all off-diagonal covariance in T(P)). The GLS method implemented in SAMMY predicts uncertainty in the cross section roughly an order of magnitude smaller than that seen from the MHMC posterior and the variance of data D , which they both assessed. The magnitude of the energydependent uncertainties in T (P) is related to the values of Λ: as Λ → ∞, it approaches the uncertainty propagated from the prior, and as Λ → 0, it approaches that of the GLS method. The drastic nature of the reduction in uncertainty can perhaps be seen more clearly in the resonance parameter PDFs shown in Fig. 2. As each of the PDFs plotted along the diagonal of the plot-matrix are normalized to integrate to one, the height of the posterior GLS computed by the MHMC algorithm (blue) are quite close to those computed by SAMMY (red) using the GLS method. The energydependent uncertainty computed by these two methods, however, differs by an order of magnitude (lower plot). Energy-dependent uncertainty is determined by the square root of the variance of the posterior models of P. These expectation values and uncertainties are compared to those of the prior data by Guber et al. [12] (black) and to the prior expectation value of models described in ENDF/B-VIII.0 [11]. PDF must be cut off to even see the prior and MHMC posterior PDFs. Sometimes the GLS and MHMC mean values shift the same direction from the prior, but the GLS posteriors are consistently and drastically more narrow because of the inherent assumptions made by GLS, the most significant assumption being that D = T (P) for posterior D and P.

Criticality for 233 U
One advantage of storing explicit PDFs of posterior resonance parameters is that they can be directly used to calculate PDFs of quantities of interest in neutron transport calculations (e.g., k e f f ), effectually propagating uncertainties from resonance parameters directly to applications.
Continuing the 233 U example above, random samples are drawn from the posterior PDFs of resonance parameters from the generalized BMC and GLS evaluations, and 233 U cross sections are reconstructed each time and used to simulate a criticality benchmark for 233 U (U233-SOL-INTER-001 in Ref. [13]). This approach involves reconstructing 1,000 cross sections from each PDF with AMPX [14] and then running 1,000 eigenvalue calculations with the KENO-VI code within the SCALE code system [14]. The PDF of k e f f , assuming the posterior PDF of P from the GLS evaluation, has a standard deviation of ∼93 pcm. The PDF from the MHMC evaluation, in which the evaluator has incorporated explicit uncertainty information in the model-data residual, is more conservative with a standard deviation of ∼490 pcm. The experimental uncertainty in the benchmark is ∼850 pcm.

Conclusions
Conventional methods of nuclear data evaluation make inherent assumptions about the a priori and a posteriori PDFs, occasionally leading to underestimation of uncertainties or the need for ad hoc methods to properly estimate uncertainty. The framework presented here attempts to remove these assumptions and provides a method which can explicitly define the evaluator's belief in the likelihood function. The most deleterious assumption is that of perfection in model-data PDFs, which in turn makes the removal of this assumption with a new likelihood function the most impactful contribution of this work. In practice, this means that evaluators can define model defects, data defects, and other expert judgements in a welldocumented, repeatable way. This information is essential for creating nuclear data libraries that properly estimate covariance.