Conception and software implementation of a nuclear data evaluation pipeline
Georg Schnabel, Henrik Sjöstrand, Joachim Hansson, Dimitri Rochman, Arjan Koning, Roberto Capote
CConception and software implementation of a nuclear data evaluation pipeline
G. Schnabel, ∗ H. Sjöstrand, J. Hansson, D. Rochman, A. Koning, and R. Capote NAPC–Nuclear Data Section, International Atomic Energy Agency, Vienna, Austria Uppsala University, Department of Physics and Astronomy, Uppsala, Sweden Paul Scherrer Institut, 5232 Villigen, Switzerland
We discuss the design and software implementation of a nuclear data evaluation pipeline appliedfor a fully reproducible evaluation of neutron-induced cross sections of Fe above the resolved reso-nance region using the nuclear model code TALYS combined with relevant experimental data. Theemphasis of this paper is on the mathematical and technical aspects of the pipeline and not on theevaluation of CONTENTS
I. INTRODUCTION 2II. MATHEMATICAL BUILDING BLOCKS 4A. Multivariate normal model 4B. Representation of covariance matrices 5C. Gaussian processes 6D. Levenberg-Marquardt algorithm with prior 8E. Taylor expansion of logarithmized pdf 9III. IT BUILDING BLOCKS 10A. Handling of information in the Bayesiancontext 111. Experimental data and statisticaluncertainties 122. Systematic components 123. Mapping of systematic components toexperimental data 12 ∗ [email protected]; previously affiliated to Uppsala University
4. Correlated systematic components 13B. Retrieval and handling of experimental data 151. Searching experimental data 152. Retrieving experimental data 16C. Mapping of model predictions to EXFOR 18D. Parallelization of nuclear model invocations 19IV. THE NUCLEAR DATA EVALUATIONPIPELINE 21A. Retrieval of experimental data 22B. Generation of predictions based on a referencecalculation 22C. Rule-based correction of experimentaluncertainties 23D. Correction of experimental uncertainties basedon statistics 24E. Evaluation of the Jacobian associated with thereference calculation 26F. Setup of Gaussian processes forenergy-dependent model parameters 26G. Optimization of TALYS parameters using theLM algorithm 27H. Calculation of a MVN approximation of the a r X i v : . [ phy s i c s . d a t a - a n ] S e p ..Nuclear Data Evaluation Pipeline G. Schnabel et al. posterior pdf 27I. Generation of representative random files 28V. DISCUSSION OF PRACTICALITIES ANDRESULTS 28VI. SUMMARY AND CONCLUSIONS 32References 33A. GENERAL STATISTICS 371. Properties of variance and covarianceoperator 37B. LINEAR ALGEBRA IDENTITIES TO SPEEDUP COMPUTATIONS 371. Derivative of an inverse matrix 372. Derivative of a logarithmized determinant 373. Matrix determinant lemma 374. Woodbury identity 385. Products involving the inverse of sparsecovariance or other positive-definite matrices 38C. BAYESIAN INFERENCE 381. Details on modified Levenberg-Marquardtalgorithm 382. Prior on second-derivative as smoothnessprior 393. Necessity of second-order Taylorapproximation of posterior pdf 404. Computation of the posterior covariancematrix 415. Computation of the posterior expectation ofinsensitive paramaters 42D. DOWNLOADING THE PIPELINE ANDPACKAGES 421. About identification strings 432. Catalogue of pipeline modules 43 I. INTRODUCTION
The ability to reproduce complex scientific workflowsbecomes more and more an important concern. Increas-ing awareness that many past studies are impossible toreproduce has now even led to the term reproducibilitycrisis [1]. Nuclear data evaluation is not exempt of suchconcerns. Many choices related to the selection of the ex-periments, corrections applied to the experimental data,and the choice of evaluation method are often not welldocumented enough to enable the reproduction of an eval-uation. The lack of traceability regarding decisions takenby human experts to produce an evaluation is a seriousobstacle to building on the efforts from the past. Timeand money needs to be spent on re-analyzing data and in-formation in articles, which could have been better spentotherwise. However, we also want to stress that availableknowledge in the form of evaluations, even if they arenot reproducible, can still be beneficially used as priorknowledge in future evaluations.Reproducibility and transparency is the standard prac-tice in everyday neutronics simulations performed in thenuclear industry. It may also be regarded as especiallyimportant in the evaluation of standard reactions, suchas the neutron standards [2–4]. These reactions serve asa reference in many nuclear experiments to translate rel-ative cross section measurements to absolute ones. Forthis reason, much effort coordinated by the IAEA andvarious nuclear data evaluation projects has been putinto the evaluation of these reactions. The documenta-tion of these efforts is important, not only in the form ofreports but also as scripts to enable the reproduction ofthese evaluations. The availability of the evaluations asa well-structured sequence of scripts facilitates studyingthe impact of specific assumptions on the evaluation andprovides a solid basis for improved evaluations, such asthose of neutron data standards, in the future.Awareness about the need of storing nuclear data in anaccessible and transparent way and the need for automa-tion of evaluation and verification processes is rapidly in-creasing as indicated by the ongoing efforts of the WPECsubgroup 49 on reproducibility in nuclear data evalua-tion and the recently approved WPEC subgroup 50 on acurated experimental database with a focus on program-matic readability. Another encouraging effort is under-taken by WPEC subgroup 45 to improve the data man-agement and processes to validate nuclear data libraries.Also established solutions exist that can facilitate thecreation of transparent and reproducible evaluations, suchas the correction system attached to the EXFOR web re-trieval system [5] and the
EXFOR-CINDA for application package available on the IAEA-NDS website enabling pro-grammatic access to information in the EXFOR library [6].Yet another effort to facilitate the retrieval and manage-ment of experimental data from the EXFOR library isthe conversion of the EXFOR library to a NoSQL JSONdatabase [7]. An example of the benefits of reproducibleevaluations is the TENDL library project [8]. Its focus2..Nuclear Data Evaluation Pipeline G. Schnabel et al. on automation and reproducibility from the beginninghas been a catalysator in the creation of a comprehen-sive nuclear data library, setting an example for futureevaluations.Recognizing the need for reproducibility, automationand transparency in nuclear data evaluation and moti-vated by the benefits of existing software and projects sup-porting those needs, we present a prototype of a pipelinefor nuclear data evaluation, which has been employed foran example evaluation of neutron-induced cross sectionsof Fe above the resolved resonance region.Even though the pipeline (and therefore also its prod-ucts) has still to be regarded as a prototype, it combinesseveral advancements of evaluation methodology achievedduring the last years. Examples are Gaussian process pri-ors on energy-dependent model parameters to addressmodel defects [9], the treatment of inconsistent experi-mental data by marginal likelihood optimization [10] andthe adjustment of model parameters via a prior-awareLevenberg-Marquardt algorithm [11] to take into accountnon-linearities of the nuclear physics model. All of theseapproaches are mathematically related and they can beused in different ways. Therefore it was important to im-plement these essential mathematical concepts and algo-rithms in a generic way so that at first glance seemingly un-related task settings can be dealt with in a unified manner.For instance, the adjustment of uncertainties of inconsis-tent experimental data and the adjustment of so-called hy-perparameters of Gaussian processes imposed on energy-dependent model parameters rely on the same mathemat-ical approach—marginal likelihood optimization.An additional important design aspect was to enablethe evaluation of a large collection of experimental data.By adopting an intuitive perspective on experimentaldata, model parameters and their relations in the formof a simple Bayesian network [12], we are able to applythe mathematical techniques on a large collection of ex-perimental data without the need for data reduction tech-niques as a preparatory step. We refer to the mathematicalconcepts and techniques used in the pipeline as mathe-matical building blocks and discuss them first.Besides the mathematical aspects, modern nuclear dataevaluation is also about the management of information,either in the form of experimental data or scientific knowl-edge distilled to nuclear models. Therefore aspects associ-ated with information management need to be addressedas well, such as the handling of experimental data and theparallelization of nuclear model calculations on a cluster.We refer to these needs as IT (=information technology)building blocks and they are discussed after the mathe-matical building blocks.After the presentation of the mathematical and IT build-ing blocks, we discuss the pipeline itself, which is a se-quence of scripts employing the building blocks to imple-ment the evaluation. We go through the individual steps,starting from the retrieval of experimental data from theEXFOR library, over the construction of experimentaluncertainty information using a rule-based approach, sub- sequent automated adjustment of systematic uncertain-ties associated with inconsistent experimental data, andending with the adjustment of TALYS model parametersand the generation of covariance and random files using amodified version of the TASMAN code. TALYS and TAS-MAN are part of the evaluation code system called T6 [13].Throughout the description of the individual steps, we dis-cuss features of the pipeline and possible directions forimprovement in the future. Practicalities, such as compu-tation time, and results are presented in a separate sectionafter the discussion of the steps in the pipeline.Importantly, the pipeline is open-source and has beenmade available online [14] and a Dockerfile is also avail-able [15] to facilitate the installation of all components ofthe pipeline as one bundle. The availability of the pipelineand the ability to launch it at the push of a button enablethe incremental implementation of more sophisticated al-gorithms in the future and the systematic study of theimpact of assumptions codified in the pipeline.The pipeline can therefore be regarded as a skeleton forfuture evaluations, which can be extended according to thescope of the evaluation. Going towards more global andautomated evaluations, manual choices regarding data se-lection and uncertainty specifications need to be removedand statistical algorithms made more robust. Going to-wards high-quality evaluations of specific isotopes, theselection of experimental data need to be more carefullyundertaken and assumptions regarding uncertainties andprior assumptions more scrutinized.In the actual implementation of the pipeline, we had tomake many specific choices concerning selection of exper-imental data, algorithms, and nuclear models. We stressthat many of the choices in our current implementationof the pipeline, such as the use of the nuclear modelscode TALYS, are not of fundamental significance for theconcept of the pipeline. Due to the modular design ofthe pipeline, specific algorithms or the model for fittingmay be replaced by something more pertinent in a cer-tain evaluation context. For this reason, we stress thatthis paper is as much about the concept of a pipeline andgeneral design considerations as it is about the specificimplementation details, which are not crucially linked tothe conceptual aspects of the pipeline.The paper can be read selectively depending on thereader. A reader interested in the mathematical frame-work and how its algorithms are applied in the evaluationpipeline may only read section II on mathematical build-ing blocks and section IV on the sequence of steps in thepipeline. A reader who is interested in the mathematicalframework and its translation to data structures to enablelarge scale computations may start with section II and sec-tion III A. Finally, a reader who is foremostly interestedin our ideas on the efficient management of experimen-tal data and model calculations in terms of informationtechnology may directly jump to section III.3..Nuclear Data Evaluation Pipeline G. Schnabel et al.
II. MATHEMATICAL BUILDING BLOCKS
Concerning the representation of information and al-gorithms, the pipeline relies on the functionality of theR package nucdataBaynet published on GitHub [16]. Inthis section we review the multivariate normal distribu-tion as it is the mathematical core concept of the pipeline,elaborate on the representation of information to enablelarge scale evaluations, explain the integration of Gaussianprocesses into this framework, elaborate on a prior-awareLevenberg-Marquardt algorithm for the optimization ofnon-linear nuclear physics models, and finally discuss anapproximation to the posterior distribution based on themultivariate normal distribution. All of these different con-cepts and algorithms have a clear link to the multivariatenormal distribution, which we strive to emphasize in thefollowing exposition.
A. Multivariate normal model
The mathematical core concept used at various placesthroughout the pipeline is the multivariate normal (MVN)distribution. If d denotes the number of variables, theMVN distribution is characterized by a center vector (cid:126)µ ofsize d and a covariance matrix Σ of dimension d × d . Theprobability density function (pdf) of a MVN distributionis given by N ( (cid:126)x | (cid:126)µ, Σ ) = 1 (cid:112) (2 π ) d det Σ × exp (cid:20) −
12 ( (cid:126)x − (cid:126)µ ) T Σ − ( (cid:126)x − (cid:126)µ ) (cid:21) . (1)Noteworthy, a MVN distribution can exactly capture lin-ear relationships as well as uncertainties, which are bothreflected in the specification of the covariance matrix.The variables associated with elements in (cid:126)x can beamong other things model parameters, measured valuesin an experiment, and statistical and systematic errorsassociated with measured values. Some of these variables,such as measured values are evidently observable, whereasothers, such as model parameters, can only be indirectlyinferred from variables whose values were observed.Before we move on, we stress the important distinc-tion between errors and uncertainties here: An error is aquantity contributing to the difference between the truthand the measurement. As the value of an error cannot bemeasured, the likelihood of potential values is typicallymodeled by a normal distribution and the standard de-viation of this distribution is referred to as uncertainty.In the pipeline, we model measurement errors explicitlybecause it introduces the possibility to check results forconsistency and enables a mathematical representation ofinformation to speed up computations.To discuss how the values of unobservable variables canbe inferred from variables whose values were observed, wepartition both the center vector and the covariance matrix into blocks corresponding to observable variables denotedby the index ’obs‘ and unobservable variables denoted bythe index ’hid‘ (for hidden): (cid:126)µ = (cid:18) (cid:126)µ obs (cid:126)µ hid (cid:19) , Σ = (cid:18) Σ obs Σ obs,hid Σ T obs,hid Σ hid (cid:19) (2)In the covariance matrix, the off-diagonal blocks Σ obs,hid encode linear relationships between observed and unob-served variables.Given that values of a subset of the variables have beenobserved, i.e., (cid:126)x obs = (cid:126)o , the pdf of the unobserved variablesis modified due to the propagation of information gainedabout the observed variables. The mathematical opera-tion to propagate this information is called conditioningwhere we go from the joint distribution π ( (cid:126)x obs , (cid:126)x hid ) tothe conditional distribution π ( (cid:126)x hid | (cid:126)x obs = (cid:126)o ) . For a MVNdistribution using the partitioned form of the center vec-tor and covariance matrix in eq. (2), the conditional pdfis again multivariate normal. The center vector and co-variance matrix are given by (cid:126)µ (cid:48) hid = (cid:126)µ hid + Σ T obs,hid Σ − obs ( (cid:126)o − (cid:126)µ obs ) (3) Σ (cid:48) hid = Σ hid − Σ T obs,hid Σ − obs Σ obs,hid . (4)These formulas (or mathematical equivalent ones) areknown as the Generalized Least Squares (GLS) methodin the field of nuclear data evaluation, e.g., [17, 18], imple-mented in various codes, e.g., [19–21]. Usually, unobservedvariables are model parameters or the values of cross sec-tions on a predefined energy mesh and observed variablescorrespond to measurements made in experiments. How-ever, we stress that any other type of variable can beadded as well. For instance, model parameters and sys-tematic errors can both be present at the same time inthe set of unobserved variables. We further elaborate onthis possibility in the next section, as it is used in theevaluation pipeline.In eqs. (3) and (4), the covariance matrix is assumedto be perfectly known, which may not always be the casein practice. As an example, if the calibration uncertaintyof a detector is not reported, only a reasonable range ofthis uncertainty may be available by taking into accountexperiences from similar experiments. In such cases, wecan set the values of missing uncertainties in a way thatis consistent to the information from other experiments,the model for fitting, and respecting reasonable limits.Regarding the consistency with the model for fitting,we should be aware that by our decision of employing acertain nuclear physics model in a conventional statisticalinference procedure, such as the generalized least squares(GLS) method, e.g., [17], we express our absolute confi-dence in the ability of the model to represent the truthgiven that values of model parameters are set correctly.Due to the complexity of nuclear processes, we usuallydo not believe that nuclear physics models can perfectly reproduce the truth. Therefore, simpler and very flexiblemathematical models, such as a piecewise linear function,4..Nuclear Data Evaluation Pipeline G. Schnabel et al. are sometimes used instead to avoid the distortion of re-sults due to an imperfect model. This path was takenin the evaluation of neutron data standards [2–4]. An al-ternative is to endow a nuclear physics model with moreflexibility by adding a simpler mathematical model ontop that captures the trends in the experimental data thenuclear physics model cannot reproduce. This approachis associated with the idea of model defects in the field ofnuclear data evaluation, e.g., [11, 22–24].One approach to fit missing or uncertain elements ofthe covariance matrix to enforce consistency between ex-perimental datasets and the model for fitting is marginallikelihood optimization (MLO), e.g., [25, sec. 5.4.1]. It re-lies on the marginal distribution of the observed values π ( (cid:126)x obs ) obtained from π ( (cid:126)x obs , (cid:126)x hid ) by integrating overthe unobservable variables (cid:126)x hid . Marginal distributions ofa MVN distribution are again multivariate normal. Weobtain them by removing the blocks from the center vec-tor and covariance matrix in eq. (2) corresponding to theunobserved variables, hence the marginal pdf is given by N ( (cid:126)x obs | (cid:126)µ obs , Σ obs ) .Marginal likelihood maximization seeks to find valuesfor the unknown uncertainties that maximize the proba-bility density for the observed values (cid:126)x obs = (cid:126)o . More tech-nically, missing or uncertain elements in Σ obs are adjustedso that N ( (cid:126)o | (cid:126)µ obs , Σ obs ) attains the largest possible value.Writing out the logarithmized pdf helps in understandingwhat the maximization actually means: N ( (cid:126)x obs | (cid:126)µ obs , Σ obs ) = − d ln(2 π ) − ln det Σ obs − ( (cid:126)x obs − (cid:126)µ obs ) T Σ − obs ( (cid:126)x obs − (cid:126)µ obs ) . (5)The second term on the right-hand side is proportional tothe differential entropy of the multivariate normal distri-bution. The third term is the generalized χ -value, whichtakes in comparison to the conventional χ -value the fullcovariance matrix into account. Maximizing the marginallikelihood therefore amounts to simultaneously minimiz-ing the information entropy and the generalized χ -value.Because the information entropy is a measure of modelcomplexity, higher entropy means a more complex model,the maximization process aims to increase the quality ofthe fit while avoiding too complex models.The availability of the partial derivatives of the loga-rithmized determinant and the inverse covariance matrixwith respect to uncertainties in analytic form, see appen-dices B 1 and B 2, lead to analytic expressions of partialderivatives of eq. (5). Quasi-newton methods, such as theL-BFGS algorithm [26], can leverage this information inthe optimization process. B. Representation of covariance matrices
The computation of the distribution parameters of theconditional pdf in eqs. (3) and (4) and the evaluation ofthe value of the marginal pdf in eq. (5) is problematic for a large number of observable variables due to the needto invert the associated covariance matrix Σ obs and thecalculation of its determinant. Making the pipeline appli-cable to a large number of observables was a design goalfrom the very beginning. In this section we outline therepresentation of covariance matrices to overcome compu-tational obstacles, which was already suggested by [27].For the following we assume a linear model of the form M lin ( (cid:126)p ) = (cid:126)σ ref + S mod ( (cid:126)p − (cid:126)p ref ) . (6)The term ‘model’ is intentionally left open in the currentdiscussion. It can be a linear approximation of a nuclearphysics model or a piecewise linear function. Both optionsmentioned are employed at different stages in the pipeline.Other possibilities include wavelets, Fourier series, smooth-ing splines and (linear approximations of) neutron trans-port models. Also combinations of such models are ofinterest, e.g., when both differential and integral data areconsidered in an evaluation. The treatment of non-linearmodels, such as implemented in the nuclear model codeTALYS is discussed in section II D.The relationship between observable variables (cid:126)σ exp andunobservable variables can be expressed as (cid:126)σ exp = M lin ( (cid:126)p ) + S sys (cid:126)ε sys + (cid:126)ε stat (7)where S sys denotes the mapping responsible to distributethe systematic errors (cid:126)ε sys to the affected observable vari-ables. The statistical errors (cid:126)ε stat are mutually independentand therefore directly added to the sum. In principle, theterm associated with systematic errors could be split intoseveral terms reflecting different types of errors, such asthose related to the detector calibration or to the samplethickness.The notational distinction between model parametersand systematic errors can be considered irrelevant becausethey can be dealt with in a unified way. Therefore wecombine the parameter vector and the systematic errors toone vector (cid:126)u and do the same for the associated mappingmatrices: S = (cid:0) S mod , S sys (cid:1) , (cid:126)u = (cid:18) (cid:126)p(cid:126)ε sys (cid:19) , (cid:126)u ref = (cid:18) (cid:126)p ref (cid:126) (cid:19) (8)Using these combined variables, eq. (7) can be rewrittenas (cid:126)σ exp = (cid:126)σ ref + S ( (cid:126)u − (cid:126)u ref ) + (cid:126)ε stat . (9)Going one step further, we denote the sum of the constants (cid:126)σ ref and − S (cid:126)u ref as (cid:126)σ and write (cid:126)σ exp = (cid:126)σ + S (cid:126)u + (cid:126)ε stat . (10)This equation serves as the basis to determine the co-variance matrix blocks in eq. (2). To that end, we ex-ploit the bilinearity of the covariance operator Cov [ (cid:126)x, (cid:126)y ] and therefrom induced properties of the variance opera-tor Var [ (cid:126)x ] := Cov [ (cid:126)x, (cid:126)x ] . The properties of these operators5..Nuclear Data Evaluation Pipeline G. Schnabel et al. are briefly reviewed in appendix A 1. Noteworthy, Var [ (cid:126)x ] yields the covariance matrix associated with (cid:126)x . Thus thecovariance matrix Σ obs associated with (cid:126)σ exp is given by Σ obs = Var [ (cid:126)σ exp ] = S Var [ (cid:126)u ] S T + Var [ (cid:126)ε stat ] . (11)Please note that this covariance matrix is different fromwhat is typical understood under an experimental covari-ance matrix. An experimental covariance matrix Σ exp reflects potential deviations of the experimental measure-ments from an assumed ‘truth’, i.e., N ( (cid:126)σ exp | (cid:126)µ true , Σ exp ) .In contrast to that, the covariance matrix Σ obs also in-corporates the covariance matrix of the unknown modelparameters, i.e., it describes fluctuations of potential mea-surements with respect to the prior expectation E [ (cid:126)σ exp ] = (cid:126)σ + S E [ (cid:126)u ] . Therefore this covariance matrix capturesthe fluctuations of hypothetical measurements obtainedby the following sampling procedure: Model parametersare sampled from the prior distribution, correspondingpredictions are calculated, and finally a perturbation isadded to the predictions, obtained by drawing a samplefrom the experimental covariance matrix Σ exp .The covariance matrix between observable and unob-servable variables is given by Σ obs,hid = Cov [ (cid:126)σ exp , (cid:126)u ] = S Var [ (cid:126)u ] . (12)Both Σ obs and Σ obs,hid are thus completely determined bythe mapping matrix S , the prior covariance matrix U =Var [ (cid:126)u ] associated with model parameters and systematicerrors, and the prior covariance matrix D = Var [ (cid:126)ε stat ] associated with statistical errors.The covariance matrix associated with observations istherefore amenable to the decomposition Σ obs = SUS T + D . (13)We exploit this representation of Σ obs in the computa-tion of the inverse and the determinant. To express thefollowing formulas concisely, we introduce the abbrevia-tion Z = U − + S T D − S . (14)Using the Woodbury matrix identity, see appendix B 4,the inverse appearing in eqs. (3) to (5) can be written as Σ − obs = D − − D − SZ − S T D − . (15)We usually do not need to compute Σ obs nor its inversebut rather products of these matrices with other vectors ormatrices. For instance, a matrix-vector product Σ − obs (cid:126)x isperformed sequentally from right to left starting with theevaluation of D − (cid:126)x yielding a vector. More informationon the efficient computation of matrix products involvinginverse covariance matrices is in appendix B 5.For the computation of the logarithmized determinantrequired in eq. (5), a similar relationship holds, ln det( Σ obs ) = ln det( D ) + ln det( U ) + ln det( Z ) . (16)A derivation of this identity is provided in appendix B 3.Working with the logarithm is essential for the compu-tation of determinants as they become easily so large orsmall that the representational capabilities of double pre-cision floats are exceeded.At first glance, the evaluation of the inverse and deter-minant of Σ obs according to the prescriptions in eqs. (15)and (16) seems to be more involved than the direct com-putation. The advantage of these formulas is owed to thestructure of the matrices. The covariance matrix D as-sociated with statistical errors is diagonal thanks to themutual independence of statistical errors, hence its inver-sion trivial. Its diagonal contains the squared statisticaluncertainties. The number of model parameters and sys-tematic uncertainties is typically much lower than thenumber of observable variables, hence the inversion of thesmaller dimensional matrix Z can be performed fasterthan of Σ obs in eq. (13). Moreover, model parameters arealmost always assumed to be independent of systematicerrors, and systematic errors are mutually independentof each other. In effect, this leads to U being sparse andblock diagonal. The matrix S is very rectangular, i.e., thenumber of rows equalling the number of observable vari-ables is usually much larger than the number of modelparameters in typical evaluation scenarios. The blocks in S corresponding to systematic experimental errors, suchas normalization errors, contain only a very small numberof non-zero elements, hence S exhibits a large degree ofsparseness, which can be exploited in matrix multiplica-tions and the calculation of inverses and determinants.More information on how the sparsity is exploited in com-putations can be found in appendix B 5. C. Gaussian processes
A MVN distribution is a model for the joint distributionof a finite number of variables. A Gaussian process isthe generalization of a MVN distribution to an infinitenumber of variables, i.e., functions. Gaussian processesappear in the modeling of systematic errors that varyas a function of incident energy, such as the calibrationerror of a detector, and they can be used to incorporatethe notion of model defects into nuclear data evaluations[11, 22–24, 28] and more generally in the calibration ofsurrogate models for computer experiments, e.g. [29].In the pipeline they are used to inject more flexibil-ity in energy-dependent model parameters of TALYS toaddress model defects. This approach was discussed andexplored in [9]. The available functionality in the nuc-dataBaynet package [16] goes beyond the assignment ofGaussian processes to model parameters. Any group ofelements present in the vector (cid:126)u in eq. (8) can be associ-ated with a Gaussian process. A good and comprehensiveintroduction to Gaussian processes with pointers to rele-vant literature can be found in [25]. In the following we6..Nuclear Data Evaluation Pipeline G. Schnabel et al. provide a brief discussion tailored to their application inthe pipeline.A Gaussian process is defined by a mean function m ( E ) ,which provides the apriori expectation of the functionvalues f ( E ) for all energies E in the admissible domain,and a covariance function κ ( E, E (cid:48) ) , which provides thecovariance between function values f ( E ) and f ( E (cid:48) ) forall admissible pairs of energies E and E (cid:48) : m ( E ) = E [ f ( E )] and (17) κ ( E, E (cid:48) ) = Cov [ f ( E ) , f ( E (cid:48) )] . (18)A defining property of a Gaussian process is that forany finite set of energies { E i } i =1 ..N , the function values { f ( E i ) } i =1 ..N are governed by a MVN distribution. There-fore, if we combine those function values to a vector (cid:126)f ,where f i = f ( E i ) , and construct the associated centervector (cid:126)m , where m i = m ( E i ) , and covariance matrix K ,where K ij = κ ( E i , E j ) , the pdf of the function values isgiven by N (cid:16) (cid:126)f | (cid:126)m, K (cid:17) . A typical choice for the covariancefunction is the so-called squared exponential, κ ( E, E (cid:48) ) = δ exp (cid:20) − λ ( E − E (cid:48) ) (cid:21) , (19)which is employed in the implementation of the pipelineto effect the evaluation of Fe. The amplitude δ repre-sents the apriori standard deviation of all function values f ( E ) and is therefore a measure of the range functionvalues are expected to cover. The length-scale λ deter-mines the correlation of function values at different en-ergies, hence is a measure of the expected similarity offunction values associated with similar energies. Otherchoices of the covariance function are possible. The thesisof Duvenaud [30, chap. 2] contains a pertinent catalogueof them and discusses their properties. The possibility ofenergy-dependent length-scales and amplitudes and theirdetermination on the basis of a collection of neighboringisotopes is explored in [31].To predict the function values at energies of interest { E pred i } i =1 ..M based on the function values observed atenergies { E obs i } i =1 ..N , we use eqs. (3) and (4) to get theconditional pdf. The elements of required vectors andmatrices are computed as follows: ( (cid:126)µ hid ) i = m ( E pred i ) , ( Σ hid ) ij = κ ( E pred i , E pred j ) (20) ( (cid:126)µ obs ) i = m ( E obs i ) , ( Σ obs ) ij = κ ( E obs i , E obs j ) (21) ( (cid:126)o ) i = f ( E obs i ) , ( Σ obs,hid ) ij = κ ( E obs i , E hid j ) (22)However, there are two issues with the use of Gaus-sian processes in our situation. First, the inference basedon Gaussian processes can become prohibitively compu-tationally expensive with an increasing number N of ob-served function values due to the needed inversion of Σ obs . The time required for an inversion scales with the cube ofthe number of observed function values. Second, typicalchoices of covariance functions define a Gaussian processthat corresponds to a fitting function with an infinitenumber of parameters. Thus, imposing a Gaussian pro-cess on any block of variables in (cid:126)u in eq. (8), e.g., energy-dependent nuclear model parameters that are propagatedto observations via the mapping matrix S , would be im-possible because an infinite number of energy dependentvariables needed to be stored in (cid:126)u .Both issues can be dealt with by introducing a finitedimensional approximation to Gaussian processes. Thereis a variety of sparse Gaussian process approximations,such as in [32], which was also explored for the determina-tion of model defects in nuclear data evaluation in [33]. Areview of various sparse GP approximations and unifiedframework for comparison is presented in [34].Instead of using any of those approximations, we de-cided to construct a low-dimensional approximation toa GP based on linear interpolation. The reason beingthat the sparse approximations referenced in the previousparagraph are associated with dense mapping matrices S whereas this matrix is very sparse if using linear interpo-lation. Consult, e.g., [33, sec. 2.2], for a brief discussion ofsparse GPs in terms of the associated mapping matrix. Inthe remainder of this section, we introduce the essentialformulas for a sparse GP approximation based on linearinterpolation.Given two mesh points associated with energies E i and E i +1 , where E i < E i +1 and function values y i and y i +1 ,respectively, values at intermediate energies can be deter-mined by linear interpolation: g ( E ) = (cid:18) E i +1 − EE i +1 − E i (cid:19) y i + (cid:18) E − E i E i +1 − E i (cid:19) y i +1 if E i ≤ E < E i +1 (23)To state the formulas in a general way for a completemesh of energies, we introduce the abbreviation c i ( E ) = (cid:40) E i +1 − EE i +1 − E i if E i ≤ E < E i +1 otherwise (24)and d i ( E ) = (cid:40) E − E i − E i − E i − if E i − ≤ E < E i otherwise (25)as well as their sum, f i ( E ) = c i ( E ) + d i ( E ) . (26)Piecewise linear interpolation, i.e., locating for an energyof interest E the enclosing energies on the mesh and thenperforming linear interpolation according to eq. (23), cannow be concisely written as g ( E ) = M (cid:88) i =1 f i ( E ) y i . (27)7..Nuclear Data Evaluation Pipeline G. Schnabel et al. Thanks to the bilinearity of the covariance operator, thecovariance between function values at arbitrary energiescan be computed by
Cov [ g ( E ) , g ( E (cid:48) )] = M (cid:88) i =1 M (cid:88) j =1 f i ( E ) f j ( E (cid:48) ) Cov [ y i , y j ] . (28)This formula shows that the covariance matrix ˜ K with ˜ K ij = Cov [ y i , y j ] associated with a finite number of vari-ables M , in combination with linear interpolation enablesthe computation of covariances for arbitrary pairs of en-ergies E, E (cid:48) . In other words, it is a valid specification ofa covariance function κ ( E, E (cid:48) ) and together with a meanfunction m ( E ) completely characterizes a Gaussian pro-cess.Furthermore, we can identify a mapping matrix S whichallows us to incorporate Gaussian processes in a princi-pled way into the representation presented in eq. (13).For example, let us assume that we have an energy-dependent systematic error assigned to a certain experi-mental dataset that contains measured values at energies { E obs i } i =1 ..N . Now using a reduced energy mesh contain-ing the energies { E hid j } j =1 ..M as a basis for linear interpo-lation, the partial derivatives of eq. (27) with respect tothe y i , i.e., the values at the knot points, constitute theelements of the mapping matrix: S ij = ∂g ( E obs i ) ∂y j = f j ( E obs i ) . (29)Noteworthy, for a given energy E obs i , there are only twonon-zero values of f j ( E obs i ) , hence the mapping matrix isvery sparse with only two non-zero elements in each row.For the sake of an example, if Σ obs were only made upof an energy-dependent systematic uncertainty given as aGP and uncorrelated statistical uncertainties reflected in D , we would have Σ obs = S ˜ KS T + D (30) Σ obs,hid = S ˜ K (31) Σ hid = ˜ K (32)with the elements of ˜ K defined below eq. (28) and thoseof S in eq. (29). D. Levenberg-Marquardt algorithm with prior
We assumed a linear link between model parametersand observables in the discussion of section II B. How-ever, the model parameters of a nuclear physics modelare typically linked to observables in a non-linear way. Insome evaluations, the model parameters are cross sectionsand adjusted according to measurements of ratios of thesecross sections. This would also result in a non-linear link between parameters and observables. Irrespective of thecause of the non-linearity, the Jacobian matrix S mod be-comes a function of the parameter vector, S mod ( (cid:126)p ref ) . Forsimplicity, we refer from now on to a nuclear model but itshould be understood that the discussion is valid for anynon-linear link between parameters and observations.One possible solution is to apply the GLS formulasiteratively and to use an updated linear approximationof the nuclear model in each iteration. For instance, thisapproach is one option offered by the CONRAD evalua-tion code [35, 36] and SAMMY [37] and briefly discussedin [35] and in more detail in [37, sec. IV.A.3]. The modifiedLevenberg-Marquardt algorithm outlined in this sectionimproves upon this approach by the introduction of anadaptive damping term.In the following we abbreviate J = S mod . A linearapproximation of a nuclear model as in eq. (6), which wereiterate here for convenience, M lin ( (cid:126)p ) = (cid:126)σ ref + J ( (cid:126)p ref ) ( (cid:126)p − (cid:126)p ref ) , (33)is only reliable in vicinity of the expansion point (cid:126)p ref . Thechoice to pick the apriori best guess (cid:126)p as reference point (cid:126)p ref is plausible but an updated best guess (cid:126)p obtained bythe conditioning formula in eq. (3) can be misleading dueto the fact that the linear map J ( (cid:126)p ) is not representativeanymore for the region of the parameter space where (cid:126)p resides.The Levenberg-Marquardt algorithm [38, 39] extendsthe idea of the linear least squares method to a non-linear link between model parameters and predictions. Weuse a modified version of the Levenberg-Marquardt algo-rithm [11] to take into account an experimental covariancematrix and a prior parameter covariance matrix.To prepare the discussion of the modified Levenberg-Marquardt algorithm, we write the pdf of (cid:126)p conditionedon the observable variables (cid:126)σ exp in the factorized form π ( (cid:126)p | (cid:126)σ exp ) = 1 π ( (cid:126)σ exp ) (cid:96) ( (cid:126)σ exp | (cid:126)p ) π ( (cid:126)p ) . (34)The first factor /π ( (cid:126)σ exp ) is independent of (cid:126)p and hencerepresents a normalization constants to ensure the inte-gral of the conditional pdf being unity. For a non-linearmodel link to map from (cid:126)p to the observables in (cid:126)σ exp , thisnormalization constant can usually not be obtained an-alytically and one either needs (1) to take recourse to aMonte Carlo approach to estimate it or (2) approximatethe conditional pdf by a simpler functional form for whichit can be calculated analytically. Later in the discussion,we take the second approach, but for the moment it is notnecessary to have the conditional pdf π ( (cid:126)p | (cid:126)σ exp ) correctlynormalized.We continue by specifying for both pdfs, i.e., the so-called likelihood (cid:96) ( (cid:126)σ exp | (cid:126)p ) and prior π ( (cid:126)p ) , multivariatenormal distributions, (cid:96) ( (cid:126)σ exp | (cid:126)p ) = N ( (cid:126)σ exp | M ( (cid:126)p ); Σ exp ) , (35) π ( (cid:126)p ) = N ( (cid:126)p | (cid:126)p ; P ) . (36)8..Nuclear Data Evaluation Pipeline G. Schnabel et al. If the model were linear, i.e., could be written in theform of eq. (33), the resulting conditional pdf π ( (cid:126)p | (cid:126)σ exp ) in eq. (34) would be multivariate normal. Its maximumwould then be associated with the updated vector (cid:126)p ob-tained by the formula for the conditional mean vector ineq. (3) if making the identifications: (cid:126)µ hid = (cid:126)p , Σ hid = P (37) (cid:126)µ obs = M lin ( (cid:126)p ) , Σ obs = JP J T + Σ exp (38) (cid:126)o = (cid:126)σ exp , Σ obs,hid = JP (39)However, the nuclear model represented by the func-tion M ( (cid:126)p ) is non-linear and does in general not have ananalytical form because the nuclear model code may usenumerical methods to solve differential equations or im-plements a stochastic simulation. Consequently, there isno simple update formula to locate the parameter vector (cid:126)p associated with the maximum of the conditional pdf π ( (cid:126)p | (cid:126)σ exp ) .The idea of the modified Levenberg-Marquardt (LM)algorithm is to still rely on the formula for the conditionalmean eq. (3) and exploit a linear approximation of thenuclear model to locate a parameter vector closer to themaximum. Introducing the quantities, A = J T Σ − exp J + P − + λ I , (40) (cid:126)b = J T Σ − exp ( (cid:126)σ exp − (cid:126)σ ref ) + P − ( (cid:126)p − (cid:126)p ref ) , (41)the conditional mean vector is given by (if λ = 0 ): (cid:126)p prop = (cid:126)p ref + A − (cid:126)b . (42)Consult appendix C 1 for a derivation of this update for-mula. We discuss the purpose of the term λ I in a moment.It is not evident by visual inspection that this formula(with λ = 0 ) is equivalent to eq. (3) with the appropri-ate substitutions of eq. (39) but they can be transformedinto each other by making use of the Woodbury matrixidentity [40]. The Woodbury identity is briefly discussedin appendix B 4.The LM algorithm is an iterative algorithm. In eachiteration it relies on a linear approximation of the non-linear model constructed at the obtained proposal vector (cid:126)p prop of the previous iteration denoted by (cid:126)p prev . Thereforethe linear approximation employed in the current iterationis characterized by (cid:126)p ref = (cid:126)p prev , (cid:126)σ ref = M ( (cid:126)p prev ) and J = J ( (cid:126)p prev ) .In the case of a mildly non-linear model, the proposal(with λ = 0 ) in eq. (42) is close to the true maximum ofthe pdf. For stronger degrees of non-linearity, the currentlinear approximation can only be trusted in close vicin-ity of the current reference parameter vector (cid:126)p ref and theproposed parameter vector (cid:126)p prop may lie outside this re-gion of trust. This would also invalidate the trust in thepertinence of (cid:126)p prop .The innovation of the LM algorithm is the introductionof the additional term λ I with λ > and I a diagonalmatrix. If the proposal parameter vector is too far away from the reference parameter vector to rely on the linearapproximation, the value of λ is increased. Consequentlythe matrix A becomes more diagonal and so its inverse.As the vector (cid:126)b is the gradient of ln π ( (cid:126)p | (cid:126)σ exp ) evaluated at (cid:126)p ref , i.e., b i = ∂ (ln π ( (cid:126)p | (cid:126)σ exp )) /∂p i , the parameter updatedegenerates gradually to the gradient ascent method withan increasing value of λ and at the same time the stepsize becomes smaller. The adjustment of λ is thereforea mechanism to smoothly transition between an updateusing the formula for the conditional mean of the MVNdistribution if the nuclear model is only mildly non-linearand the gradient ascent method for larger degrees of non-linearity.The adjustment of λ is effected based on the com-parison between the expected increase of the value of ln π ( (cid:126)p | (cid:126)σ exp ) according to the linear approximation andthe actual increase by evaluating the nuclear physicsmodel. Let f ex ( (cid:126)p ) = ln π ( (cid:126)p | (cid:126)σ exp ) be the value based onthe exact nuclear physics model (applied in eq. (35)) and f lin ( (cid:126)p ) = ln ˜ π ( (cid:126)p | (cid:126)σ exp ) the value based on the linear ap-proximation constructed at (cid:126)p ref . The criterion for adjust-ment is the gain, see, e.g., [41], ρ = f ex ( (cid:126)p prop ) − f ex ( (cid:126)p ref ) f lin ( (cid:126)p prop ) − f lin ( (cid:126)p ref ) . (43)The strategy proposed by Marquardt [39] is to select aftereach iteration a new value of lambda according to theprescription λ (cid:48) = λ if ρ < . λ if . ≤ ρ < . λ/ if . ≤ ρ . (44)Noteworthy, if a proposed parameter vector (cid:126)p prop is asso-ciated with a lower value of the posterior pdf in eq. (34),it is rejected and a new proposal computed on the basisof the updated value of λ .In the case of a large number of model parameters,it is computationally attractive to optimize only the pa-rameters which are sensitive to the experimental data.Afterwards, due to prior correlations between sensitiveand insensitive parameters, also insensitive parametersmust be updated. This update procedure is described inappendix C 5. E. Taylor expansion of logarithmized pdf
Once the vector (cid:126)p associated with the maximum ofthe posterior distribution in eq. (13) is located using theLM algorithm, one may be tempted to apply the condi-tioning formulas eqs. (3) and (4) with the substitutionsin eq. (39) for one further improvement of the posteriorcenter vector (cid:126)p and to obtain the posterior covariance ma-trix P . This strategy only works if the statistical modelis consistent with the experimental data, i.e., the remain-ing differences between the posterior model prediction9..Nuclear Data Evaluation Pipeline G. Schnabel et al. and the experimental data can be explained by the sta-tistical and systematic uncertainties associated with theexperiments. Otherwise, in the case of the experimentaldata and/or the model being deficient, the results are mis-leading. For instance, the update of (cid:126)p could undergo acomparatively large adjustment at odds with the foundsolution by the LM algorithm, and the associated valueof π ( (cid:126)p | (cid:126)σ exp ) would be drastically reduced. The reason forthis behavior is explained in appendix C 3.If experimental data or the nuclear model have deficien-cies, a better solution to obtain the posterior covariancematrix P is to perform a second-order Taylor expansionat the posterior maximum. Misspecified experimental dataor models remain still problematic but the covariance ma-trix obtained in this way is at least consistent with thevector (cid:126)p found by the LM algorithm. The form of asecond-order Taylor approximation of the logarithmizedposterior pdf is given by ln π ( (cid:126)p | (cid:126)σ exp ) = π + G ( (cid:126)p − (cid:126)p ref )+ 12 ( (cid:126)p − (cid:126)p ref ) T H ( (cid:126)p − (cid:126)p ref ) . (45)Noteworthy, the gradient G is not associated with thenuclear model but with the posterior pdf, i.e., G i = ∂ (ln π ( (cid:126)p | (cid:126)σ exp )) /∂p i . Ideally, if the LM algorithm has con-verged and (cid:126)p ref = (cid:126)p , i.e., the expansion is done exactlyat the parameter vector corresponding to the posteriormaximum, the gradient vanishes. We want to extract theposterior covariance matrix P from this expression. Com-paring the expansion in eq. (45) with the logarithmizedpdf of a multivariate normal distribution in eq. (5), wecan identify the posterior covariance matrix P as P = − H − . (46)Therefore in order to determine the posterior covariancematrix, one needs to numerically determine the Hessianmatrix. Its elements are given by the second-order deriva-tives with respect to the model parameters, H ij = ∂ (ln π ( (cid:126)p | (cid:126)σ exp )) ∂p i ∂p j (cid:12)(cid:12)(cid:12)(cid:12) (cid:126)p = (cid:126)p ref . (47)For the numerical evaluation of the full Hessian matrix,the number of required nuclear model invocation scalesquadratically with the number of adjustable model param-eters. As the aim of the pipeline is large scale evaluationstaking into account potentially hundreds of model param-eters, the full determination of the Hessian matrix maybe impractical or even infeasible in the case of a compu-tationally expensive nuclear model.Therefore we implemented an approximation of the fullHessian matrix, which exploits the mathematical struc-ture of the pdf in eq. (34) when prior and likelihood aremultivariate normal distributions. In this case, the fullHessian matrix can be written as (see appendix C 4) H = − ( U + J T QJ + R ) , (48) with J being the Jacobian matrix of the nuclear modelappearing in eq. (33), Q being the inverse experimentalcovariance matrix and R being the inverse prior parame-ter covariance matrix. The elements of the matrix U aredetermined by a product involving second-order deriva-tives of the nuclear model, elements of the experimentalcovariance matrix and differences between model predic-tions and experimental data. The approximation consistsin replacing the exact matrix U by an approximation thatneglects second-order derivatives involving two differentmodel parameters. This approximation is explained inmore detail in appendix C 4. III. IT BUILDING BLOCKS
In our opinion, the implementation of an evaluation asa sequence of scripts in a pipeline should be formulatedat a high-level of abstraction. Scripts should be conciseso that important decisions taken in the evaluation canbe easily understood. For instance, the pipeline shouldnot contain the code for fitting parameters with all thegory details of the Generalized Least Squares (GLS) orsimilar fitting method but instead invoke a function withan emblematic name to effect such fitting. This idea doesnot mean that the code implementing lower-level function-ality, such as the GLS approach, should be inaccessible. Itjust means that such code should not clutter the scripts inthe pipeline to draw away the attention from importantaspects, such as the selection of experimental datasets forthe evaluation.Therefore this section is dedicated to—in lack of a bet-ter name—IT building blocks. We discuss functionalitythat should in our opinion be available out-of-the-box toevaluators who wish to design an evaluation pipeline. Weoutline our approach to the compartmentalization of suchfunctionality and the design considerations that guidedus. Specifically, we elaborate on the following:•
Handling of information in the Bayesian con-text:
In the previous section we outlined the math-ematical machinery relied upon in our prototypeof an evaluation pipeline. We believe that thinkingabout an evaluation in terms of, e.g., experimentalmeasurements vectors (cid:126)σ exp , experimental covariancematrices Σ exp , and mapping matrices S may be nottoo intuitive. Instead, the construction of such math-ematical objects should be automatically derivedfrom a more intuitive representation of informationin terms of specific types of systematic componentsand their association with experiments.• Handling of experimental data:
The combina-tion of information of various experiments is at theheart of nuclear data evaluation. Therefore it isan essential requirement to have flexible means tosearch and retrieve relevant experimental data, inparticular data in the EXFOR library. Besides anincreased convenience in manual evaluation work, a10..Nuclear Data Evaluation Pipeline G. Schnabel et al. better handling of experimental data facilitates therenormalization of experimental data according tothe latest evaluations of monitor reactions and opensnew possibilities for uncertainty quantification (UQ)and quality assurance based on advanced statisticalmethods and machine learning. In the course of anevaluation, it may also be helpful to augment exper-imental data with additional information derivedfrom the raw data, e.g., a covariance matrix derivedfrom the individual systematic uncertainty compo-nents. Adding such information should be facilitatedas much as possible from the perspective of an evalu-ator. Also facilitating linking together relevant data,such as evaluations of standard reactions with ex-perimental data where they have been employed asmonitor reactions, can have a positive impact on thequality and transparency of evaluations.•
Mapping of model predictions to EXFOR:
Not all measurements documented in the EXFOR li-brary can be directly related to the output of nuclearmodel codes, such as TALYS employed in our pro-totype of the pipeline. As an example, many exper-iments measure ratios of nuclear reactions whereasnuclear model codes provide the individual reactionsas output. The translation of results from nuclearmodel codes to predictions comparable with experi-mental data is an essential requirement of any eval-uation and should therefore be easy to achieve byan evaluator.•
Parallelization of nuclear model invocations:
With the availability of ever increasing computingpower reflected in growing computational demandsof evaluation approaches, cluster computing is notany longer only reserved to large research projectsbut can be a viable option even for a single evaluator.Ideally, evaluators should not be required to careabout how computations are run in parallel butinstead empowered to focus completely on what theywant to compute.The following subsections provide more details on de-sign considerations to implement the functionality of eachIT building block. These building blocks have been imple-mented as
R packages , i.e., modules that can be used inthe programming language R [42]. The restriction to oneprogramming language is certainly a disadvantage butother languages, such as Python, provide equivalent func-tions and data types, which makes the discussion of designconsiderations still worthwhile. A valuable improvementin the future would be an implementation of the IT build-ing blocks so that they can be conveniently employed inany popular programming language.
A. Handling of information in the Bayesian context
In section II B we explained that it is computationallyadvantageous to avoid the explicit construction of covari- ance matrices and instead maintain a representation ofthe form Σ obs = SUS T + D (49)with the covariance matrix U associated with systematicuncertainties and model parameter uncertainties and thediagonal covariance matrix D associated with statisticaluncertainties.Besides the aspect of computational efficiency, there isanother advantage of this representation in terms of in-formation management. It is always possible to constructthe full covariance matrix on the basis of uncertainty com-ponents but without further knowledge impossible to de-compose Σ obs back into its individual components.As explained in section II B, the statistical model under-lying this representation of the covariance matrix couplingthe vector of experimental measurements (cid:126)σ exp with thevector (cid:126)u of systematic error components and potentiallymodel parameters is given by (cid:126)σ exp = (cid:126)σ + S (cid:126)u + (cid:126)ε stat . (50)The matrix S reflects how systematic error componentsand model parameters—summarized under the term sys-tematic components —are summed up and assigned to ex-perimental data points.The R package named nucdataBaynet [16]—an abbrevi-ation of nuclear data in a Bayesian network —provides thefunctions to do Bayesian inference exploiting the decom-position in eq. (50) for efficient computation. To discussthe data structures employed by this package, we consideran experiment with three measured data points being af-fected by the the following systematic error components:(1) a detector calibration error, (2) sample thickness error,and (3) background noise. There can be other systematicerrors but it suffices to consider these three in our dis-cussion. These systematic errors can usually be assumedindependent of each other. For instance, the measurementdevice to determine the thickness of the specimen is notcoupled by any physical process to the measurement proce-dure to determine the calibration of the detector. The dia-gram in fig. 1 depicts the situation as a so-called Bayesiannetwork , e.g. [12]. Each node is associated with a value.An arrow from a node A to another node B indicates thata certain part of the value associated with B is given bya linear transformation of the value of A . In the graph,systematic errors contribute to the distortion of all threedata points whereas each statistical error is connected toa single measured data point. As an important remark,systematic components do not have to be connected to alldata points, e.g., when considering data points stemmingfrom several independent experiments.The depicted structure of the relation between mea-sured data points, systematic errors and statistical errorsis very typical and we second the recommendation in [27]that it may be regarded as the default structure for therepresentation of information of a nuclear experiment.For this reason, we decided to base the nucdataBaynet package implementing the statistical functionality out-11..Nuclear Data Evaluation Pipeline G. Schnabel et al. detector calibration errorsample thickness errorbackground error data point 1data point 2data point 3 statistical error 1statistical error 2statistical error 3Figure 1. Bayesian network reflecting the typical relation between measured data points and associated systematic and statisticalerrors. Systematic errors are assumed to be mutually independent but each systematic error can possibly affect all data datapoints. Statistical errors are also mutually independent and each of them is associated with exactly one data point. lined in the previous section on this specific representation.In the remainder of this section, we sketch how this repre-sentation has been translated to data structures and howthese data structures are used to establish the mappingof systematic components to measured data points.
1. Experimental data and statistical uncertainties
The experimental data points together with their statis-tical uncertainties are stored in a tabular data structureknown as data frame in Python and R. For their handling,we use the R package data.table [43] which offers supportfor efficient searching, sorting and grouping of rows in dataframes. This data structure with some example contentis depicted in table I.
IDX EXPID DIDX REAC L1 DATA UNC380 13764002 573 (N,TOT) 2.00 3255.80 127.05381 13764002 574 (N,TOT) 1.99 3659.20 116.40382 13764002 575 (N,TOT) 1.99 3246.60 128.05383 13764002 576 (N,TOT) 1.99 3408.70 123.96384 13764002 577 (N,TOT) 1.98 3359.80 124.44385 13764002 578 (N,TOT) 1.98 2449.20 162.79Table I. Structure of the data frame containing experimentaldata as used by package nucdataBaynet . The essential columns are
IDX , DATA and
UNC . Thecolumn
IDX indicates the position in (cid:126)σ exp , the colum
DATA , the measured value at that position and
UNC theassociated statistical uncertainty. Depending on the typeof systematic component, e.g., systematic error or modelparameter, some of the other columns may be relevantto map the systematic components to the measured datapoints.
2. Systematic components
Systematic components are also stored in a data frame.An example is provided in table II. These components are enclosed in a blue dotted rectangle in fig. 1. The es-
IDX EXPID DIDX ERRTYPE DATA UNC1 REACEXP-N,INL-01 1 pw 0.00 100.002 REACEXP-N,INL-01 2 pw 0.00 100.003 REACEXP-N,INL-01 3 pw 0.00 2.004 EXPID-13764002 1 sys-rel 0.00 0.105 EXPID-22316003 1 sys-abs 0.00 100.006 EXPID-23134005 1 sys-rel 0.00 0.107 EXPID-32201002 1 sys-rel 0.00 0.068 EXPID-40532014 1 sys-rel 0.00 0.10Table II. Structure of the data frame containing systematiccomponents as used by package nucdataBaynet . Column
GP-TYPE has been omitted for better display. sential columns are
IDX , ERRTYPE , DATA and
UNC .The column
IDX indicates the position of the systematiccomponent in (cid:126)u , the column
DATA stores its value, and
UNC the associated uncertainty. Typically, systematic er-rors are a priori expected to be absent, hence the zerosin the
DATA column. This is the default assumption inthe GLS method as employed for nuclear data evaluation,even though systematic errors are usually not made ex-plicit there. The non-zero uncertainty assigned in column
UNC allow for non-zero posterior expectations.The important assumption taken is that all systematiccomponents are mutually independent, which allows tostore the uncertainties in a column of this data frameinstead of managing a full covariance matrix. This inde-pendence assumption must sometimes be abandoned, e.g.,when dealing with energy dependent systematic uncer-tainties. We are going to discuss this case later. For themoment we assume the independence assumption to hold.
3. Mapping of systematic components to experimental data
Having introduced the two data structures for experi-mental data and systematic components, we can discussthe mapping of systematic components to experimentaldata points, illustrated by the arrows in fig. 1. There can12..Nuclear Data Evaluation Pipeline G. Schnabel et al. be many types of systematic components, such as relativesystematic errors, absolute systematic errors, and parame-ters characterizing the shape of a piecewise linear function.The type of systematic component is indicated in the col-umn
ERRTYPE . In the example in table II, the identifier pw denotes a piecewise linear function and the correspond-ing elements in the DATA and
UNC column contain theprior expectation and uncertainty, respectively, of the val-ues at the mesh points. Analogously, the identifiers sys-rel and sys-abs refer to relative and absolute systematicerrors, respectively. It is crucial to be able to map system-atic components to measured data points. Consideringeqs. (49) and (50), this means to be able to construct thematrix S effecting the mapping.During the conception of the pipeline, it was evidentfrom the beginning that we cannot implement the map-pings for all conceivable types of systematic componentsand therefore our objective was to design an extensiblearchitecture to make the package as future-proof as possi-ble. The final implementation corresponds directly to thefollowing picture: Each type of systematic component isassociated with a handler which knows how to deal withthis particular type. A manager module keeps track of theavailable handlers. A request to the manager to constructthe mapping matrix S for a given pair of a data framewith experimental data expDt and another one with sys-tematic components sysDt , the manager splits the dataframe sysDt in chunks where each cunk corresponds to aspecific type. For each chunk, the manager asks in turnthe available handlers if they know how to map this typeof systematic component. The first handler answering inthe affirmative is tasked with the construction of the map-ping matrix for that chunk. At the very end, the managercombines all the results of the different handlers to a finalresult, which is the matrix S to map all systematic com-ponents at once to the experimental data. In the diagramshowing the Bayesian network, the matrix S defines theconnections between the systematic components and themeasured data points.In the language of informatics, each handler is a moduleproviding the two functions:• getErrTypes() : informs about the types of system-atic components the handler module is able to map.• map(expDt, sysDt) : returns the matrix implement-ing the mapping of systematic components in sysDt to the experimental data in expDt . The data frame expDt contains all experimental data whereas sysDt is a chunk of the full data frame containing only the ERRTYPE the handler knows how to map.These functions are the required interface to interact withthe manager. Handlers can and usually do possess morefunctions for their setup. To make available a new typeof systematic component, one has to implement a handlercontaining the two functions mentioned above.Beside the essential columns
IDX , DATA and
UNC indata frame expDt and
IDX , ERRTYPE , DATA and
UNC in data frame sysDt , individual handlers typically rely on some of the other columns to establish the mapping.Columns not needed by a particular handler are ignoredby that handler. This flexibility is very helpful as newcolumns can be introduced to implement the functionalityof a new handler without interfering with the functioningof already existing ones.To make the functioning of a handler more comprehen-sible, we outline the handler implementing the mappingof normalization errors. This handler receives from themanager the subset of rows in data frame sysDt where thecolumn
ERRTYPE contains the string sys-rel . We takeas an example the row in table II where the column
IDX is 4. The string
EXPID-13764002 in the
EXPID columnof sysDt informs the handler that this systematic errorcomponent is associated with rows in expDt that containthe string in the column
EXPID . The handlerthen identifies these rows in expDt and retrieves the in-dices in
IDX . Taking table I as example, the first indexthere in combination with the index in sysDt leadsto the specification S , = 3255 . where the value isthe cross section to have the mapping of a relative system-atic error component. The other elements of S associatedwith other rows of expDt are assigned analogously. Thisexample also shows that the resulting matrix S is sparse.The mapping of a normalization error requires only onenon-zero element per row. This sparsity of the mappingmatrix is very typical for systematic error componentsand exploited in the nucdataBaynet package by relyingon the support of sparse matrices provided by the Matrix package [44].Noteworthy, the functionality of the handler dealingwith normalization errors is more general. For instance, anormalization error can be assigned to experimental datapoints based on the string in the
REAC column. Therebyeach experimental data point of a specific reaction channelis associated with the same normalization error. Such anormalization error can be introduced in addition to nor-malization errors at the level of datasets. The next sectionon the handling of experimental data elaborates on theflexibility in data retrieval from the EXFOR library. Forinstance, it is straight-forward to add to the data frame expDt a column with the institution where the experi-ment was performed or a column indicating the detectortype. Normalization errors could also be associated withdata points based on those attributes, which enables thecreation of elaborated error structures while maintainingcomputational tractability thanks to the decompositionof the covariance matrix explained before.
4. Correlated systematic components
The assumption of mutual independence of the system-atic components does not always hold. Some systematicerrors of the experiment may be fully or partially corre-lated over energy. Also Gaussian processes imposed onenergy-dependent TALYS parameters, i.e., functions ofenergy, represent a prior that correlates the parameter13..Nuclear Data Evaluation Pipeline G. Schnabel et al. values at nearby energies. From a technical point of view,there is no difference between a Gaussian process on anenergy-dependent parameter and a systematic experimen-tal error correlated over energy. Both cases are addressedby the introduction of off-diagonal elements in the co-variance matrix U . As the data structure to representsystematic components outlined in section III A 2 only al-lows the storage of uncertainties, i.e., assuming U to bediagonal, another data structure is needed to account forthe correlations defined by a Gaussian process.To discuss this case, we show in table III an example of adata frame for systematic components which also containsa Gaussian process specification for energy-dependent sys-tematic components. EXPID ERRTYPE GPTYPE DATA UNC ENEXPID-32201002 sys-rel NA 0.00 0.06 NAEXPID-40532014 sys-rel NA 0.00 0.10 NATALYS-d1adjust n endep sqrexp 1.00 0.10 1.00TALYS-d1adjust n endep sqrexp 1.00 0.10 2.00Table III. Example of a data frame for systematic components sysDt with a Gaussian process specification. Some columnsnot relevant in the current in the current discussion are omit-ted in this table for better display. The
ERRTYPE named talyspar_endep is abbreviated. NA stands for not available and is a special value in R.
The assignment of the GP is indicated by the string sqr-exp in the column
GPTYPE . Here the GP is imposed onthe energy-dependent TALYS parameter associated withthe input keyword d1adjust n but GPs can also be as-signed to systematic experimental errors. As an aside, weimplemented a specific handler for the mapping (or in thiscontext more typical diction: propagation) of model pa-rameters to measured data points. The handler is responsi-ble for rows with
ERRTYPE being either talyspar_endep or talyspar .The important column for correlated systematic com-ponents is GPTYPE . This column indicates the type ofGaussian process imposed on the systematic componentsin order to introduce correlations. As for the mapping ofsystematic components to measured data points, special-ized GP handlers are responsible for the different Gaussianprocess types. The implementation of a GP handler fora specific
GPTYPE can rely on any column in the sysDt data frame to construct the parts of the covariance matrix U associated with this GP type. The GP handler for the GPTYPE identified by sqrexp , in particular, requires thecolumn EN , which indicates the energy. This informationis used in the construction of the covariance matrix onthe basis of the covariance function defined in eq. (19).A GP handler needs to implement the following twofunctions:• getGPTypes() : returns a list of strings of GP typessupported by the handler.• cov(sysDt, gpDt) : returns the covariance elementsof systematic components sysDt of systematic com- ponents associated with a Gaussian process typesupported by the handler.Besides these two functions, GP handlers usually haveadditional functions for their setup. A manager moduleaware of the available handlers splits the data frame of sys-tematic components in chunks according to the GPTYPE and delegates the individual chunks to the appropriateGP handlers to compute the parts of the covariance ma-trix. This manager module is the same that handles themapping from systematic components to measured datapoints discussed in section III A 3. In a request to the man-ager module to return the complete covariance matrix U ,the manager merges the covariance elements returned bythe GP handlers with the diagonal elements correspond-ing to the uncertainties in the UNC column in the dataframe sysDt , see table III. In the current implementation,values in the
UNC column in sysDt are ignored if a GPis assigned to the systematic component.Covariance functions usually depend on hyperparame-ters, such as the amplitude δ and the length scale λ inthe case of the squared exponential covariance function.These specifications are stored in a data frame, which wecall in the following gpDt . An example of such a dataframe is given in table IV. This data frame is required by IDX EXPID GPTYPE PARNAME PARVAL1 TALYS-d1adjust n sqrexp sigma 0.102 TALYS-d1adjust n sqrexp len 2.003 TALYS-d1adjust n sqrexp nugget 0.00Table IV. Example of a data frame with the specification ofhyperparameters of the Gaussian processes. the GP handlers to compute the covariance elements in U . For this reason, it must also be passed as argument tothe cov function listed above.As a final remark regarding the nucdataBaynet pack-age, the structure imposed on the experimental data bystoring them in the outlined data structures is exploitedin the computation of the inverses and determinants ap-pearing the GLS formulas eq. (4) and marginal likelihoodformula eq. (5) by relying on the matrix identities ineqs. (15) and (16). The same identities are also employedin the analytical computation of the derivatives of thelogarithmized marginal likelihood with respect to hyper-parameters in Gaussian processes or any other type ofuncertainty associated with a systematic component. Thederivatives of inverse matrices and logarithmized deter-minants appearing in this representation of the logarith-mized marginal likelihood can be evaluated analytically,see appendices B 1 and B 2. These derivatives are takeninto account in the adjustment of hyperparameters anduncertainties by locating the maximum of the marginallikelihood using the L-BFGS algorithm [26].14..Nuclear Data Evaluation Pipeline G. Schnabel et al. B. Retrieval and handling of experimental data
The management of experimental data is an essentialrequirement for the creation of evaluation pipelines witha focus on automation and reproducibility. Both the re-trieval of data using flexible capabilities of search and thehandling of retrieved data are important.The
EXFOR library [6] maintained and continuouslyextended by the International Network of Nuclear Reac-tion Data Centres (NRDC) contains a large collection ofnuclear cross section data and related quantities. It re-lies on the EXFOR format [45] to store and organize theinformation associated with experiments.The
EXFOR format is a text-format, hence readableby humans and hierarchically structured with the top twolevels being so-called entries and subentries . It exhibitsa rigid structure and at times complex syntactical rules,such as those for the designation of the reaction system.A detailed description of this format is available in [45].Unfortunately, it is sometimes perceived as impracticalto extract relevant data. A possible reason may be thatpopular programming languages, such as Java, Python,Matlab and R, do not natively support the extraction ofinformation stored in the EXFOR format.Several development activities were concerned with eas-ier access to the EXFOR library and data in the EXFORformat. The EXFOR web retrieval system [5] enables theuser to search data based on a comprehensive catalogueof criteria and retrieve the result in a variety of formatsincluding JSON and XML, which are better suited fornumerical processing than the original EXFOR format.Quick plotting in the web browser is also supported. An-other initiative to facilitate working with the data in EX-FOR was the development of the
X4toC4 code [46], whichallows the translation of the EXFOR format to a tabu-lar format named C4 . The x4i code by D.A. Brown anda pure Python 3 implementation, x4i3 [47], provide aninterface to the EXFOR library.Our philosophy on the management of experimentaldata in the EXFOR library was to convert the data fromthe EXFOR format to the JSON format first. The JSONformat enables the storage of hierarchically organized dataand supports strings and numbers. Thanks to these fea-tures, the logical structure of the EXFOR format and alldata stored therein can be perfectly preserved. Therefore,considering only the format itself, there is no advantageof the JSON over the EXFOR format to store the dataof nuclear experiments as they are equivalent in termsof logical structure and information. However, the JSONformat is a standard format commonly employed for theexchange of information between web servers and webbrowsers. Consequently, it is supported by a wide varietyof programming languages—in contrast to the EXFORformat.We used the R package exforParser [48] to translatethe EXFOR to the JSON format. As this parser does notmodify the logical organization of the data, it can probablycope robustly with future updates and extensions of the { "$and": [ { "BIB.REACTION": { "$regex": "\\(26-FE-56\\(N,[^)]+\\)[^,]*,,SIG\\)", "$options": "" }}, { "DATA.TABLE.DATA": { "$exists": true }}, { "DATA.TABLE.EN": { "$exists": true }} ]} Listing 1. Example of a search query in the MongoDB querylanguage. This query was used in the pipeline to retrieveneutron-induced cross sections of Fe.
EXFOR format, keeping the need to update the parseritself to a minimum.Besides the convenient extraction and manipulation ofexperimental data facilitated by the JSON format, it isalso important to have flexible capabilities of search tolocate relevant data. At the time of writing, more than23k experiments (=number of entries) and more than 150kreactions (=number of subentries) are recorded in theEXFOR library.Our approach to address this need was to create adocument-oriented database containing the JSON ob-jects with EXFOR data. The fundamental concept of adocument-oriented databases is a document, a bundle ofinformation belonging to one conceived entity. As an ex-ample, a document could store information about a person,such as their hair color, eye color and age. For comparison,a SQL database is based on the notion of tables where in-formation about an entity may be distributed over severaltables and one table usually contains partial informationof several entities.One of the popular document-oriented databases is
MongoDB , which enables the direct storage of
JSON ob-jects. It provides a flexible API (=application interface)to search, retrieve and modify data, which can be usedfrom a large variety of programming languages and fromthe command line.A script to convert the experimental data in the EX-FOR format to the JSON format and feed it into a Mon-goDB database has been made available in the repository createExforDb on GitHub [49]. A Dockerfile automatingthe full installation of the MondoDB EXFOR databaseincluding a manual is also available online [15]. More de-tails about the translation of the EXFOR library to aMongoDB database are available in [7].
1. Searching experimental data
After this general discussion on how we prepared theexperimental data in the EXFOR library, we provide abrief walk-through of searching and retrieving data toexemplify the convenient handling of experimental data.An example search query using the MondoDB querylanguage is shown in listing 1. Search queries are expressedas JSON objects. Search criteria are defined for specific15..Nuclear Data Evaluation Pipeline G. Schnabel et al. fields, such as the
REACTION field in the example. Inthe hierarchical representation of EXFOR, some fields aresubfields of other fields, e.g., the
REACTION field is asubfield of the
BIB (=bibliography) section. Fields deeperin the hierarchy are addressed by adding their parent fieldas a prefix, such as
BIB.REACTION in the example.A large variety of elementary search filters are availableand two are used in the example. The filter in line fouris a regular expression, e.g., [50]. A regular expression de-fines a search pattern that is matched against text. Itssyntax allows the specification of sophisticated patterns.For instance, the part N,[ˆ)]+ matches a sequence start-ing with the string N, and then matches a number ofone or several arbitrary characters that are not a closingbracket. This example provides just a small glimpse at thevast possibilities to define search patterns using regularexpressions. Of course, an end-user may not be botheredwith this syntax and common queries can and we thinkshould be packaged in functions or scripts, also executableon the command line.The second criterion employed two times in line sevenand eight in the example allows to match only documentsthat contain a specific field. In the case of the example,the document must contain the fields
DATA.TABLE.EN and
DATA.TABLE.DATA , i.e., the data table must con-tain a column with incident energies and the measuredquantity. Several other elementary search filters exist. Forexample, fields containing numbers can be matched basedon comparison operators, such as larger-than and smaller-than , and fields containing arrays can be matched basedon whether they contain a specific range of values. Thissearch filter can be useful if one is interested in experimen-tal data measured in a specific range of incident energies.Such elementary search filters can be combined usingthe logical operators and and or . The listing shows theapplication of the and operator specifying that all elemen-tary search filters must match at the same time. Morecomplex combinations of logical comparison operators arepossible.The outlined features of the query language hint at theflexibility to formulate search queries. As MongoDB andits query language can cope with EXFOR data in its orig-inal nested structure and the heterogeneous appearanceof data in different subentries, e.g., fields missing in onesubentry and present in another one, we think it is a goodand future-proof way to handle information in the EX-FOR library. The conversion of the EXFOR library into aMongoDB database can be done in dozens of minutes on acontemporary laptop, hence such a database can be easilykept in sync with the continuously evolving EXFOR li-brary. As an important remark, other implementations ofdocument-oriented databases than MongoDB exist whichprovide equivalent functionality. The EXFOR library hasalso been translated to a CouchDB database [51], whichcomes with a more permissive license.
2. Retrieving experimental data
The MongoDB API is accessible from a variety of pro-gramming languages and also from the command line.As the pipeline and all supporting packages are writtenin the programming language R, we sketch with shortcode snippets how experimental data can be retrieved inthis language. Access to a MongoDB database in R ispossible via the mongolite package [52] and a small pack-age
MongoEXFOR [53] was developed providing—in ouropinion—slightly more convenient access to EXFOR data.The latter package is a thin wrapper around the formerpackage. The following code examples demonstrate theaccess to information in EXFOR via the interface of the
MongoEXFOR package.There are two options for the retrieval of experimentaldata given in EXFOR subentries:1. A result can be returned as a collection of nestedlists , each of these containing the information ofone EXFOR subentry. A nested list in R is roughlyequivalent to a JSON object. A code skeleton toobtain a result in this format is given by it = exforIterator(queryStr)while (!is.null((curSub <- it$getNext()))) {
2. A result can be returned as a data frame , which is anExcel-sheet like table, summarizing the informationof several subentries. A code skeleton to obtain thedata frame is given by result = db$find(queryStr, {
Both options require the argument queryStr in the func-tion call, a string that contains the query specificationpassed to the MongoDB API. An example of a querystring was given in listing 1.The advantage of the first option being that the fullinformation of a specific subentry is obtained which canthen be processed before it enters the evaluation. Weused this option in the rule-based correction of uncertaintyinformation in the evaluation pipeline. In this task setting,the presence and absence of certain columns in the datatable of an EXFOR subentry determines the sequenceof actions taken. For instance, the existence of a totaluncertainty component is assessed via the R expression is.null(curSub$DATA$TABLE$`ERR-T`) and if present, the consistency with statistical and sys-tematic uncertainty components assessed. This expres-sion also demonstrates the convenience of access to anyinformation in an EXFOR subentry by specifying thechain of field names separated by a dollar sign to traverse16..Nuclear Data Evaluation Pipeline G. Schnabel et al. result = db$find(queryStr,{ exforid = ID reac = BIB$REACTION detector_raw = nullToNA(BIB$DETECTOR) detector = str_match(detector_raw, '\\(([^)]+)\\)')[2] energies = DATA$TABLE$EN measurement = DATA$TABLE$DATA list(ID=exforid, REAC=reac, DET=detector, EN=energies, XS=measurement) }) Listing 4. Retrieval of EXFOR data from several subentries inthe form of a data frame the hierarchical structure of the subentry. The backticksenclosing the last field name are necessary to avoid theinterpretation as a subtraction operation.The second option is to summarize the information ofseveral subentries in a data frame , a data structure akinto a table in Excel. We use the example code snippet inlisting 4 to discuss the features of this retrieval option.The find function evaluates the code between lines twoand nine enclosed by curly braces in the context of eachmatching subentry in the MongoDB EXFOR database.Differently stated, variables such as
BIB$REACTION inline three contain the
REACTION string of the subentrycurrently processed. Any valid R code can be executedin the context of the subentry, hence anything from basicmatrix operations to advanced statistical algorithms ispossible. As an aside, the ability to execute arbitrary Rcode in the context of a subentry also enables the imple-mentation of a correction system. R code can be storedalong the data in the subentry and when the data is ac-cessed, the R code can be executed first to correct andtransform the experimental data before retrieval. In thisway, corrections recommended by experts or those used inprevious evaluations can be automatically applied beforethe data is further processed for an evaluation.The code in line four and five is an example of slightlymore elaborate preprocessing. First, in line four it is at-tempted to retrieve the character string in the
DETEC-TOR field. Due to the variable being wrapped by thefunction nullToNA , a missing
DETECTOR field leads tothe assignment of the special value NA (=not available).Functionality as provided by the function nullToNA isvery important to deal robustly with missing values be-cause many fields in EXFOR subentries are not guaran-teed to be present but we may still want to take theirinformation into account if they are.To explain line five, we have to take a brief excursion tothe conversion of the EXFOR library. As the field contentof subentries has been preserved from the translation ofthe EXFOR library to a MongoDB database, so-called EX-FOR codes and free text descriptions may occur together.EXFOR codes are abbreviations put between brackets ref-erencing, e.g., detector types and institutions. For instance,the string (SCIN) in the field
DETECTOR indicates theutilization of a scintiliation detector in an experiment. Fields often contain free text that provides additional in-formation to humans, which is often the explanation ofthe EXFOR code. For instance, the
DETECTOR field ofan EXFOR subentry may contain (NAICR) NA-I SPEC-TROMETER .Returning to the discussion of the code snippet, thefunction str_match from the stringr package [54] is usedin line five to extract only the EXFOR code from the
DETECTOR field. Irrespective of the fact that the regu-lar expression passed as second argument looks dauntingto the uninitiated, this line of code demonstrates that aconcise instruction is enough to retrieve only the EXFORcode from a field—an important requirement to convertEXFOR data into formats understood by processing codes.As an aside, the effortless manipulation of strings and datain high-level languages such as Python and R was one ofthe reasons why the logical structure and field contentsof the EXFOR library were preserved in the translationto a MongoDB database. Derived formats and databasessuited for numerical processing can be easily produced us-ing the powerful and flexible facilities of these languages.Finally, we discuss the meaning of line eight and ninein the code snippet. After the retrieval of data and theirassignment to variables, the content of the resulting dataframe needs to be specified. The last statement list(...) in the code block enclosed by curly braces defines thecolumns, their name, and content of the resulting dataframe. An excerpt from the resulting data frame usingthe code snippet is illustrated in table V.
ID REAC DET EN XS10022010 (26-FE-56(N,P) . . . ) NAICR 14.60 113.0010031005 (26-FE-56(N,P) . . . ) NAICR 14.80 109.0010037004 (26-FE-56(N,EL) . . . ) SCIN 5.05 2129.0010037004 (26-FE-56(N,EL) . . . ) SCIN 5.58 1916.0010037005 (26-FE-56(N,TOT) . . . ) SCIN 5.05 3607.00Table V. Excerpt from the resulting data frame obtained bythe code snippet in listing 4. The full EXFOR reaction speci-fication is truncated in column
REAC for better display.
The intention of this walk-through was to give thereader an idea of the retrieval options and a glimpse onthe ease and flexibility to retrieve data. The functionalityof this IT building block puts an evaluator in the positionto retrieve relevant data by writing a few lines of code,which may be considered an essential requirement for theimplementation of an evaluation pipeline. Furthermore,the convenient extraction of features associated with ex-perimental data, such as the detector type in the example,in combination with the possibility of large scale Bayesiancomputation provided by the nucdataBaynet package dis-cussed in section III A enables the consideration of corre-lations between experiments due to detectors, monitors,institutions, etc. on an unprecedented scale in the field ofnuclear data evaluation. It also opens the door to conve-niently retrieve and prepare the data in a format requiredby advanced statistical algorithms and machine learningmethods.17..Nuclear Data Evaluation Pipeline G. Schnabel et al.
C. Mapping of model predictions to EXFOR
Nuclear model codes, such as TALYS [13] and EM-PIRE [55, 56], are capable to predict a large variety ofobservables. These codes write their predictions in theform of numeric tables to files. Sometimes the resultsof model codes cannot be directly compared to the ex-perimental data in the EXFOR library. This is the casefor measurements of sums and ratios of reaction crosssections. However, the availability of predictions for allmeasured quantities is essential in evaluations where a nu-clear model is employed as fitting function. Furthermore,because our objective are large scale evaluations with thevision to take into account all available and trustworthyexperimental data of the EXFOR library, the mapping ofmodel predictions to experimental data needs to be solvedin a general and extensible manner.In this section, we sketch the design and functionalityof the R package talysExforMapping [57] for this purpose.This package addresses the need to translate predictionsof the nuclear model code TALYS to experimental data inEXFOR. Even though the package is tailored to TALYS,its general philosophy and basic structure can be adoptedfor other codes as well, such as EMPIRE.It was clear from the beginning that the large varietyof observables stored in the EXFOR library in combi-nation with the complex EXFOR format and differentoptions within the format to represent the same informa-tion makes it unfeasible to come up with a comprehensivesolution to the mapping of predictions to EXFOR data inone shot. Consequently, the implementation of mappingfunctionality as an IT block was guided by the followingdesiderata:1. Support for the mapping of new observable typescan be incrementally added whenever the needarises.2. As it is not clear in advance what information maybe relevant to determine the appropriate mappingfor an observable, functions to implement a mappingfunction should be provided with all information inan EXFOR subentry.Desideratum (1) was addressed by employing the con-cept of handler modules, as it was already done in sec-tion III A for the management of information in theBayesian context. Each handler is responsible for themapping of predictions to a specific type of observable.Support for the mapping to a new observable type isthen established by the implementation of a new handler.Desideratum (2) was addressed by passing the completeinformation of an EXFOR subentry to the functions of ahandler module.A handler should be able to tell which predictions itneeds from a model code so that they can be mappedto the types of observables in EXFOR subentries. Ob-servables of interest in the subentries are referenced in adata frame called expDt in the following. The structure of this data frame with some example content is illustratedin table VI. The essential colums are
IDX , EXPID and
IDX EXPID DIDX REAC L1 DATA1 10529004 207 (26-FE-56(N,INL) . . . ) 1.00 1043.802 10529004 208 (26-FE-56(N,INL) . . . ) 1.00 967.603 10529004 209 (26-FE-56(N,INL) . . . ) 1.00 731.204 10529004 210 (26-FE-56(N,INL) . . . ) 1.01 568.50Table VI. Structure of data frame expDt with example content.Columns L2 and L3 are omitted for better display. DIDX . The purpose of column
IDX is to impose an orderamong the observables. The column
EXPID indicates theEXFOR ID of the subentry where the measured value ofthe observable is recorded. The column
DIDX indicatesthe row in the data table of the subentry containing themeasured value. The information in the other columns isredundant because it can always be retrieved from thereferenced subentry. The meaning of L1 , L2 and L3 isdependent on the EXFOR reaction string in REAC . Forinstance, the code
SIG at the end of the reaction stringdenotes cross sections and then the column L1 containsthe associated incident energy and the other columns L2 and L3 are ignored.The predictions required to map to the observables in expDt are referenced in a data frame needsDt . Such a dataframe with example content is shown in table VII. The IDX PROJ ELEM MASS REAC L11 N FE 56 CS/EL 1.502 N FE 56 CS/EL 2.003 N FE 56 CS/REAC/100000/TOT 1.004 N FE 56 CS/REAC/100000/TOT 1.005 N FE 56 CS/REAC/100000/TOT 1.00Table VII. Structure of data frame with example content todefine required TALYS predictions. The column names
PROJ and
ELEM are abbreviations of the full column names
PRO-JECTILE and
ELEMENT for better display. Columns L2 and L3 are omitted for the same reason. column REAC combines the definition of the reactionand measured quantity. The prefix
CS/ indicates a crosssection and other prefixes exist which are associated withangular distributions, energy spectra and double differen-tial cross sections. Afterwards follows the specification ofthe reaction. The naming convention mirrors closely thenames of the TALYS output files. The meaning of thecolumns L1 , L2 , L3 depends on the quantity. In the caseof a cross section, L1 indicates the incident energy and L2 and L3 are ignored.After the introduction of the data frame expDt referenc-ing observables of interest in a list of EXFOR subentries subents and the data frame needsDt referencing the re-quired model predictions, we discuss how the handlers areemployed. There is one manager module that keeps trackof the handlers and delegates the work to the appropriatehandlers. It provides the following interface:• canMap(subents) : Informs for each subentry in the18..Nuclear Data Evaluation Pipeline G. Schnabel et al. list subents whether a handler is available thatknows how to translate predictions of TALYS tothe observables in that subentry.• needs(expDt, subents) : Returns a data frame in theform of table VII with the predictions required fromTALYS to enable the mapping to the observablesstated in data frame expDt .• map(expDt, needsDt, subents) : Translates the pre-dictions of TALYS given in the data frame needsDt to the observables referenced in expDt . The dataframe needsDt is of the form as exemplified in ta-ble VII and must possess the additional column V1 with the associated predictions.The individual handlers provide the same functions as themanager module but with different parameters as theywork on the level of a single subentry. The manager mod-ule also provides other functions for registering handlersand querying information about them. From the point ofview of mathematics, it is often useful to have the trans-lation of the predictions in needsDt to the observables in expDt as a mapping matrix S . This functionality is alsooffered by the manager module.In summary, due to the extensible architecture of thepackage talysExforMapping and the fact that it deals withthe unaltered structure and information of the EXFORsubentry, it can probably organically co-evolve with theoriginal EXFOR library. D. Parallelization of nuclear model invocations
Modern nuclear data evaluation relying on a nuclearmodel code is computationally expensive. Evaluationmethods based on Monte Carlo sampling, such as UMC-G [58, 59], UMC-B [60], BFMC [61, 62] and BMC [63]require a large number of model code invocations withvaried model parameter sets to achieve convergence. Eval-uation methods based on optimization, such as [11] em-ploying the prior-aware LM algorithm and [64] using theKalman filter, require several invocations of the nuclearmodel to compute a gradient or Jacobian matrix. A con-ceptual overview of some of these methods is presentedin [65, 66]. A quantitative comparison of some of the meth-ods, also including a study of the impact of model defects,has been performed in [67].In the prototype of the pipeline, the joint optimizationof 147 TALYS parameters requires 148 model calcula-tions in each iteration and each model invocation takestens of minutes. Considering the computational demand,the parallel execution of nuclear model calculations on amulti-processor machine or computer cluster becomes anecessity.The implementation of the functionality of this IT build-ing block was guided by the following desiderata:• The implementation of parallelization functionalityshould be such that it is easy to set it up on a wide variety of different multi-core and cluster environ-ments.• The development of functionality to be parallelizedshould be possible on the local machine, i.e., func-tions and packages can be run both locally and re-motely on the multi-core or cluster environment inthe same way, which is good for testing and qualityassurance.• There should be as few restrictions as possible onwhat can be parallelized, i.e., everything from basictext processing to the invocation of a nuclear modelscode should be doable.Several solutions exist addressing parallelization of calcu-lations but none of those we are aware of fully addressedall points of our list. Scheduling engines, such as the SunGrid Engine (SGE) [68], manage the queuing of calcula-tion requests and distribution of calculations to computenodes. Based on our experience, not all computer clustersat scientific institutions have a scheduling engine installed.Of course, if it is available it can be relied upon. Softwareframeworks such as MapReduce [69] and Spark [70] en-able the processing of structured data, text and imagesbut we felt they are not the right tool to invoke nuclearmodel codes that produce several thousand output files ineach run which can be immediately discarded after execu-tion. Also the programming language R provides supportfor the parallel execution of R functions on chunks of in-put data by means of the parallel package part of theR language [42]. It uses sockets for communication be-tween compute nodes. The functionality of this packageaddresses many points of our list. Unfortunately, the useof sockets is not possible in some scenarios due to firewallrules.Considering our use case, which is running several in-stances of the same model code in parallel, the design ofour solution was based on the following assumptions:• Individual calculation requests are bulky taking afew minutes or longer so that latency due to networkcommunication or other overhead is negligible.• Individual calculations are independent of eachother so they do not need to exchange informationbetween each other.• All computational units use Linux, provide SSH ac-cess and have access to a shared file system.The last assumption is a reflection of our experience withscientific cluster environments.Several R packages were implemented building on eachother to enable the parallelization of computations. Thepackage interactiveSSH [71] and rsyncFacility [72] enablethe execution of commands in a bash prompt on a remotemachine and the transfer of files via SSH. Building on topof these packages, the package remoteFunctionSSH [73]enables the remote execution of functions written in R ina transparent way. The package clusterSSH [74] builds19..Nuclear Data Evaluation Pipeline G. Schnabel et al. remHnd = initSSH(...) fun = function(x,y) { x + y} remfun = remHnd$createRemoteFunction(fun) remote_result = remfun(2,4) local_result = fun(2,4) Listing 5. Definition of a function and its execution on a remotemachine on top of remoteFunctionSSH to enable the parallel exe-cution of R functions on chunks of data on a multi-coremachine or computer cluster. Finally, the package cluster-TALYS [75] enables the parallelization of TALYS calcula-tions and retrieval of the output produced. These packagesrely only on the remote execution of commands via SSHand the transfer of files to exchange data and code be-tween the local machine and the remote machine, hencethey enable parallelization as long as assumption threeabove holds. As a side note, in the meantime an R pack-age named ssh [76] has become available, which appearsto be a more robust implementation of the functionalityprovided by interactiveSSH .We are not going to discuss the technicalities of thesepackages here. Instead, we first highlight the functionalityof the package remoteFunctionSSH as it was a catalysatorto implement the required functionality for the paralleliza-tion of TALYS model calculations. Afterwards, we providea concise example to run TALYS calculations in parallelusing the clusterTALYS package.An example of the usage of the remoteFunctionSSH package is shown in listing 5. The setup of the remoteconnection is done in line one. The parameters of thefunction initSSH omitted for conciseness in the code snip-pet are the host address, user credentials and a directoryon the local machine and on the remote machine, whichare used for data transfer. Line two defines a functionnamed fun in the local R session. This function is trans-ferred to the remote machine and associated with the localvariable remfun . The invocation of the function remfun inline four executes the function fun remotely and returnsthe result as usual in the local R session. For comparison,the same function is also executed on the local machinein line five.The important feature is that from the viewpoint of theuser, there is no difference between invoking fun and itsremote counterpart remfun because both take exactly thesame parameters and return the same result. The abilityto design and test functions in the local R session andthen run them on a remote machine in the same way as ifthey were local functions is tremendously helpful for therapid development of functionality. This example featureda simple sum of two numbers but any valid R code canbe executed remotely, implementing anything from simplearithematical operations to the invocation of a schedulingengine, such as SGE, or the invocation of a nuclear modelcode.Finally we highlight the package clusterTALYS whichenables the parallel execution of TALYS on a computer talysClust = initClusterTALYS(...) paramset = list(projectile = "n", element = "Fe", mass = 56, energy = c(1,2,3) ) outSpec = data.table(REAC = "CS/TOT", L1 = c(1,2,3,4), L2 = 0, L3 = 0 ) paramsetList = replicate(3, paramset, simplify=FALSE) runObj = talysClust$run(paramsetList, outSpec, pollTime=30) result = talysClust$result(runObj) Listing 6. Setup and execution of several TALYS calculationsin parallel including the retrieval of predictions. cluster. A basic example of the usage of this package isshown in listing 6. An object is initialized in line one whichbundles the parallelization functionality. The parametersof function initClusterTalys are omitted for conciseness.In line two a list with TALYS input keywords and theirvalues is defined. Noteworthy, the energy keyword takesan array of numbers, which are written to a separatefile as required by TALYS during the setup of a TALYScalculation. Similar enhancements were introduced for en-ergy dependent model parameters. Specifications, such as "rvadjust(10) n" = 1.05 with the energy in bracket aretranslated to additional files with energy/value tables asrequired by TALYS. Line seven defines a data frame withthe predictions of interest, which is of the form as exem-plified in table VII without the columns
PROJ , ELEM and
MASS . The parameter set in paramSet is stored threetimes in a list paramsetList . TALYS is then launched inline 13 for each parameter set in paramsetList . The resultretrieved in line 15 is a list of data frames where eachdata frame corresponds to a calculation with a specificparameter set in paramsetList .This example demonstrated the seamless integrationof TALYS into the programming language and statisti-cal environment of R, regarding both parallel execution ofthe model code and retrieval of generated predictions. Theconvenient interaction with R is very beneficial for the cre-ation of an evaluation pipeline, as it facilitates the formu-lation of code logic at a high-level of abstraction. Further-more, advanced methods of statistics, machine learningand optimization available in R thanks to efforts of a largenumber of contributors can be easily applied on TALYS.Because R code can be run not only on the command linebut also in an interactive computer programming envi-ronment (REPL), such as RStudio [77], on-demand andhands-on data analysis are also facilitated.20..Nuclear Data Evaluation Pipeline G. Schnabel et al.
IV. THE NUCLEAR DATA EVALUATIONPIPELINE
The mathematical and IT building blocks discussed inthe two previous sections are a means to an end, whichis to provide best estimates and uncertainty informationof nuclear observables on the basis of information fromnuclear models and experiments. Generally, the evaluationprocess can be divided into the following stages:1.
Collection of experimental data.
Data suitableto be used in the evaluation have to be identified.Ideally, they are already added to the EXFOR li-brary where they can be retrieved via a targetedsearch. For more recent measurements, they haveto be extracted from the original publication or bydirect communication with the authors of the exper-iment.2.
Correction of experimental data.
Experimen-tal data need to be renormalized according to thelatest evaluations of monitor reactions. Outliers andinconsistencies in the experimental data have to beidentified and appropriate measures taken. Thesemeasures may be the exclusion of datasets, an in-depth analysis of original experimental setups to fig-ure out the problem with the data or an automaticcorrection relying on rules of thumb or algorithms ofstatistics and machine learning. Noteworthy, thesemeasures are not mutually exclusive. For instance,data can be manually corrected by human expertsas good as possible and remaining inconsistenciesaddressed by an automated statistical treatment.3.
Fitting a model to the corrected experimen-tal data.
The model may be a nuclear physicsmodel or a simple mathematical function, such asa piecewise linear function. To have the advantagesof both, which are good extrapolation in the caseof the nuclear physics model and great flexibility tomimic the trends in the data for the mathematicalmodel, both models may be combined to a meta-model, which is exactly the idea of model defects,e.g., [11, 22–24]. Predictions of this meta-model arethe sum of the nuclear physics model and the flexi-ble mathematical model. Another popular approachis to replace the exact nuclear physics model by aMVN distribution, as done in UMC-G [58, 59]. Amean vector and covariance matrix are estimatedfrom a sample of model predictions with varied pa-rameter sets. Nowadays, there seems to be a conver-gence to Bayesian methods for model fitting. Theyenable the inclusion of prior knowledge and result-ing estimates come always with associated uncer-tainties.4.
Generation of data files.
The final results of theevaluation need to be cast into formats understand-able by verification, processing and simulation codes.The format adopted worldwide by major institutions linked to nuclear physics is the ENDF-6 format [78].An emerging alternative to the ENDF-6 format isthe GNDS format [79].Some evaluations were informed by feedback from inte-gral experiments. The above mentioned sequence of stagesremains valid if one broadens the scope of the term ex-perimental data to include integral experiments as well.From the point of Bayesian statistics, the distinction be-tween so-called differential data and integral experimen-tal data is irrelevant. Any information can be included inthe Bayesian fitting procedure as long as the link to theother pieces of information, e.g., other experimental dataor model parameters, is mathematically well defined.Having said that, at present integral and differentialdata are considered differently for technical reasons. Inthe case of integral simulations, complicated processingbeyond the simple ACE format is often needed, such asself-shielding calculations, macroscopic group cross sec-tions tabulated as a function of specific parameters andbuckling factors. Such processing often relies on approxi-mations for computational tractability or for compatibilitywith the limitations of specific file formats. Ignoring theseeffects in a Bayesian evaluation using integral data bearsthe risk of adjusting differential data too much in an at-tempt to compensate for the unaccounted errors due toapproximations during processing. First attempts to ac-count for such effects are presented in [80, 81] using theMLO approach. The long runtime of neutron transportcodes for certain benchmark experiments is another obsta-cle to the inclusion of integral data in Bayesian methodsfor nuclear data evaluation.In the following we present our specific implementationof a nuclear data evaluation pipeline passing through thefour stages. Our implementation has to be regarded as aprototype without the pretense that selected experimentaldata, employed corrections or statistical algorithms arethe ultimate answer to everything.For instance, a shortcoming in the current version of thepipeline is that experimental data retrieved from the EX-FOR library are not renormalized according to the latestevaluations of the employed monitor reaction. More gen-erally, a variety of errors may affect measurement results,discussed in detail in [82]. Ideally, an in-depth analysis ofthe experiments included in an evaluation is performedto determine the error sources and their potential im-pact. However, this has not been done for this exampleevaluation. Another issue is that a crucial assumption ofthe statistical model (in the nuclear physics sense) imple-mented in TALYS does not hold for Fe. Because theproton number in Fe is close to a magic number, thelevel density is very low and effects of its structure canbe observed up to about 6 MeV. The statistical model as-sumes these effects to be completely averaged out whichis only the case if the average level width is much largerthan the average level spacing. No-model fits, e.g., basedon Gaussian processes [83], may provide better fits to theexperimental data in this energy range.Yet, we do believe that the coherent assembly of inno-21..Nuclear Data Evaluation Pipeline G. Schnabel et al. vative approaches, such as the correction of experimentaldata based on both rules and an statistical algorithm, andthe specification of Gaussian process priors on energy-dependent model parameters can serve as inspiration forfuture nuclear data evaluation efforts. As the pipeline ismodularly designed and can be run at the push of a button,modifications in any step can be gradually implementedand their impact systematically evaluated. For instance,another model for fitting than the nuclear model codeTALYS could be used in the pipeline. In the evaluation ofneutron standards [4] one prefers a flexible mathematicalfitting function, such as a piecewise linear model, to avoidpotential distortions of the result due to the rigidity of anuclear model.In the following sections we go through each step, see ta-ble VIII, and discuss how the functionality of the math-ematical building blocks provided by the nucdataBaynet package and the IT building blocks are applied in thepipeline. Each step is associated with a script in [14] whosefilename contains the step index as a prefix followed bya description roughly matching the one in the table. Thecorrespondence to the four general stages is as follows:Step one corresponds to stage one. Steps three to six areassociated with stage two. Steps seven and eight are im-plementing the stage three. Finally, step nine correspondsto stage four.The pipeline is open-source and has been publishedon GitHub [14] along with a manual with guidance oninstallation, usage, and the structure of the pipeline.
A. Retrieval of experimental data
In the first step, relevant experimental data were re-trieved from the MongoDB EXFOR database, see sec-tion III B for the discussion of this IT building block. Tothat end, we formulated a search query implementing thefollowing criteria:1. The
REACTION field must match the regular ex-pression:(26-FE-56(N,[ˆ)]+)[ˆ,]*„SIG) , i.e., neutron inducedcross sections of .2. Both the EN and DATA column must be presentin the data table of the EXFOR subentry.We are aware that the regular expression is certainlynot the most user-friendly way to define reactions to besearched. The augmentation of the MongoDB EXFORdatabase by additional fields containing the projectile,target, reaction process and quantity separately removesthe need for a regular expression. The MongoDB databaseis augmented with such fields but they have not been usedin this query.The obtained subentries compatible with the searchcriteria above are then only considered in the evaluationif: 1. They have valid uncertainty information, i.e., theycontain any of the fields
DATA-ERR , ERR-S , and
ERR-T .2. It is known how to produce a TALYS predictionthat can be related to the data in the subentry.Of the subentries fulfilling these criteria, only data pointsbetween 2 and 30 MeV are selected for the evaluation.The validity of uncertainty information is checked by afunction in a package named exforUncertainty [84]. Morefunctionality of this package is discussed in step 3 deal-ing with the rule-based correction of uncertainties. Thedetermination whether TALYS predictions can be relatedto the measurements in the subentries is done via thepackage talysExforMapping , see section III C for an expla-nation of this IT building block.
B. Generation of predictions based on a referencecalculation
A calculation of TALYS with default parameters is per-formed to obtain reference predictions of the observablespresent in the experimental data. These predictions arethe basis to translate relative experimental uncertaintiesto absolute ones to avoid the underestimation of evaluatedquantities—a manifestation of the phenomenon known asPeelle’s pertinent puzzle (PPP) [85] in the field of nucleardata. Excerpts of the informally distributed memorandumby Peelle can be found in [86].TALYS parameters are adjusted by making use of the adjust keywords in the TALYS input file, e.g., rvadjustn 1.05 . An effective parameter ˜ p is therefore the multi-plication factor applied to the default value p ref used byTALYS: p = ˜ p × p def (51)Furthermore, each parameter ˜ p was further constrainedto the range [0 . , . by applying the transformation ˜ p = x + δ × − exp( − k ( t − x ))1 + exp( − k ( t − x )) (52)and setting x = 1 , δ = 0 . , k = 4 . Making use of thetransformed parameters t , unconstrained optimization al-gorithms can be applied and the restriction on the rangestill enforced. For the sake of conciseness, any mentionof model parameters in the following discussion includ-ing other subsections always refers to the transformedones. As a side note, a variety of other transformationsis conceivable that may be better in certain situations. Aconvenient feature of the employed transformation is thata transformed parameter value equals approximately theoriginal parameter value if ˜ p is not too far away from one.Instead of a reference calculation to avoid PPP as donein the pipeline, there are other options to remove or atleast alleviate the undesired bias. For instance, a smoothcurve can be fit to the experimental data and used as the22..Nuclear Data Evaluation Pipeline G. Schnabel et al. Step Description reference in the conversion. Statistical fluctuations are’washed out‘ and consequently the bias reduced. This op-tion may be better if the shape of the reference predictionwas not well aligned with the experimental data.Another approach is to apply an iterative fitting strat-egy where the conversion of absolute to relative uncertain-ties in the current iteration is based on the intermediatefit of the previous iteration. This approach was suggestedby Chiba and Smith [86] and followed for the evaluationof neutron standards where two iterations of the Gener-alized Least Squares method were undertaken. Relativeuncertainties were converted to absolute ones based on theresult of the first iteration. Also in [11] the importanceof adjusting the covariance matrix due to the presenceof relative uncertainties during the iterative fitting pro-cedure was emphasized. This latter paper employs theprior-aware Levenberg-Marquardt algorithm outlined insection II D which relies on an additional damping term compared to the approach of Chiba and Smith. A simu-lation study presented in [28, chap. 6.3] gives indicationthat a damping term is indeed necessary to ensure con-vergence.
C. Rule-based correction of experimentaluncertainties
The information on uncertainties of experiments in EX-FOR is often incomplete. The idea put forward in [87]is to introduce or modify uncertainty components basedon pragmatic and reasonable rules. For instance, datasetswere discarded if the specified uncertainties were deemedtoo low; and an additional uncertainty component wasintroduced if only either a statistical or systematic un-certainty component provided. The rule-based approachenables the construction of complete uncertainty informa-tion.The rules implemented in the pipeline differ from thoseof [87] and should be considered as placeholders to be sub-stituted by better ones in the future. They have been bun-dled in the R package exforUncertainty [84]. We brieflyoutline the approach adopted for the rule-based construc-tion of uncertainty information here.First, an attempt is made to extract systematic andstatistical uncertainties. Statistical uncertainties are re- trieved from either the
ERR-S or the
DATA-ERR columnin the subentry. A situation of both fields being presentis not accepted and the subentry discarded from the eval-uation. If the statistical error is given in percent in thesubentry, it is converted to an absolute uncertainty usingthe predictions of the reference calculation obtained inthe previous step of the pipeline.Systematic uncertainties are retrieved from columnsnamed
ERR-0 , ERR-1 , etc. Uncertainties in one columnare assumed to be fully correlated, and uncertainties ofdifferent columns uncorrelated. The information whethera systematic uncertainty component is absolute or rela-tive is preserved. Moreover, if an uncertainty componentdeclared as absolute is a constant proportion of the exper-imental cross section in the
DATA column, it is convertedto a relative one.Once the available information on statistical and sys-tematic uncertainties is extracted from a subentry, poten-tially missing information is reconstructed according tothe following rules:• If only the total uncertainty component in column
ERR-T is available, both the statistical and system-atic uncertainty component are assumed to be 10%.In effect, the information in
ERR-T is completelyignored.• If either the statistical or systematic uncertaintycomponent is missing and the total uncertainty isavailable in the column
ERR-T , the missing uncer-tainty component will be reconstructed. Let us de-note by δ miss the missing uncertainty component,by δ avail the available one, and by δ tot the total un-certainty component. The missing uncertainty com-ponent is then reconstructed as δ miss = (cid:113) δ tot − δ avail . (53)A negative difference under the root indicates aninconsistency in the uncertainty information of thesubentry and is penalized by assuming 10% for themissing uncertainty component.• If statistical, systematic and total uncertainty com-ponent are provided simultaneously but the quantity ∆ := δ tot − δ sys − δ stat (54)23..Nuclear Data Evaluation Pipeline G. Schnabel et al. does not vanish, the absolute value of ∆ is addedas both a statistical and a systematic uncertaintycomponent.The amount of uncertainty introduced in response tomissing or inconsistent uncertainties is ad-hoc becausesomething had to be done. A better approach would be toperform a statistical analysis of uncertainties provided insimilar experiments to determine a reasonable value forthe missing uncertainty. The automated filling in of miss-ing data by relying on available data of similar cases isknown as imputation in statistics. Yet another approachwould be to base the introduction of additional uncer-tainty on templates that contain ranges for typical uncer-tainty components, e.g., detector calibration uncertainty.Such templates have been recently created for certainexperiments [88]. We also think that the validation of un-certainties using the physial uncertainty bounds methodpresented in [89] contains relevant ideas in this context. D. Correction of experimental uncertainties basedon statistics
Nuclear experiments are complex endeavors and it mayhappen that not all sources of uncertainties are alwaysrecognized. A discussion of this possibility in differentexperimental situations and effects in evaluations canbe found in [90]. The presence of unrecognized (or atleast unreported) uncertainties can be detected if severaldatasets from different experiments are inconsistent witheach other. The difference between the measured valuesof different datasets is then too large to be explainable bythe reported uncertainties.The most desirable resolution of such a situation wouldbe to figure out the forgotten uncertainty and take itinto account in the evaluation. This manual approach canbe time-consuming and the timely success or even beingsuccessful at all is not guaranteed. Still, in evaluations au-tomated to a large extent such as those of the TENDL li-brary [8], something must be done about inconsistent data.We feel that due to the complexity of proper uncertaintyquantification, the inclusion of an automatic treatment ofinconsistent data is also valuable in manual evaluationsof specific isotopes. Even if one does not accept the auto-matic correction of uncertainties in such evaluations, analgorithm to automatically correct uncertainties can beused as a safeguard informing the evaluator of missinguncertainties and their potential magnitude.The solution implemented in the pipeline is to intro-duce additional systematic uncertainties and adjust theirvalue based on the maximization of the marginal like-lihood (MLO) in eq. (5). Another possibility could bethe exclusion of inconsistent experimental datasets us-ing an outlier detection algorithm. We think algorithmsfor outlier detection and uncertainty correction can bebeneficially applied side-by-side in future nuclear dataevaluation pipelines. The experimental covariance matrix Σ exp is the resultof the rule-based approach implemented as the previousstep of the pipeline. As an important reminder, we neverexplicitly build up this matrix but always use the repre-sentation in section II B, i.e., Σ exp = S sys U sys S T sys + D . (55)The covariance matrix U sys is associated with systematicerrors and the covariance matrix D with statistical errors.To automatically determine the additional systematicuncertainties in the case of inconsistent data, we augmentthe experimental covariance matrix in the following way: Σ augexp ( (cid:126)λ ) = Σ exp + S extra U extra S T extra . (56)The matrix U extra is diagonal and each element λ i rep-resents the square of an additional normalization uncer-tainty associated with a specific dataset. The augmentedexperimental covariance matrix Σ augexp is therefore a func-tion of (cid:126)λ , which are regarded as fitting parameters in MLO.The element e ij in the i-th row and j-th column of S extra is non-zero if the i-th element of (cid:126)σ exp belongs to the exper-imental dataset associated with the extra normalizationuncertainty λ j . In this case, the element e ij is given bythe reference cross section obtained in the second stepassociated with the i-th element of the experimental mea-surement vector (cid:126)σ exp . This construction means that each λ i is a relative and fully correlated normalization errorfor a specific dataset.The experimental covariance matrix expresses how farmeasured values are expected to be away from the un-known true values. In an evaluation and for the appli-cation of the MLO approach, we always need a modelfor what the truth may look like. For the determinationof additional experimental systematic uncertainties, wedo not want to rely on a nuclear physics model. If sucha model were deficient, the adjusted values of the addi-tional uncertainties would not reflect unrecognized errorsin the experiment but rather deficiencies of the nuclearmodel. Therefore, instead of using a physics model, we usea piecewise linear model as defined in eq. (27). As long asthe mesh is dense enough, it can very well approximateany continuous function.We recall the important statement that a covariancematrix of the values at the mesh points P in combinationwith linear interpolation induces a Gaussian process, i.e.,the covariances between values at arbitrary energies canbe evaluated according to eq. (28). The mapping of thevalues at the energy mesh of the piecewise linear functionto the experimental energies is established by the matrix S mod whose construction was described in eq. (29). Thefull observational covariance matrix Σ obs can be writtenas Σ obs = Σ augexp + S mod PS T mod . (57)Again, as explained in section II B, this matrix is never ex-plicitly stored but always dealt with in the representation24..Nuclear Data Evaluation Pipeline G. Schnabel et al. Σ obs = SUS T + D (58)and being of the following form in the current step: U = U sys extra
00 0 P and S = S sys U extra S mod T . (59)There is a strong prior expectation that cross sectioncurves are smooth above a certain incident energy. In suchregions, we believe that non-expert humans are mostlyguided by the shape of an envisaged cross section curveunderlying the experimental data to make a judgement ofwhether datasets are consistent or not. In other words, webelieve that neither the absolute values of this envisagedcross section curve nor the steepness of its slope matter.What matters is its smoothness, i.e., the rate of change ofthe slope. If too rapid variations in the cross section curveare required to be within experimental uncertainties closeto the points of the various datasets, the datasets wouldbe deemed inconsistent.One way forward to formalize this notion about reason-able shapes is to softly constrain the second derivative ofthe piecewise linear function to fit the cross section curve.We briefly sketch the essential ideas here. Appendix C 2contains a more rigorous discussion in terms of equations.In order to construct the corresponding prior covariancematrix P we consider the formula for numeric differentia-tion, d f ( E i ) d E ≈ f ( E i +1 ) − f ( E i ) E i +1 − E i , (60)where the energies E i and E i +1 are neighbors on the re-duced energy grid as introduced in the context of eq. (23)for the sparsification of the Gaussian process. The nu-merical differentiation rule defined in eq. (60) is a linearoperator and can thus be rewritten as a matrix multipli-cation. Let (cid:126)f denote the vector of function values at theenergies of the mesh and (cid:126) ∆ the vector of first derivativesat the energies E , E , . . . , E M − of the mesh. We canthen construct a matrix S diff such that (cid:126) ∆ ≈ S diff (cid:126)f . Im-portantly, according to the prescription in eq. (60), wecannot evaluate the first derivative of the function at thehighest energy of the mesh E M .The vector of second derivatives (cid:126) ∆ can be calculatedby another application of S diff . With a slight abuse ofnotation because the dimension of S diff is reduced in thesecond application, this can be written as (cid:126) ∆ = S diff S diff (cid:126)f .Due to S diff having one more column than row, this ma-trix cannot be inverted. This is a reflection of the factthat integration—the inverse operation to differentiation—is only unique up to an integration constant. To make S diff invertible, we augment it by inserting one row of theform (1 , , . . . , at the top so that S augdiff S augdiff (cid:126)f yields an augmented vector (cid:126) ∆ aug = ( f ( E ) , ∆ , (cid:126) ∆ ) T . The first ele-ment is the function value f ( E ) , the second element ∆ the first derivative at f ( E ) , and each subsequent elementat position i represents the second derivative at f ( E i − ) .Because S augdiff is a square matrix and invertible, the origi-nal parameters can be expressed by the transformed ones,i.e., (cid:126)f = R (cid:126) ∆ with R = ( S augdiff S augdiff ) − . The sandwich for-mula yields the covariance matrix P associated with theoriginal parameters as a function of the covariance matrix ˜ P of the transformed parameters, i.e., P = R ˜ PR T .For simplicity and computational efficiency, we assume ˜ P to be diagonal. We set the prior variance for f ( E ) andits first derivative to large values to express our ignoranceabout the cross section value and the slope. The determi-nation of the variances associated with second derivativeswere guided by the cross section at the peak of the curve.These assignments were, however, somewhat ad-hoc anda more principled approach for their determination meritsinvestigation in the future.Due to the fact that the prior covariance matrix of thepiecewise linear model is diagonal in the space of secondderivaties and there is a bijective and linear mappingbetween the function values at the locations of the meshand the values of the second derivaties at the locationsof the mesh, it is computationally more efficient to usethe second derivatives as model parameters. Therefore, ineq. (59), the block P in the covariance matrix U needs tobe substituted by ˜ P and the block S mod in the Jacobian S mod needs to be replaced by S mod ( S augdiff S augdiff ) − .At this stage we have discussed the construction of allparts of the observational covariance matrix Σ obs specifiedin eq. (57). With the additional specifications (cid:126)µ obs = (cid:126) and (cid:126)x obs = (cid:126)σ exp , we can evaluate the marginal likelihoodgiven in eq. (5) and adjust the additional normalizationuncertainties in (cid:126)λ (which determine U extra ) to maximizeit. Importantly, the systematic errors of datasets in dif-ferent reaction channels were assumed to be independent.This assumption enabled the application of MLO to eachreaction channel in isolation.The important result of this step is the term in Σ obs that corresponds to the augmented experimental covari-ance matrix Σ augexp , see eq. (56), based on the adjustedadditional systematic uncertainties (cid:126)λ .Finally, we want to hint at possible improvements ofthis approach for uncertainty correction. Instead of man-ually designing a Gaussian process using an intuition ofhow humans would identify problematic data, one couldalso infer the properties of the Gaussian process in a data-driven way from a collection of reaction systems. This ideawas explored in [31] with the aim to determine a modeldefect prior for reaction channels, which can then be usedin an evaluation procedure. The Gaussian approach con-struction described there could also be a way forward toextend the idea of MLO uncertainty correction to energyranges where resonant-like structures start to emerge.Another idea hinting at a potential improvement waspresented in the context of validation [91, 92]. There, var-25..Nuclear Data Evaluation Pipeline G. Schnabel et al. ious quantities derived of cross section curves, such as thepeak cross section, were regarded as functions of quanti-ties associated with properties of the nuclei, such as theasymmetry. These functions were fitted globally to a col-lection of experimental data and the empirical distributionof deviations from these derived quantities constructed.Problematic reaction systems were then identified as out-liers with respect to this empirical distribution. We thinkthat such an approach could be turned around and notonly used for validation but also for uncertainty correctionin the future. E. Evaluation of the Jacobian associated with thereference calculation
At this stage of the pipeline, we have a fully specifiedexperimental covariance matrix Σ augexp required for the sta-tistical inference carried out in subsequent steps. Thisstatistical inference comprises the determination of hy-perparameters associated with the Gaussian processes forenergy-dependent model parameters in the next step andthe adjustment of TALYS model parameters thereafter.Both steps mentioned rely on the Jacobian matrix J of thenuclear model whose elements are given by J ij = ∂σ i /∂p j with σ i being the prediction associated with the i-th ele-ment of the experimental measurement vector (cid:126)σ exp and p j being the j-th TALYS model parameter. More detailson how the Jacobian matrix is going to be used in thesesteps are explained then.The numerical evaluation of elements in the Jacobianmatrix is based on the formula ∂σ i ∂p j = σ i ( (cid:126)p ref + hp i (cid:126) ∆ i ) − σ i ( (cid:126)p ) hp i . (61)The only non-zero element of the vector (cid:126) ∆ i at index i is one. Using this specification, the step size h of thefinite difference approximation denotes the magnitude ofthe perturbation as proportion of the reference parametervalue. Too large values of h cause the numerical derivativeto be a worse approximation to the exact derivative. Toosmall values of h on the other hand make the result of thenumerical evaluation more prone to round-off errors asthe ratio of two very small numbers is taken. The valuethat has been employed in the pipeline is h = 0 . , i.e., aone percent perturbation of each TALYS parameter.Please note a more accurate formula for the numericdifferentiation would be ∂σ i ∂p j = σ i ( (cid:126)p ref + h(cid:126) ∆ i ) − σ i ( (cid:126)p − h(cid:126) ∆ i )2 hp i (62)because the lowest-order term of the error cancels out.However, this formula requires approximately twice asmany evaluations of the physics model as the formulain eq. (61). Because about thousand model parameterswere considered adjustable in our example evaluationof neutron-induced cross sections of Fe , we relied on eq. (61) to cut down on computational resources needed.Despite this reduction of the computational cost, we stillhad to rely on a scientific cluster with about 150 process-ing cores to evaluate the Jacobian matrix, which took afew hours. The parallelization of model calculations wasdone with the help of the clusterTALYS package discussedas an IT building block in section III D.
F. Setup of Gaussian processes forenergy-dependent model parameters
Nuclear models are often if not most of the times defi-cient. This statement probably sounds rather rude to theears of theoretical physicists who work hard to developmodels with more accurate descriptions of physics pro-cesses. We acknowledge these efforts and deeply respectthese people pushing the limits, yet we are uncompromis-ing concerning our statement.A model is potentially deficient if we have reason to be-lieve that its predictions cannot perfectly reproduce thetruth. Therefore even a very good model with an excel-lent capability of extrapolation and predictions within afew percent of trustworthy data is potentially deficient ifthe differences between predictions and data cannot beexplained away by experimental uncertainties.Work to incorporate model defects into nuclear dataevaluations was presented in [22, 23, 93, 94] and laterthe connection to Gaussian processes identified and for-malized [24, 28]. Some work dealing with the notion ofmodel defects in a more general setting outside nucleardata predates the efforts in the field of nuclear data, suchas [29, 95].We implemented an idea proposed and explored in [9]to use the capability of TALYS to allow variations of pa-rameters as function of incident energy to simulate thetreatment of model defects. Effectively, this measure in-jects more flexibility into the nuclear physics model and itbecomes capable of capturing a broader range of shapes.The advantage of this approach being that physics con-straints are automatically taken care of by the nuclearmodel code TALYS.The idea is the following. Nuclear model parameters,such as the real radius r v of the optical model, by de-fault assumed constant for a specific isotope, are re-garded as functions of incident energy, e.g., r v ( E ) . Thesefunctions are assumed to be piecewise linear functions,i.e., the function is characterized by the function val-ues on a mesh of energies. Intermediate values are ob-tained by linear interpolation. Therefore, each energy-dependent model parameter in the model parameter vec-tor (cid:126)p is represented as a collection of model parame-ters, i.e., the values of the paramter at a specific in-cident energy. In other words, we expand a parame-ter p k in (cid:126)p = ( . . . , p k , . . . ) T ) assumed to be energy-dependent as (cid:126)p = ( . . . , p k ( E ) , p k ( E ) , . . . , p k ( E M ) , . . . ) with E , E , . . . , E M being the energies of the predefinedmesh. As a reminder, the model parameters p k here are26..Nuclear Data Evaluation Pipeline G. Schnabel et al. the transformed ones defined in eq. (52).The covariance of parameter values at different meshpoints are defined by a squared exponential covariancefunction, κ k ( E i , E j ) = Cov [ p k ( E i ) , p k ( E j )] = δ k exp (cid:18) − λ k ( E i − E j ) (cid:19) + η ij . (63)The amplitude δ k represents the apriori standard devi-ation of the parameter values and is considered to beindependent of incident energy. The length-scale λ k is ameasure of similarity of a parameter at nearby incidentenergies. The so-called nugget parameter η ij is assigned asmall numerical value, e.g., − , if E i and E j are equaland zero otherwise. It is introduced for numerical stabilityand its impact on the result is usually negligible. Theseparameters of the covariance function are called hyper-parameters. In the evaluation of Fe in the pipeline, weallow the amplitude and length-scale to be different foreach energy-dependent model parameter.At this point it becomes clear that energy-dependentTALYS parameters are not conceptually different fromenergy-independent ones from an evaluation point of view.They are represented by values in the parameter vector (cid:126)p and the covariance functions are employed to computethe associated blocks in the prior covariance matrix P .The determination of the hyperparameters δ k and λ k isdone via marginal likelihood optimization—the same ap-proach that has already been employed in the determina-tion of additional systematic uncertainties of experimentsin step four discussed in section IV D. Therefore, we dealagain with the observational covariance matrix in the rep-resentation of eq. (58). The important difference is that P is now the prior covariance matrix associated with nuclearmodel parameters and not the function values of a piece-wise linear function at the mesh points as it was in stepfour. Analogously, each element of the Jacobian matrix S mod represents now the derivative of a prediction σ i withrespect to a model parameter p j , i.e., ( S mod ) ij = ∂σ i /∂p j .This Jacobian matrix has already been computed in stepfive.In more detail the MLO approach is carried out inthe following way: We adjust the hyperparameters in thecovariance functions in eq. (63) which alter the prior co-variance matrix P in eq. (59) and consequently the valueof the logarithmized marginal likelihood in eq. (5). To de-termine the values of the hyperparamters to maximize themarginal likelihood, we use the L-BFGS algorithm [26].The package nucdataBaynet supports these operations. G. Optimization of TALYS parameters using theLM algorithm
At this stage of the pipeline, all covariance matrices re-quired for an evaluation are determined. The experimentalcovariance matrix was constructed using a rule-based ap-proach in step three and additional relative systematic uncertainties introduced in step four to account for unrec-ognized or unreported systematic uncertainties. The finalexperimental covariance matrix was denoted as Σ augexp . Inthe previous step—step six—we determined the hyperpa-rameters of the Gaussian processes imposed on energy-dependent model parameters, which were necessary tocompute the associated elements of the prior covariancematrix P .The posterior distribution π ( (cid:126)p | (cid:126)σ exp ) ∝ (cid:96) ( (cid:126)σ exp | M ( (cid:126)p )) π ( (cid:126)p ) (64)with (cid:96) ( (cid:126)σ exp | M ( (cid:126)p )) = N ( (cid:126)σ exp | M ( (cid:126)p ); Σ augexp ) (65) π ( (cid:126)p ) = N ( (cid:126)p | (cid:126)p ; P ) (66)is therefore fully determined. The predictions of the nu-clear model TALYS are denoted by M ( (cid:126)p ) and this modellink is non-linear.In order to locate the parameter vector (cid:126)p max associ-ated with the maximum of π ( (cid:126)p | (cid:126)σ exp ) , we employed themodified Levenberg-Marquardt algorithm outlined in sec-tion II D.Importantly, parameters for adjustment were selectedbased on their sensitivity to experimental data. Each col-umn in the sensitivity matrix S reflects the impact of anadjustment of a specific parameter on the predicted ob-servables. A parameter was only taken into account foradjustment if the maximum of the values in the corre-sponding column S squared was larger than one. Thiscriterion reduced the number of parameters from roughlythousand to about 150 parameters.As a final remark, an alternative formulation of theparameter exclusion criterion in terms of relative changesof cross sections induced by relative changes of modelparameters may be easier to work with and to develop anintuition about it. H. Calculation of a MVN approximation of theposterior pdf
After the localisation of the model parameter vectorthat corresponds to the maximum of the posterior distri-bution in the previous step, associated uncertainty infor-mation needs to be determined as well. Due to the modellink between parameters and observables being non-linear,neither the posterior distribution of parameters nor modelpredictions has to be multivariate normal. To deal witha general distribution, Monte Carlo methods, such as im-portance sampling [96] and Markov Chain Monte Carlo(MCMC) [97], can be used to obtain representative sam-ples of it. Alternatively, the shape of the posterior distri-bution can be approximated by an analytic expression.We opted for the latter option and approximated the pos-terior distribution by a multivariate normal distributionas described in section II E. Even though this decision wasborn out of the necessity of computational tractability, wethink it may be defensible on the following grounds:27..Nuclear Data Evaluation Pipeline G. Schnabel et al. • Whenever experimental data strongly constrain theparameters, the model link is in good approximationlinear within the ranges of parameter values allowedby the posterior uncertainties. Due to the assump-tion of a multivariate normal prior and likelihood,the posterior is then also in good approximationmultivariate normal.• The posterior distribution of model parametersweakly constrained by the experimental data is closeto the prior distribution, which we specified as mul-tivariate normal.• Stronger deviations of the approximative posteriorparameter distribution from the exact distributionhave to be expected for intermediate cases. If theobjective is not unbiased parameter estimation butuncertainty quantification on the observable side,which is typically the case in nuclear data evaluation,the distortion of the distribution on the observableside is not critical as long as posterior uncertaintiesare somewhat compatible with those associated withthe exact distribution.The last point may be considered an optimistic intuition.Certainly mathematical toy scenarios can be constructedwhich invalidate this assumption. It would be interestingand pertinent to study possible distortions in evaluationsdue to such an approximation in typical evaluation sce-narios. However, such studies are outside the scope of thispaper.As described in section II E, the computation of the fullHessian matrix required to obtain the posterior parametercovariance matrix can be unpractical. Its computation forabout 150 model parameters requires more than ten thou-sand model invocations, each lasting tens of minutes. Toimprove the time complexity from N with N being thenumber of model parameters, we exploited the approxi-mation briefly described in section II E and in more detailin appendix C 5 and appendix C 4.Consequently, the second order Taylor approximationdescribed in section II E leads to the following approxima-tion of the posterior distribution, π ( (cid:126)p | (cid:126)σ exp ) = N ( (cid:126)p | (cid:126)p max , ˜ P ) , (67)with (cid:126)p max being the parameter vector associated with themaximum of the posterior distribution and ˜ P an approx-imate posterior covariance matrix constructed accordingto the prescription detailed in appendix C 4. I. Generation of representative random files
Two options are popular to perform uncertainty quan-tification in nuclear applications. One option is to propa-gate the covariance matrices associated with observablesin evaluated nuclear data files (ENDF files) to the quanti-ties of interest in the specific application. The other optionfollowed, e.g., in the Total Monte Carlo approach [98] is to rely on a representative sample of ENDF files where eachfile contains the values of observables obtained by a drawfrom the posterior distribution. Each ENDF file is thenemployed in the simulation of the nuclear application andthe resulting empirical distribution of the quantities ofinterest analyzed, e.g., by plotting histograms or by com-puting summary statistics such as the mean vector andstandard deviation.The pipeline produces both a collection of random filesand a best ENDFF file with covariance matrices. For thispurpose, a sample of parameter vectors is drawn from theapproximative posterior distribution in eq. (67) obtainedin the previous step. We stress here that the probabilitydistribution constructed in step eight to draw samplesand to extract a covariance matrix is different from theone used in the TMC approach applied for TENDL.TALYS calculations are then performed for each of theseparameter sets. The TEFAL code converts the output ofeach calculation to a corresponding ENDF file. A modifiedversion of the TASMAN code loops over the outputs of theindividual calculations and generates a ‘best’ ENDF filethat also includes covariance matrices. Our modificationof TASMAN [99] enabled the use of precalculated TALYSoutputs instead of letting TASMAN generate a sample ofparameter sets and invoke TALYS. TALYS, TEFAL andTASMAN are part of the T6 evaluation system [13].One may question the necessity of sampling for the con-struction of an ENDF file. For the multivariate normalapproximation of the posterior distribution, we arguedthat this approximation is justified if posterior uncertain-ties constrain the model parameters to a domain wherethe nuclear model is approximately linear. In this case lin-ear error propagation would be sufficient. However, evenif the model link between observables of the experimentaldata and the selected parameters used in the LM fittingis linear, the model link between predictions of reactionchannels without experimental data and model parame-ters may not be linear.Irrespective of this theoretical justification, the primarymotivation to rely on this procedure for the generationof an ENDF file with covariance information was con-venience. In the view of a constrained time budget, wewanted to avoid the need to implement functions to aug-ment an ENDF file with covariance information. The TAS-MAN code can take care of this but derives the covariancematrix from a sample of model calculations. We are op-timistic that the emerging GNDS format [79] is going toimprove the management of evaluated nuclear data.
V. DISCUSSION OF PRACTICALITIES ANDRESULTS
Even though the focus of the paper is on the conceptuallevel, we want to present results to convey an idea ofthe working of the algorithms and features of the finalevaluation. We also discuss practicalities to demonstratethe feasibility of the employed algorithms.28..Nuclear Data Evaluation Pipeline G. Schnabel et al.
For the evaluation of neutron-induced reactions of Fe,the search of experimental data in the EXFOR libraryaccording to the search specification in step one, see sec-tion IV A, yielded data for nine reaction channels, whichare (n,2n), (n,a), (n,d), (n,el), (n,inl), (n,n+p), (n,p), (n,t)and (n,tot), amounting to a total number of 4333 datapoints. Most data points (3895) are associated with (n,tot)and most of the remaining data points (324) are associatedwith (n,p).The effect of the automatic addition of normalizationuncertainties to experimental datasets via the MLO ap-proach performed in step four (section IV D) is shown atthe example of the (n,2n) reaction channel in fig. 2. Toobtain this result, we fixed the dataset with EXFOR ident-fication number 23171003 to have a statistical uncertaintyof 5% and a normalization uncertainty of 0.5%. With thisimposed assumption, the MLO approach introduces anadditional normalization uncertainty to one dataset thatalso visually appears as an outlier.Without this manual fix of the uncertainty of onedataset, the MLO approach would add the additional nor-malization uncertainties to the datasets we intuitivelyregarded as consistent, see fig. 3. Moreover, in the regionaround 15 MeV very good and almost mono-energetic neu-tron sources are available so that typical uncertainties ofactivation measurements are between two and five per-cent. Such a large increase of uncertainties is thereforenot plausible from an experimental point of view.This behavior of the MLO approach certainly needs fur-ther exploration and a more thorough study in order to gotowards more automated evaluations and less human in-volvement. Regarding the current implementation of theMLO approach, the more uncertainty information is al-
Figure 2. Example of the automatic correction of experimentalsystematic uncertainties based on the MLO approach in stepfour of the pipeline. The error bars attached to the data pointsare total uncertainties, i.e., including statistical and systematicuncertainties, after the rule-based correction. The extensionof the error bars of the blue dataset indicated by the greenextension line represents the increase of the total uncertaintydetermined by the MLO approach due to the inconsistencybetween datasets. Figure 3. Example of the automatic correction of experimentalsystematic uncertainties based on the MLO approach if notfixing the uncertainty of dataset with EXFOR identificationnumber 23171003. The extension of the error bars by the greenlines represents the increase of the total uncertainty due to anadditional normalization uncertainty determined by the MLOapproach. ready precisely provided by a human expert, the higher thechance that the MLO approach determines an assignmentof extra uncertainties compatible with our expectations.Therefore, the codification of choices imposed by a hu-man and the use of statistical algorithms to correct dataare not mutually exclusive. In such cases, an importantrequirement of a statistical algorithm is that extra uncer-tainties should only be introduced if the data are indeedinconsistent, i.e., the human expert has missed sources ofuncertainty.Irrespective of the difference in the results with andwithout fixing the uncertainty of a specific dataset, theMLO approach seems in both cases to be able to cor-rect uncertainties to establish consistency among datasets.This statement holds in general if the parametrization ofthe covariance matrix is flexible enough. Specifically inthe current implementation of the pipeline, the adjustableparameters associated with the covariance matrix allowfor the introduction of additional normalization uncertain-ties. If the inconsistency in the data is not a mere normal-ization issue, the solution found by the MLO approachwill therefore not be ideal. Noteworthy, as demonstratedin [10], the MLO approach is not limited to the adjust-ment of normalization uncertainties but any other type ofuncertainty, e.g., systematic uncertainties partially corre-lated over energy, can be determined as well. Despite thispossibility, in the case of a severely deficient dataset, itis probably better to remove it, potentially automaticallyusing an outlier detection algorithm.After the determination of extra experimental normal-ization uncertainties, the Gaussian process hyperparam-eters of energy-dependent model parameters were deter-mined via the MLO approach in step six, see section IV F.The optical model parameters v , d , w , v so , w so , r C forneutron, proton, deuteron, tritium, Helium-3, and alpha29..Nuclear Data Evaluation Pipeline G. Schnabel et al. are regarded as energy dependent. For further explana-tion of these parameters consult the manual accompany-ing the nuclear reaction code TALYS. The energy mesh ofeach of these 36 energy-dependent parameters consistedof 16 energies uniformly spaced between 0 and 30 MeVamounting to a total number of 576 effective parametersassociated with energy-dependent parameters. The adjust-ment of these parameters to maximize the logarithmizedmarginal likelihood in eq. (5) requires the Jacobian matrix(=sensitivity matrix) between model parameters and ob-servables measured in the experiments, computed in stepfive, see section IV E. To obtain the Jacobian matrix of di-mension × due to 4333 experimental data pointsand 925 model parameters, 926 model calculations hadto be carried out, each lasting about 20 minutes. On thecomputer cluster at the division of applied nuclear physicsat Uppsala university, we performed about 150 calcula-tions in parallel, hence the execution time was around twohours.Once the Jacobian had been computed, the adjustmentof the 576 parameters associated with energy-dependentparameters to maximize the logarithmized marginal like-lihood in eq. (5) took a bit more than a minute (withoutparallelization) using the L-BFGS-B algorithm [26]. Weforced all length-scales λ k to be in the range between2 and 50 MeV and all amplitudes δ k to be in the rangebetween 0.1 and 0.5. In other words, the uncertainty asso-ciated with energy-dependent parameters was bound tobe between 10% and 50%. We imposed the minimum onthe length-scales to make the optimization process robustagainst potential kinks in the experimental data. The min-imum imposed on the amplitudes prevented the completeremoval of uncertainties associated with energy-dependentparameters insensitive to the experimental data. This min-imum mirrors the a priori uncertainty of 10% imposed onenergy-independent TALYS model parameters. We wantto stress here that these adjusted uncertainties have tobe regarded as prior uncertainties even though they wereinformed by the experimental data. The determinationof hyperparameters by using the experimental data toobtain a fully specified prior is at the heart of empiricalBayesian methods [100, 101] where it is interpreted as anapproximation to a full Bayesian treatment.The short runtime of a bit more than a minute to deter-mine the 72 Gaussian process hyperparameters by opti-mizing the logarithmized marginal likelihood despite theoccuring inversion of a covariance matrix of dimension × in each iteration is possible thanks to thedecomposition of the covariance matrix presented in sec-tion II B and the exploitation of the matrix identitiesoutlined in appendix B. Especially, the availability of ananalytical gradient of the logarithmized marginal likeli-hood thanks to the identities in appendices B 1 and B 2significantly speeds up the computation compared to anumeric evaluation of the gradient. All of this function-ality has been made available in the R package nuc-dataBaynet [16] whose mathematical underpinnings werediscussed in section II and technical implementation as- Figure 4. Posterior expectations and uncertainties of theenergy- independent model parameters most constrained bythe experimental data. Parameters are sorted according toposterior uncertainty. pects in section III A.Before the invocation of the LM algorithm to adjustTALYS model parameters based on the prior knowledgeabout experimental data and parameters, the Jacobianmatrix obtained in step five, see section IV E, was usedto exclude insensitive parameters from adjustment. A pa-rameter was excluded if the maximal value in the cor-responding column in the Jacobian matrix was smallerthan one. This criterion reduced the number of param-eters from 925 to 147 parameters. The LM algorithmrequired about 15 iterations to obtain optimized valuesfor the 147 model parameters. The Jacobian associatedwith these parameters needs to be recomputed in everyiteration, which amounts to 148 model invocations periteration. These model invocations were performed in par-allel on the computer cluster. Due to a time requirementof about twenty minutes per model calculation, the LMalgorithm terminated after about five hours.The posterior approximation evaluated in step eight,see section IV H, relies on the evaluation of the diagonalelements of the Hessian matrix associated with the loga-rithmized marginal likelihood. The numeric evaluation of30..Nuclear Data Evaluation Pipeline G. Schnabel et al.
Figure 5. Posterior expectations and uncertainties of theenergy- dependent model parameters most constrained by theexperimental data. the second-derivatives of the marginal likelihood with re-spect to the 147 model parameters required × modelinvocations. Once these second-derivatives are available,the computation time to obtain the approximate posteriordistribution for all 925 TALYS model parameters accord-ing to the approach in appendix C 4 and appendix C 5 isnegligible.The resulting posterior expectations and uncertain-ties for energy-independent model parameters most con-strained by the experimental data are shown in fig. 4. Wesee that the posterior uncertainty of many parameters,e.g., gpadjust 25 55 , is not significantly reduced comparedto the 10% prior uncertainty. The posterior expectationsand uncertainties of the energy-dependent model param-eters most constrained by the data are shown in fig. 5.Among the 36 energy-dependent parameters, only a smallproportion is significantly informed by the experimentaldata. The visible spike in d1adjust n at about 5 MeV isa reflection of the parameter transformation explainedin section IV B to constrain the parameter values to theinterval [0 . , . .Regarding the posterior expectations and uncertaintiesof the model parameters most informed by the experimen- tal data, the inclusion of 147 model parameters in the LMoptimization appears to be very cautious and in futurelarge scale evaluations the criterion to exclude insensitivemodel parameters can be designed to be more aggressive.The posterior for the observables is depicted in fig. 6.The fits appear to be reasonable and follow well the trendsin the experimental data. Many of the posterior uncer-tainties appear to be very small, maybe too small to beaccepted in a real evaluation. This may be indication thatGaussian processes on energy-dependent parameters arenot always enough to effectively address the deficiency ofnuclear models, and the introduction of model defects onthe observable side needs to be considered as well.From a physics point of view, it is questionable whetherthe measurements in the EXFOR library tagged as(n,n+p) and (n,d) really measure different reaction pro-cesses. As the deuteron is a weakly bound particle, theQ-value for the (n,n+p) and (n,d) reaction are not verydifferent and consequently these reactions cannot be dis-cerned using the activation method. The extensible designof the IT building block discussed in section III C enablesthe introduction of more sophisticated handlers that takeinto account the measurement method and the reactionstring to determine the appropriate mapping from modelpredictions to the experimental data. In the example ofthe (n,n+p) and (n,d) reaction, the appropriate mappingwould then be to use the corresponding residual produc-tion cross section instead of the exclusive reaction crosssections if the experiment relied on activation.The experimental data included in the example eval-uation spanned the energy range from 2 to 30 MeV, anumber of 4333 data points. To get a better grasp of whatis possible in future evaluations, we also launched thepipeline with an extended collection of experimental databy including also data in the range from 0.1 to 2 MeVincreasing the number of data points to 21868. The ex-perimental covariance matrix (never explicitly computed)contains two dense blocks, each larger than × ,due to the normalization uncertainty of the associateddatasets. The EXFOR identifications for these datasetsare 13511003 and 22316003.Clearly, the nuclear models in TALYS are not expectedto describe the resonances appearing in this energy range.The sole purpose of this exercise was to determine whetherthe algorithms in the pipeline can cope with this amountof data. All steps were performed without problems andthe total runtime of the pipeline not significantly increased.This can be explained by the fact that the runtime is al-most exclusively determined by the runtime of the nuclearmodels code TALYS. Each invocation of a TALYS calcula-tion lasted more than ten minutes whereas the statisticalalgorithms despite the increased number of data pointstook one order of magnitude less time thanks to the speedups explained in appendix B.31..Nuclear Data Evaluation Pipeline G. Schnabel et al. Figure 6. Posterior expectations and uncertainties of the evaluated cross section curves corresponding to the posterior of themodel parameters.
VI. SUMMARY AND CONCLUSIONS
We presented a prototype of a nuclear data evaluationpipeline which covers all aspects of the nuclear data evalu-ation process, starting from the retrieval of experimentaldata and their preparation over fitting a nuclear-modeland the generation of ENDF files. The pipeline is an ex-ample of a fully reproducible evaluation. It has been madeavailable online [14] as well as a Dockerfile [102] to facili-tate the installation of all required components as a singlebundle.Several innovations in evaluation methodology of thelast years found their way into the pipeline:• Rule-based construction and correction of experi-mental uncertainties• The introduction of extra systematic uncertaintiesof experimental data based on marginal likelihoodoptimization• Gaussian process priors on energy-dependent modelparameters• The exact treatment of the non-linear model byusing a modified Levenberg-Marquardt algorithm,which takes the experimental covariance matrix andprior parameter covariance matrix into account• An approximation of the posterior distribution that also incorporates information on second-orderderivatives of the nuclear modelThe latter four items benefited from a unified represen-tation of statistical errors, systematic errors and model pa-rameters in the form of a simple Bayesian network, whichenables to take into account a large collection of experi-mental datasets without the need of data reduction as apreparatory step. The ability to deal with a large amountof correlated experimental data implies the promise todeal with exact experimental covariance matrices on thescale of EXFOR. Both Monte Carlo and deterministicevaluation methods, such as GLS, are equally amenableto the advantages of this representation of experimentaldata.The pipeline was also discussed in terms of buildingblocks of information technology. We emphasized the needto make the following functionality available to creatorsof nuclear data evaluation pipelines:• Convenient programmatic access to experimentaldata• The convenient execution of model calculations inparallel in multi-processor environments or on sci-entific clusters and the retrieval of results• The convenient generation of model predictions forthe observables in the experimental datasets.32..Nuclear Data Evaluation Pipeline G. Schnabel et al.
We outlined design considerations of our implementationof these IT building blocks and provided short usage exam-ples to make the potential advantages for data treatmentand the creation of pipelines more tangible. In partic-ular, the conversion of the EXFOR library to a JSONdatabase using the MongoDB database engine greatly fa-cilitated the programmatic access to any information inthe EXFOR library. We anticipate that the ease of accessto any information in the EXFOR library will not onlyspeed up the creation of evaluation pipelines but makenuclear data more accessible to advanced statistical treat-ment and machine learning methods, such as exploredin [10, 33, 103–105].As another example, one important focus of the ITbuilding block implementing the parallelization of modelcalculations was the functioning on a variety of multi-processor or scientific cluster configurations, which is es-pecially relevant in the scientific environment where theconfiguration of scientific clusters can be quite differentacross different institutions.In the transformation or mapping of the output fromnuclear model codes to the observables measured in the ex-periments, we emphasized the need to have an extensiblecode architecture that enables the incremental addition ofsupport for new types of observables. We argued that func-tions responsible for the mapping should have access toall information—including meta information—associatedwith experimental data because any information may bepotentially relevant for a proper mapping. Therefore, func-tions for mapping always received the full information ofan EXFOR subentry in our implementation of this ITbuilding block.Several evaluation methods exist, either based on MonteCarlo sampling or deterministic update formulas, and ad-vancements in the field of Bayesian statistics and machinelearning will gradually also lead to improvements of meth-ods for nuclear data evaluation. We therefore stress thepoint that the value of the evaluation pipeline is not aboutthe specific choices of methods or data but as a demonstra-tion of an automated and reproducible evaluation involv-ing advanced statistical methods and parallelization on acomputer cluster. During the discussion of the pipeline, we also pointed out possible improvements, such as:• The introduction of computational fields in the Mon-goDB JSON database to facilitate the search andaccess of information in the EXFOR library• The automatic renormalization of experimental datain the EXFOR library according to the newest eval-uations of monitor reactions• The enhancement of the rule-based approach for theconstruction of uncertainty information by relyingon statistical data imputation and uncertainty tem-plates• The addition of outlier detection algorithms as acomplement to automatic uncertainty correction viathe MLO approach• The introduction of model defects on the observableside to deal with deficient models if the increasedflexibility achieved by energy-dependent model pa-rameters is not enoughBesides these potential improvements, there is anotherone deserving special attention in the future. Similar toEXFOR, existing evaluated data also represents historicknowledge which could be taken into account in the nu-clear data pipeline for the creation of a new evaluation.Rather than insisting on re-evaluating everything fromscratch, the release of new evaluations could be broughtforward in time. Hence, ease of access to existing evalu-ated data, similar to EXFOR in JSON format, is requiredto make such evaluation processes efficient. The GNDSformat is an important step in this direction.There are certainly many other possible improvements.The modular design of the presented evaluation pipelinewhich can be run at the push of a button enables rapidtesting of new algorithms and modified assumptions. Thescripts constituting the pipeline [14] as well as all support-ing packages have been made available online. A Dock-erfile [102] also available online greatly facilitates the in-stallation of all required components of the pipeline as asingle bundle. [1] F. Fidler and J. Wilcox, “Reproducibility of scientific re-sults,” in
The Stanford Encyclopedia of Philosophy (E. N.Zalta, ed.), Metaphysics Research Lab, Stanford Univer-sity, winter 2018 ed., 2018.[2] A. Carlson, V. Pronyaev, R. Capote, G. Hale, Z.-P. Chen,I. Duran, F.-J. Hambsch, S. Kunieda, W. Mannhart,B. Marcinkevicius, R. Nelson, D. Neudecker, G. Noguere,M. Paris, S. Simakov, P. Schillebeeckx, D. Smith, X. Tao,A. Trkov, A. Wallner, and W. Wang, “Evaluation of theNeutron Data Standards,”
Nuclear Data Sheets , vol. 148,pp. 143–188, Feb. 2018.[3] A. Carlson, V. Pronyaev, R. Capote, G. Hale, Z.-P. Chen,I. Duran, F.-J. Hambsch, S. Kunieda, W. Mannhart,B. Marcinkevicius, R. Nelson, D. Neudecker, G. Noguere, M. Paris, S. Simakov, P. Schillebeeckx, D. Smith, X. Tao,A. Trkov, A. Wallner, and W. Wang, “Corrigendum to“Evaluation of the Neutron Data Standards” [Nucl. DataSheets 148, p. 143 (2018)],”
Nuclear Data Sheets , vol. 163,pp. 280–281, Jan. 2020.[4] A. Carlson, V. Pronyaev, D. Smith, N. Larson, Z. Chen,G. Hale, F.-J. Hambsch, E. Gai, S.-Y. Oh, S. Badikov,T. Kawano, H. Hofmann, H. Vonach, and S. Tagesen,“International Evaluation of Neutron Cross Section Stan-dards,”
Nuclear Data Sheets , vol. 110, pp. 3215–3324,Dec. 2009.[5] V. Zerkin and B. Pritychenko, “The experimental nu-clear reaction data (EXFOR): Extended computerdatabase and Web retrieval system,”
Nuclear Instru- et al. ments and Methods in Physics Research Section A: Accel-erators, Spectrometers, Detectors and Associated Equip-ment , vol. 888, pp. 31–43, Apr. 2018.[6] N. Otuka, E. Dupont, V. Semkova, B. Pritychenko,A. Blokhin, M. Aikawa, S. Babykina, M. Bossant,G. Chen, S. Dunaeva, R. Forrest, T. Fukahori, N. Fu-rutachi, S. Ganesan, Z. Ge, O. Gritzay, M. Herman,S. Hlavač, K. Kat¯o, B. Lalremruata, Y. Lee, A. Mak-inaga, K. Matsumoto, M. Mikhaylyukova, G. Pikulina,V. Pronyaev, A. Saxena, O. Schwerer, S. Simakov, N. Sop-pera, R. Suzuki, S. Takács, X. Tao, S. Taova, F. Tárkányi,V. Varlamov, J. Wang, S. Yang, V. Zerkin, and Y. Zhuang,“Towards a More Complete and Accurate ExperimentalNuclear Reaction Data Library (EXFOR): InternationalCollaboration Between Nuclear Reaction Data Centres(NRDC),”
Nuclear Data Sheets , vol. 120, pp. 272–276,June 2014.[7] G. Schnabel, “A computational EXFOR database,” arXiv:1908.00209 [nucl-ex, physics:nucl-th,physics:physics] , May 2020. arXiv: 1908.00209.[8] A. Koning, D. Rochman, J.-C. Sublet, N. Dzysiuk,M. Fleming, and S. van der Marck, “TENDL: CompleteNuclear Data Library for Innovative Nuclear Science andTechnology,”
Nuclear Data Sheets , vol. 155, pp. 1–55, Jan.2019.[9] P. Helgesson and H. Sjöstrand, “Treating model defectsby fitting smoothly varying model parameters: Energydependence in nuclear data evaluation,”
Annals of Nu-clear Energy , vol. 120, pp. 35–47, Oct. 2018.[10] G. Schnabel, “Fitting and Analysis Technique for In-consistent Nuclear Data,” arXiv:1803.00960 [nucl-ex,physics:nucl-th, physics:physics] , Mar. 2018. arXiv:1803.00960.[11] P. Helgesson and H. Sjöstrand, “Fitting a defect non-linear model with or without prior, distinguishing nuclearreaction products as an example,”
Review of ScientificInstruments , vol. 88, p. 115114, Nov. 2017.[12] J. Pearl,
Causality: models, reasoning, and inference .Cambridge, U.K. ; New York: Cambridge UniversityPress, 2000.[13] A. Koning and D. Rochman, “Modern Nuclear Data Eval-uation with the TALYS Code System,”
Nuclear DataSheets , vol. 113, pp. 2841–2934, Dec. 2012.[14]
Prototype of evaluation pipeline .[15] G. Schnabel,
Dockerfile to set up EXFOR MongoDBdatabase . See appendix D.[16] G. Schnabel,
R package: nucdataBaynet . See appendix D.[17] D. W. Muir, “Evaluation of correlated data using par-titioned least squares: a minimum-variance derivation,”
Nuclear Science and Engineering , vol. 101, no. 1, pp. 88–93, 1989.[18] D. L. Smith, “A least-squares computational "tool kit",”Tech. Rep. ANL/NDM-128, Argonne National Labora-tory, Argonne, Illinois, Apr. 1993.[19] D. M. Hetrick and C. Y. Fu, “GLUCS: A generalized least-squares program for updating cross section evaluationswith correlated data sets.” Unknown, Oct. 1980.[20] D. W. Muir, A. Trkov, I. Kodeli, R. Capote, andV. Zerkin, “The Global Assessment of Nuclear Data,GANDR,” EDP Sciences, 2007.[21] D. Muir, “Global Assessment of Nuclear Data Require-ments (GANDR project),” 2007.[22] H. Leeb, D. Neudecker, and T. Srdinko, “Consistent Pro-cedure for Nuclear Data Evaluation Based on Modeling,”
Nuclear Data Sheets , vol. 109, pp. 2762–2767, Dec. 2008.[23] D. Neudecker, R. Capote, and H. Leeb, “Impact of modeldefect and experimental uncertainties on evaluated out-put,”
Nuclear Instruments and Methods in Physics Re-search Section A: Accelerators, Spectrometers, Detectorsand Associated Equipment , vol. 723, pp. 163–172, Sept.2013.[24] G. Schnabel and H. Leeb, “Differential Cross Sections andthe Impact of Model Defects in Nuclear Data Evaluation,”
EPJ Web of Conferences , vol. 111, p. 09001, 2016.[25] C. E. Rasmussen and C. K. I. Williams,
Gaussian Pro-cesses for Machine Learning . Cambridge, Mass.: MITPress, 2006.[26] R. H. Byrd, P. Lu, J. Nocedal, and C. Zhu, “A LimitedMemory Algorithm for Bound Constrained Optimiza-tion,”
SIAM Journal on Scientific Computing , vol. 16,pp. 1190–1208, Sept. 1995.[27] N. M. Larson, “A Concise Method for Storing and Com-municating the Data Covariance Matrix,” Tech. Rep.ORNL/TM-2008/104, 941045, Oak Ridge National Lab-oratory, Oct. 2008.[28] G. Schnabel,
Large scale Bayesian nuclear data evalua-tion with consistent model defects . PhD thesis, Technis-che Universität Wien, Vienna, 2015.[29] M. C. Kennedy and A. O’Hagan, “Bayesian calibrationof computer models,”
Journal of the Royal Statistical So-ciety: Series B (Statistical Methodology) , vol. 63, pp. 425–464, Aug. 2001.[30] D. Duvenaud,
Automatic model construction with Gaus-sian processes . PhD thesis, University of Cambridge,2014.[31] G. Schnabel and H. Sjöstrand, “A first sketch: Con-struction of model defect priors inspired by dynamictime warping,” arXiv:1811.03874 [nucl-ex, physics:nucl-th, physics:physics] , Nov. 2018. arXiv: 1811.03874.[32] E. Snelson and Z. Ghahramani, “Sparse Gaussian pro-cesses using pseudo-inputs,” in
Advances in neural infor-mation processing systems , pp. 1257–1264, 2006.[33] G. Schnabel, “Estimating model bias over the completenuclide chart with sparse Gaussian processes at the ex-ample of INCL/ABLA and double-differential neutronspectra,”
EPJ Nuclear Sciences & Technologies , vol. 4,p. 33, 2018.[34] J. Quiñonero-Candela and C. E. Rasmussen, “A unifyingview of sparse approximate Gaussian process regression,”
Journal of Machine Learning Research , vol. 6, no. Dec,pp. 1939–1959, 2005.[35] C. D. S. Jean, B. Habert, P. Archier, G. Noguere,D. Bernard, J. Tommasi, and P. Blaise, “Uncertainty Eval-uation of Nuclear Reaction Model Parameters using Inte-gral and Microscopic Measurements with the CONRADCode,”
Journal of the Korean Physical Society , vol. 59,pp. 1276–1279, Aug. 2011.[36] P. Archier, C. De Saint Jean, O. Litaize, G. Noguère,L. Berge, E. Privas, and P. Tamagno, “CONRAD Evalu-ation Code: Development Status and Perspectives,”
Nu-clear Data Sheets , vol. 118, pp. 488–490, Apr. 2014.[37] N. M. Larson, “Updated user’s guide for SAMMY: Multi-level R-Matrix Fits to Neutron Data Using Bayes’ Equa-tion, Revision 4,” Tech. Rep. ORNL/TM-9179/R8, OakRidge National Laboratory, Oak Ridge, 2008.[38] K. Levenberg, “A method for the solution of certain non-linear problems in least squares,”
Quarterly of AppliedMathematics , vol. 2, pp. 164–168, July 1944. et al. [39] D. W. Marquardt, “An Algorithm for Least-Squares Esti-mation of Nonlinear Parameters,”
Journal of the Societyfor Industrial and Applied Mathematics , vol. 11, pp. 431–441, June 1963.[40] M. A. Woodbury,
Inverting Modified Matrices . Statis-tical Research Group, Memo. Rep. no. 42, PrincetonUniversity, Princeton, N. J., 1950.[41] K. Madsen, H. Nielsen, and O. Tingleff, “Methods fornon-linear least squares problems (2nd ed.),” p. 60, 012004.[42] R Core Team,
R: A Language and Environment for Sta-tistical Computing . R Foundation for Statistical Com-puting, Vienna, Austria, 2017.[43] M. Dowle and A. Srinivasan, data.table: Extension of‘data.frame‘ , 2019. R package version 1.12.8.[44] D. Bates and M. Maechler,
Matrix: Sparse and DenseMatrix Classes and Methods , 2019. R package version1.2-18.[45] O. Schwerer, “EXFOR Formats Manual,” Tech. Rep.IAEA-NDS-207, IAEA, Vienna, Aug. 2015.[46] D. Cullen and A. Trkov, “PROGRAM X4TOC4 (Ver-sion 2001-3) Translation of Experimental Data from theEXFOR Format to a Computation Format,” Tech. Rep.IAEA-NDS-80, IAEA, Vienna, Mar. 2001.[47] A. Fedynitch, “afedynitch/x4i3,” May 2020.[48] G. Schnabel,
R package: exfor Parser . See appendix D.[49] G. Schnabel,
Script to add EXFOR masterfiles to Mon-goDB database . See appendix D.[50] “regular expression | Encyclopedia.com.”[51] G. Schnabel,
Dockerfile to set up EXFOR CouchDBdatabase . See appendix D.[52] J. Ooms, “The jsonlite package: A practical and con-sistent mapping between json data and r objects,” arXiv:1403.2805 [stat.CO] , 2014.[53] G. Schnabel,
R package: MongoEXFOR . See ap-pendix D.[54] H. Wickham, stringr: Simple, Consistent Wrappers forCommon String Operations , 2019. R package version1.4.0.[55] M. Herman, R. Capote, B. Carlson, P. Obložinský, M. Sin,A. Trkov, H. Wienke, and V. Zerkin, “EMPIRE: NuclearReaction Model Code System for Data Evaluation,”
Nu-clear Data Sheets , vol. 108, pp. 2655–2715, Dec. 2007.[56] M. Herman, R. Capote, M. Sin, A. Trkov, B. V. Carlson,P. Oblozinsky, C. Mattoon, H. Wienke, S. Hoblit, Y. Cho,G. Nobre, V. Plujko, and V. Zerkin, “EMPIRE-3.2 Malta-Modular system for nuclear reaction calculations andnuclear data evaluation,” Tech. Rep. INDC (NDS)-0603,IAEA, Vienna, Aug. 2013.[57] G. Schnabel,
R package: talysExforMapping . See ap-pendix D.[58] D. L. Smith, “A unified Monte Carlo approach to fastneutron cross section data evaluation,”
Proc. of the 8 thInternational Topical Mtg. on Nucl. Applics. and Util.of Accelerators, Pocatello, July , p. 736, 2008.[59] R. Capote and D. L. Smith, “An Investigation of the Per-formance of the Unified Monte Carlo Method of NeutronCross Section Data Evaluation,”
Nuclear Data Sheets ,vol. 109, pp. 2768–2773, Dec. 2008.[60] R. Capote, D. L. Smith, A. Trkov, and M. Meghzifene, “ANew Formulation of the Unified Monte Carlo Approach(UMC-B) and Cross-Section Evaluation for the Dosime-try Reaction Mn(n,g) Mn,”
Journal of ASTM In-ternational , vol. 9, pp. 179–196, Mar. 2012. [61] E. Bauge, S. Hilaire, and P. Dossantos-Uzarralde, “Eval-uation of the Covariance Matrix of Neutronic Cross Sec-tions with the Backward-Forward Monte Carlo Method,”EDP Sciences, 2007.[62] E. Bauge and P. Dossantos-Uzarralde, “Evaluation of theCovariance Matrix of
Pu Neutronic Cross Sectionsin the Continuum Using the Backward-Forward Monte-Carlo Method,”
Journal of the Korean Physical Society ,vol. 59, p. 1218, Aug. 2011.[63] A. J. Koning, “Bayesian Monte Carlo method for nu-clear data evaluation,”
The European Physical JournalA , vol. 51, p. 184, Dec. 2015.[64] T. Kawano and K. Shibata, “Covariance evaluation sys-tem,” Tech. Rep. JAERI-Data/Code 97-037, JAERI,Japan, 1997.[65] R. Capote, D. Smith, and A. Trkov, “Nuclear data eval-uation methodology including estimates of covariances,”
EPJ Web of Conferences , vol. 8, p. 04001, 2010.[66] “Covariance data in the fast neutron region,” Tech. Rep.NEA/NSC/WPEC/DOC(2010)427, OECD-NEA, Oct.2011.[67] P. Helgesson, D. Neudecker, H. Sjöstrand, M. Grosskopf,D. L. Smith, and R. Capote, “Assessment of Novel Tech-niques for Nuclear Data Evaluation,” in
Reactor Dosime-try: 16th International Symposium (M. H. Sparks, K. R.DePriest, and D. W. Vehar, eds.), pp. 105–116, 100 BarrHarbor Drive, PO Box C700, West Conshohocken, PA19428-2959: ASTM International, Nov. 2018.[68] W. Gentzsch, “Sun grid engine: towards creating a com-pute power grid,”
Proceedings First IEEE/ACM Inter-national Symposium on Cluster Computing and the Grid ,pp. 35–36, 2001.[69] J. Dean and S. Ghemawat, “MapReduce: simplified dataprocessing on large clusters,”
Communications of theACM , vol. 51, pp. 107–113, Jan. 2008.[70] M. Zaharia, R. S. Xin, P. Wendell, T. Das, M. Armbrust,A. Dave, X. Meng, J. Rosen, S. Venkataraman, M. J.Franklin, A. Ghodsi, J. Gonzalez, S. Shenker, and I. Sto-ica, “Apache Spark: a unified engine for big data process-ing,”
Communications of the ACM , vol. 59, pp. 56–65,Oct. 2016.[71] G. Schnabel,
R package: interactiveSSH . See appendix D.[72] G. Schnabel,
R package: rsyncFacility . See appendix D.[73] G. Schnabel,
R package: remoteFunctionSSH . See ap-pendix D.[74] G. Schnabel,
R package: clusterSSH . See appendix D.[75] G. Schnabel,
R package: clusterTALYS . See appendix D.[76] J. Ooms, ssh: Secure Shell (SSH) Client for R , 2019. Rpackage version 0.6.[77] RStudio Team,
RStudio: Integrated Development Envi-ronment for R . RStudio, PBC., Boston, MA, 2020.[78] A. Trkov, M. Herman, and D. Brown, “ENDF-6 For-mats Manual,” Tech. Rep. BNL-203218-2018-INRE,Brookhaven National Laboratory, Upton, NY 11973-5000, Feb. 2018.[79] “Specifications for the Generalised Nuclear DatabaseStructure,” Tech. Rep. NEA No. 7519, OECD NuclearEnergy Agency, 2020.[80] H. Sjöstrand and G. Schnabel, “Monte Carlo integral ad-justment of nuclear data libraries – experimental covari-ances and inconsistent data,”
EPJ Web of Conferences ,vol. 211, p. 07007, 2019.[81] D. Siefman, M. Hursin, H. Sjostrand, G. Schnabel,D. Rochman, and A. Pautz, “Data assimilation of post- et al. irradiation examination data for fission yields from GEF,”
EPJ Nuclear Sciences & Technologies , vol. 6, p. 52, 2020.[82] D. Smith and N. Otuka, “Experimental Nuclear ReactionData Uncertainties: Basic Concepts and Documentation,”
Nuclear Data Sheets , vol. 113, pp. 3006–3053, Dec. 2012.[83] H. Iwamoto, “Generation of nuclear data using Gaus-sian process regression,”
Journal of Nuclear Science andTechnology , vol. 57, pp. 932–938, Aug. 2020.[84] G. Schnabel,
R package: exforUncertainty . See ap-pendix D.[85] R. Peelle, “Peelle’s pertinent puzzle,”
Informal memoran-dum dated October , vol. 13, 1987.[86] S. Chiba and D. Smith, “A suggested procedure for resolv-ing an anomaly in least-squares data analysis known as“Peelle‘s Pertinent Puzzle“ and the general implicationsfor nuclear data evaluation,” Tech. Rep. ANL/NDM–121,10121367, Argonne National Lab., IL (USA), Sept. 1991.[87] P. Helgesson, H. Sjöstrand, A. Koning, J. Rydén,D. Rochman, E. Alhassan, and S. Pomp, “CombiningTotal Monte Carlo and Unified Monte Carlo: Bayesiannuclear data uncertainty quantification from auto-generated experimental covariances,”
Progress in NuclearEnergy , vol. 96, pp. 76–96, Apr. 2017.[88] D. Neudecker, B. Hejnal, F. Tovesson, M. White,D. Smith, D. Vaughan, and R. Capote, “Template forestimating uncertainties of measured neutron-inducedfission cross-sections,”
EPJ Nuclear Sciences & Tech-nologies , vol. 4, p. 21, 01 2018.[89] D. Neudecker, M. C. White, D. E. Vaughan, and G. Srini-vasan, “Validating nuclear data uncertainties obtainedfrom a statistical analysis of experimental data with the“Physical Uncertainty Bounds” method,”
EPJ NuclearSciences & Technologies , vol. 6, p. 19, 2020.[90] R. Capote, S. Badikov, A. Carlson, I. Duran, F. Gunsing,D. Neudecker, V. Pronyaev, P. Schillebeeckx, G. Schn-abel, D. Smith, and A. Wallner, “Unrecognized Sourcesof Uncertainties (USU) in Experimental Nuclear Data,”
Nuclear Data Sheets , vol. 163, pp. 191–227, Jan. 2020.[91] R. Forrest and J. Kopecky, “Statistical analysis of crosssections—A new tool for data validation,”
Fusion Engi-neering and Design , vol. 82, pp. 73–90, Jan. 2007.[92] R. Forrest, J. Kopecky, and A. Koning, “Detailed analysisof (n,p) and (n,alpha) cross sections in the EAF-2007and TALYS-generated libraries,”
Fusion Engineering andDesign , vol. 83, pp. 634–643, May 2008.[93] M. T. Pigni and H. Leeb, “Uncertainty Estimates of Eval-uated Fe Cross Sections Based on Extensive Modellingat Energies Beyond 20 MeV,” in
Proc. Int. Workshop onNuclear Data for the Transmutation of Nuclear Waste.GSI-Darmstadt, Germany , 2003.[94] H. Leeb, “Covariances for Evaluations based on Exten-sive Modelling,” in
AIP Conference Proceedings , vol. 769,(Santa Fe, New Mexico (USA)), pp. 161–164, AIP, 2005.ISSN: 0094243X.[95] B. J. N. Blight and L. Ott, “A Bayesian Ap-proach to Model Inadequacy for Polynomial Regression,”
Biometrika , vol. 62, p. 79, Apr. 1975.[96] A. B. Owen,
Monte Carlo theory, methods and examples .2013.[97] S. Brooks, A. Gelman, G. L. Jones, and X.-L. Meng, eds.,
Handbook for Markov chain Monte Carlo . Boca Raton:Taylor & Francis, 2011.[98] A. J. Koning and D. Rochman, “Towards Sustainable Nu-clear Energy: Putting Nuclear Physics to Work,”
Annals of Nuclear Energy , vol. 35, pp. 2024–2030, Nov. 2008.[99] G. Schnabel,
Patch for TASMAN . See appendix D.[100] H. Robbins, “An empirical bayes approach to statis-tics,” in
Proceedings of the Third Berkeley Symposium onMathematical Statistics and Probability, Volume 1: Con-tributions to the Theory of Statistics , (Berkeley, Calif.),pp. 157–163, University of California Press, 1956.[101] J. S. Maritz,
Empirical Bayes methods . Place of publica-tion not identified: Routledge, 2018. OCLC: 1036741238.[102] G. Schnabel,
Dockerfile to set up evaluation pipeline . Seeappendix D.[103] J. A. Hirdt and D. A. Brown, “Data mining the EXFORdatabase using network theory,” arXiv:1312.6200 [nucl-th, physics:physics] , Dec. 2013. arXiv: 1312.6200.[104] N. R. Dwivedi, “Trees and islands – machine learningapproach to nuclear physics,” 2019.[105] B. Whewell, M. Grosskopf, and D. Neudecker, “Evaluat-ing 239pu(n,f) cross sections via machine learning usingexperimental data, covariances, and measurement fea-tures,”
Nuclear Instruments and Methods in Physics Re-search Section A: Accelerators, Spectrometers, Detectorsand Associated Equipment , vol. 978, p. 164305, 2020.[106] D. A. Harville,
Matrix Algebra From a Statistician’s Per-spective . New York, NY: Springer New York, 1997.[107] A. George and J. W. Liu, “The Evolution of the Mini-mum Degree Ordering Algorithm,”
SIAM Review , vol. 31,pp. 1–19, Mar. 1989.[108] T. A. Davis,
Direct methods for sparse linear systems .Fundamentals of algorithms, Philadelphia: Society forIndustrial and Applied Mathematics, 2006.[109] P. Labs, “IPFS Powers the Distributed Web.” LibraryCatalog: ipfs.io.[110] M. Stevens, E. Bursztein, P. Karpman, A. Albertini,and Y. Markov, “The First Collision for Full SHA-1,”in
Advances in Cryptology – CRYPTO 2017 (J. Katzand H. Shacham, eds.), vol. 10401, pp. 570–596, Cham:Springer International Publishing, 2017. Series Title:Lecture Notes in Computer Science.[111] H. Wickham, P. Danenberg, G. Csárdi, and M. Eugster, roxygen2: In-Line Documentation for R , 2020. R packageversion 7.1.1. et al.
Appendix A: GENERAL STATISTICS1. Properties of variance and covariance operator
Especially when working with linear combinations ofrandom vectors governed by multivariate normal distribu-tion, using the properties of the variance and covarianceoperator is very convenient to derive distribution param-eters.Given random vectors (cid:126)x and (cid:126)y , we denote by
Cov [ (cid:126)x, (cid:126)y ] the matrix containing the covariances between elementsof (cid:126)x and (cid:126)y , i.e., (Cov [ (cid:126)x, (cid:126)y ]) ij = Cov [ x i , y j ] .With α being a real scalar and S a matrix, the followingproperties hold: Cov [ (cid:126)x, (cid:126)y ] = (Cov [ (cid:126)y, (cid:126)x ]) T (A1) Cov [ α (cid:126)x, (cid:126)y ] = α Cov [ (cid:126)x, (cid:126)y ] (A2) Cov [ S (cid:126)x, (cid:126)y ] = S Cov [ (cid:126)x, (cid:126)y ] (A3) Cov [ (cid:126)x + (cid:126)y, (cid:126)z ] = Cov [ (cid:126)x, (cid:126)z ] + Cov [ (cid:126)y, (cid:126)z ] (A4)The variance operator is given by Var [ (cid:126)x ] := Cov [ (cid:126)x, (cid:126)x ] .Therefore it always yields a symmetric matrix and has inaddition the properties Var [ α (cid:126)x ] = α Var [ (cid:126)x ] , (A5) Var [ S (cid:126)x ] = S Var [ (cid:126)x ] S T , (A6) Var [ (cid:126)x + (cid:126)y ] = Var [ (cid:126)x ] + Var [ (cid:126)y ] + 2 Cov [ (cid:126)x, (cid:126)y ] . (A7)For independent random variables, the last term contain-ing the covariances between (cid:126)x and (cid:126)y vanishes. Therefore,the variance of a sum of independent random variables isgiven by Var (cid:34) N (cid:88) i =1 (cid:126)x i (cid:35) = N (cid:88) i =1 Var [ (cid:126)x i ] . (A8) Appendix B: LINEAR ALGEBRA IDENTITIES TOSPEED UP COMPUTATIONS1. Derivative of an inverse matrix
Assume a matrix M ( (cid:126)κ ) which is a function of poten-tially several parameters κ l summarized in (cid:126)κ . We denoteits inverse by ˜M ( (cid:126)κ ) . The relation between these two ma-trices in terms of their components is given by (cid:88) j M ij ( (cid:126)κ ) ˜M jk ( (cid:126)κ ) = δ ij , (B1)with δ ij being one if i = j and zero otherwise. Taking thepartial derivative with respect to an element κ l of (cid:126)κ gives (cid:88) j (cid:32) ∂ M ij ( (cid:126)κ ) ∂κ l ˜M jk ( (cid:126)κ ) + M ij ( (cid:126)κ ) ∂ ˜M jk ( (cid:126)κ ) ∂κ l (cid:33) = 0 . (B2)This relation can be expressed in terms of matrix products, M ( (cid:126)κ ) ∂ ˜M ( (cid:126)κ ) ∂κ l = − ∂ M ( (cid:126)κ ) ∂κ l ˜M ( (cid:126)κ ) . (B3)Multiplying by ˜M ( (cid:126)κ ) from the left yields the final result ∂ ˜ M ( (cid:126)κ ) ∂κ l = − ˜ M ( (cid:126)κ ) ∂ M ( (cid:126)κ ) ∂κ l ˜ M ( (cid:126)κ ) . (B4)
2. Derivative of a logarithmized determinant
In the following we assume that a matrix M is a func-tion of potentially several parameters κ , κ , . . . summa-rized in a vector (cid:126)κ . In the pipeline, (cid:126)κ can contain system-atic and statistical uncertainties and hyperparameters ofGaussian processes. The logarithmized probability densityfunction of a multivariate normal distribution contains thelogarithmized determinant of a covariance matrix. Theavailability of the gradient in analytic form of this loga-rithmized determinant with respect to a set of parameters (cid:126)κ is extremely useful in gradient based optimization algo-rithms and is indeed computable.Using the chain rule, the derivate of ln det M ( (cid:126)κ ) can bewritten as ∂ ln det M ( (cid:126)κ ) ∂κ l = 1det M ( (cid:126)κ ) ∂ det M ( (cid:126)κ ) ∂κ l . (B5)Jacobi’s formula [106, p. 305] provides us with the deriva-tive of the determinant, ∂ det M ( (cid:126)κ ) ∂κ l = Tr (cid:18) adj (cid:0) M ( (cid:126)κ ) (cid:1) ∂ M ( (cid:126)κ ) ∂κ l (cid:19) . (B6)The adjugate matrix appearing in this expression is de-fined by [106, p. 192] M adj ( M ) = det( M ) ⇒ adj ( M ) = det( M ) M − (B7)Inserting eq. (B7) into eq. (B6) and the resulting expres-sion into eq. (B5) yields the final result: ∂ ln det M ( (cid:126)κ ) ∂κ l = Tr (cid:18) ( M ( (cid:126)κ )) − ∂ M ( (cid:126)κ ) ∂κ l (cid:19) . (B8)Please note that only the diagonal elements of the matrixproduct have to be computed to evaluate the trace. Fora matrix of size N × N , this reduces the time complexityfrom N to N .
3. Matrix determinant lemma
For the derivation of the matrix determinant lemma inthe version used in this paper, note that det (cid:0) A + UV T (cid:1) = det (cid:0) A (cid:0) + A − UV T (cid:1)(cid:1) = det( A ) det (cid:0) + A − UV T (cid:1) . et al. (B9)The application of Sylvester’s determinant identity [106,p. 416], i.e. det( + AB ) = det( + BA ) , yields det (cid:0) A + UV T (cid:1) = det( A ) det (cid:0) + V T A − U (cid:1) . (B10)Now replace U by the matrix product UW and extract W to obtain det (cid:0) A + UWV T (cid:1) =det( A ) det (cid:0)(cid:0) W − + V T A − U (cid:1) W (cid:1) . (B11)Making use of det( AB ) = det( A ) det( B ) , we get the finalresult det (cid:0) A + UWV T (cid:1) =det( A ) det( W ) det (cid:0) W − + V T A − U (cid:1) . (B12)
4. Woodbury identity
The Woodbury matrix identity [40], [106, p. 424] enablesthe inversion of a matrix of a specific structure in analternative way: ( A + UCV ) − = A − − A − U ( C − + VA − U ) − VA − (B13)with A , U , C , and V being matrices of appropriate dimen-sion.The validity of the formula can be verified by multiply-ing both sides by ( A + UCV ) from the left. As the leftside equals the identity matrix, it suffices to evaluate theright hand side of eq. (B13), ( A + UCV ) × [ A − − A − U ( C − + VA − U ) − VA − ]= + UCVA − − ( U + UCVA − U )( C − + VA − U ) − VA − = + UCVA − − UC ( C − + VA − U )( C − + VA − U ) − VA − = + UCVA − − UCVA − = . (B14)
5. Products involving the inverse of sparsecovariance or other positive-definite matrices
In section II B we had a matrix of the form Z = U − + S T D − S whose inverse needs to be multiplied with vectors or matrices. To carry out such products, inthe nucdataBaynet package we compute the Cholesky de-composition, e.g., [106], of covariance matrices and otherpositive-definite matrices first, Z = LL T , (B15)with L being a lower triagonal matrix.To evaluate a matrix product of the form R = XZ − X T , note that R = X ( LL T ) − X T = (cid:2) X ( L T ) − (cid:3) (cid:2) L − X T (cid:3) = MM T (B16)Therefore it suffices to evaluate M = X ( L T ) − and thencompute the final result by R = MM T . The computationof M can be formulated as the problem to solve the setof linear equations L T M = X , (B17)which can be efficiently done by back-substitution due to L T being upper triagonal.Also the computation of matrix products of the form S = XZ − can be restated as the solution to a set oflinear equations in two stages, L T M = X then LS = M , (B18)where the first equation can be solved by back-substitutionand the second one by forward-substitution.In addition to being positive-definite, we can expect Z mentioned at the beginning of this section to be sparsebecause:• The covariance matrix U associated with systematiccomponents is diagonal or at least block-diagonaland so is its inverse.• The sparseness of the mapping matrix S in com-bination with the diagonal structure of D leads to S T D − S being sparse, too.The Cholesky factor of sparse matrices is in general notsparse, hence the benefit of sparsity of Z would be lost forsolving the systems of linear equations outlined above.Therefore a reordering of rows of columns of the originalmatrix can be performed before computing the Choleskydecomposition so that resulting matrices are sparse. Apopular algorithm to determine the permutations are vari-ants of the minimum degree ordering algorithm [107, 108].For the operations with sparse matrices including matrixproducts, solving systems of linear equations with sparsematrices and sparse Cholesky decompositions, we rely onthe functionality of the R package Matrix [44].
Appendix C: BAYESIAN INFERENCE1. Details on modified Levenberg-Marquardtalgorithm
The LM algorithm [38, 39] extends the linear-leastsquares method to be applicable to non-linear models.38..Nuclear Data Evaluation Pipeline G. Schnabel et al.
We use a modification of the LM algorithm to account forthe experimental covariance matrix and prior parametercovariance matrix.The objective is to locate the maximum of the posteriorpdf π ( (cid:126)p | (cid:126)σ exp ) = 1 π ( (cid:126)σ exp ) (cid:96) ( (cid:126)σ exp | (cid:126)p ) π ( (cid:126)p ) . (C1)The likelihood and prior parameter distribution are givenby (cid:96) ( (cid:126)σ exp | (cid:126)p ) = N ( (cid:126)σ exp | M ( (cid:126)p ); Σ exp ) , (C2) π ( (cid:126)p ) = N ( (cid:126)p | (cid:126)p ; P ) . (C3)The model link M ( (cid:126)p ) is non-linear. The LM algorithm isan iterative approach and it relies in each iteration on alinear approximation of the non-linear model, M lin ( (cid:126)p ) = (cid:126)p ref + J ( (cid:126)p − (cid:126)p ref ) . (C4)The introduction of δ(cid:126)p = (cid:126)p − (cid:126)p ref , (cid:126)d = (cid:126)σ exp − (cid:126)σ ref . allowsus to write ( −
2) ln π ( (cid:126)p | (cid:126)σ exp ) = (cid:16) (cid:126)d − J δ(cid:126)p (cid:17) T Σ − exp (cid:16) (cid:126)d − J δ(cid:126)p (cid:17) + ( δ(cid:126)p + (cid:126)p ref − (cid:126)p ) T P − ( δ(cid:126)p + (cid:126)p ref − (cid:126)p ) + C (C5)with the constant C absorbing everything which is inde-pendent of δ(cid:126)p . The constant is of no significance as we aregoing to take derivatives with respect to elements in δ(cid:126)p .Now we isolate terms containing δ(cid:126)p and regroup theexpression according to their order, δ(cid:126)p T (cid:0) J T Σ − exp J + P − (cid:1) δ(cid:126)p − δ(cid:126)p T (cid:16) J T Σ − exp (cid:126)d + P − ( (cid:126)p − (cid:126)p ref ) (cid:17) + T (C6)where T denotes the transpose of the terms explicitlywritten out. To determine the vector δ(cid:126)p that minimizesthis expression, we calculate the gradient of eq. (C6) andrequire it to vanish: (cid:0) J T Σ − exp J + P − (cid:1) δ(cid:126)p − (cid:16) J T Σ − exp (cid:126)d + P − ( (cid:126)p − (cid:126)p ref ) (cid:17) = (cid:126) . (C7)Rearranging yields (cid:0) J T Σ − exp J + P − (cid:1) δ(cid:126)p = (cid:16) J T Σ − exp (cid:126)d + P − ( (cid:126)p − (cid:126)p ref ) (cid:17) . (C8)Please note that the right-hand side is proportional tothe gradient of eq. (C6) evaluated at the reference pa-rameter vector, i.e., δ(cid:126)p = 0 , used as expansion point forthe construction of the linear approximation. Therefore itis also proportional to the gradient of the logarithmizedposterior distribution ln π ( (cid:126)p | (cid:126)σ exp ) evaluated at (cid:126)p = (cid:126)p ref . We introduce the following abbreviations: A = J T Σ − exp J + P − + λ I (C9) (cid:126)b = J T Σ − exp ( (cid:126)σ exp − (cid:126)σ ref ) + P − ( (cid:126)p − (cid:126)p ref ) . (C10)The matrix I is the identity matrix. The introduction ofthe term λ I allows to make A more diagonal by increasingthe value of λ . The purpose of this mechanism will becomeapparent in a moment.The update equation for δ(cid:126)p can be written as A δ(cid:126)p = (cid:126)b (C11) (cid:126)p prop = (cid:126)p ref + δ(cid:126)p (C12)With increasing λ , the matrix A becomes more similar tothe identity matrix. Consequently, the update degeneratesgradually to the gradient descent method. For λ = 0 , theupdate is equivalent to the GLS method.
2. Prior on second-derivative as smoothness prior
Assume that we have tuples of cross values y i and as-sociated energies x i , i.e. { ( x i , y i ) } i =1 ..N , in some reactionchannel. Let the energies x i be sorted in ascending order,i.e. x i < x i +1 . We can use a polygon chain to obtain crosssections for any intermediate energies x , f ( x ) = N (cid:88) i =1 I i ( x ) (cid:0) u i ( x ) y i + v i +1 ( x ) y i +1 (cid:1) (C13)with I i ( x ) = 1 if x i < x < x i +1 , and otherwise. Wehave introduced the following abbreviatios which will beuseful later on: u i ( x ) = x i +1 − xx i +1 − x i and v i ( x ) = x − x i − x i − x i − (C14)Given that the grid spanned by the x i ’s is sufficientlydense, any continuous function can be approximated withsufficient precision.The prediction for cross sections ˜ y = (˜ y , · · · , ˜ y M ) atlocations ˜ x = (˜ x , · · · , ˜ x M ) T using eq. (C13) can be alsowritten in matrix notation. To that end, we introduce thesensitivity matrix S ji := (cid:40) u i (˜ x j ) if x i < ˜ x j < x i +1 v i (˜ x j ) if x i − < ˜ x j < x i . (C15)Equation (C13) now becomes ˜ y = Sy . (C16)In principle, we could introduce a prior on the valuesin the vector y . However, sometimes we would like toexpress our knowledge in terms of allowed differences or39..Nuclear Data Evaluation Pipeline G. Schnabel et al. even changes of differences. This can be achieved by areparametrization: y i = (cid:40) y if i = 1 y + (cid:80) i − k =1 ∆ y k if i ≥ (C17)We can summarize the new parameters in a vector u =( y , ∆ y , · · · , ∆ y N − ) T and define the following matrix T ij = (cid:40) if i ≥ j otherwise (C18)in order to identically rewrite eq. (C17) in matrix notationas y = T u . (C19)To get a parametrization in terms of second derivatives,we apply once again a similiar transformation, ∆ y i = (cid:40) ∆ y if i = 1∆ y + (cid:80) i − k =1 ∆ y k if i ≥ (C20)We can summarize the new parameters in a vector v =( y , ∆ y , ∆ y , · · · , ∆ y N − ) T and transform them to u defined above using the following matrix R ij = if i = j = 11 if i ≥ and i ≥ j otherwise (C21)in order to obtain u = Rv . (C22)The concatenation of all the matrices yields finally theprediction of the polygon chain in terms of the param-eters corresponding to the (finite version) of the secondderivatives. ˜ y = Zv with Z = ST R (C23)The prior can be formulated for the elements in v and itcan be assumed that no correlations between elements in v exist, π ( v ) = N ( v | (cid:126) , D ) . (C24)The vector of zeros indicates that the best guess for thecurvature of the cross section curve is zero, i.e., changesof the slope are penalized.Even though the prior covariance matrix D for the val-ues of the second-order derivative is diagonal, correlationsare induced between the function values on the mesh dueto the particular form of the mapping matrix Z .
3. Necessity of second-order Taylor approximationof posterior pdf
After the localisation of the posterior maximum by theLM algorithm, one may be tempted to apply the GLSformulas to compute the posterior covariance matrix asit only requires the Jacobian matrix of the nuclear model,i.e., first-order derivates of the nuclear model. Here weexplain why this approach is potentially flawed and alsosecond-order derivatives of the model need to be takeninto account.The logarithmized posterior pdf is given by ln π ( (cid:126)p | (cid:126)σ ) = −
12 ( (cid:126)σ − M ( (cid:126)p )) T Q ( (cid:126)σ − M ( (cid:126)p )) −
12 ( (cid:126)p − (cid:126)s ) T R ( (cid:126)p − (cid:126)s ) + C . (C25)The matrix Q is the inverse of the experimental covariancematrix Σ exp , the vector (cid:126)σ contains the measurements, thematrix R is the inverse of the prior parameter covariancematrix, the vector (cid:126)s denotes the prior expectation of (cid:126)p and M yields the vector with model predictions based onparameter set (cid:126)p . The argument of M will be dropped fromnow on for notational convenience. The term C summa-rizes the terms independent of (cid:126)p and ensures the propernormalization of the pdf.We rewrite the right hand side in element-wise notation: − N (cid:88) i =1 N (cid:88) j =1 ( σ i − M i ) Q ij ( σ j − M j ) − M (cid:88) m =1 M (cid:88) n =1 ( p m − s m ) R mn ( p n − s n ) + C . (C26)The upper limit N of the first two summations is thenumber of experimental data points and the upper limit M of the last two summations is the number of modelparameters.Taking the first derivative with respect to a model pa-rameter p k yields G k := N (cid:88) i =1 N (cid:88) j =1 ( ∂ k M i ) Q ij ( σ j − M j ) − M (cid:88) n =1 R kn ( p n − s n ) , (C27)where we used the short-hand notation ∂ k := ∂/∂p k . Wedenote this expression as G k because it is the k-th elementof the gradient of the logarithmized posterior pdf withrespect to the model parameters.We take once again the derivative with respect to a40..Nuclear Data Evaluation Pipeline G. Schnabel et al. model parameter p t H kt := N (cid:88) i =1 N (cid:88) j =1 ( ∂ k ∂ t M i ) Q ij ( σ j − M j ) − N (cid:88) i =1 N (cid:88) j =1 ( ∂ k M i ) Q ij ( ∂ t M j ) − R kt . (C28)This expression states the elements H kt of the Hes-sian matrix of the logarithmized posterior pdf. Werewrite eq. (C28) in matrix notation, H = − ( U + J T QJ + R ) , (C29)with the elements of the matrix U being of the form U kt = − N (cid:88) i =1 N (cid:88) j =1 ( ∂ k ∂ t M i ) Q ij ( σ j − M j ) (C30)and the elements of the Jacobian matrix J being J ik := ∂ k M i .For a linear model, the matrix U in eq. (C29) van-ishes and recalling that for a linear model the relationbetween the Hessian and the covariance matrix is givenby P = − H − , we recover the GLS solution for the pos-terior covariance matrix, P = ( J T QJ + R ) − .For non-linear models, the neglect of the first term con-taining the second-order derivative of the model can beproblematic for two reasons:• The posterior is not contracted enough so that thelinear approximation of the model is not sufficientin the domain of the parameter space admissible bythe posterior distribution.• The model is deficient in a statistical sense, i.e., theresiduals ( σ j − M j ) are large and/or not about asoften positive as negative, i.e., the model system-atically over- or underestimates the data. Even ifthe non-linearity of the model is mild, the first termmay get inflated due to the large residuals so thatits contribution to the Hessian cannot be neglectedanymore.By implication, the generalized least squares method hasa higher chance of producing a valid posterior covariancematrix for a non-linear model if the model is not deficient.The modelling of model defects by Gaussian processescould improve the situation for a nuclear model with de-ficiencies as the elements Q ij in the inverse experimentalcovariance matrix are decreased, hence also reducing themagnitude of the first term in eq. (C28).
4. Computation of the posterior covariance matrix
In order to reduce the computational requirement ofthe LM algorithm, one may skip the computation of the elements in the Jacobian matrix associated with modelparameters insensitive to the experimental data and setthese elements to zero.After the optimization by the LM algorithm, the poste-rior covariance matrix needs to be computed for all param-eters, including those insensitive to the experimental data.The reason being that correlations between parametersimposed by the prior lead to a propagation of informa-tion from the sensitive parameters to the insensitive ones.For instance, Gaussian processes introduce correlationsbetween energy-dependent optical model parameters, e.g.,between the values of r v at 5 MeV and 10 MeV. If r v at10 MeV gets excluded from the LM optimization, we havea situation where a sensitive parameter is already a prioricorrelated to an insensitive one.In this section we explain the calculation of the fullposterior covariance. For the derivation we will take intoaccount the specific structure of the Hessian matrix dueto the presence of insensitive parameters and introducean approximation to make the computation tractable forcomputationally expensive models, such as TALYS.To determine the posterior covariance matrix of all pa-rameters, we consider a second-order Taylor approxima-tion of the logarithmized posterior pdf constructed at thereference parameter vector (cid:126)p ref : ln π ( (cid:126)p | (cid:126)σ exp ) =ln π + G ( (cid:126)p − (cid:126)p ref ) + 12 ( (cid:126)p − (cid:126)p ref ) T H ( (cid:126)p − (cid:126)p ref ) . (C31)The normalization constant is given by ln π =ln π ( (cid:126)p ref | (cid:126)σ exp ) and G is the gradient and H the Hessianmatrix of the logarithmized posterior pdf evaluated at (cid:126)p ref . The form of the gradient was given in eq. (C27) andthat of the Hessian in eq. (C29). Given the nuclear modelcan be well approximated by a second-order Taylor expan-sion in the parameter domain admissible by the posteriordistribution, the posterior covariance matrix can be ap-proximated as P = − H − .We split the parameter vector in blocks (cid:126)p = ( (cid:126)p T , (cid:126)p T ) T with (cid:126)p containing the sensitive parameters and (cid:126)p con-taining the insensitive parameters. The Jacobian matrixis divided into blocks accordingly, i.e., J = ( J , ) . Insensi-tive parameters do not have any impact on the observablesunder consideration, which leads to the block of zeros inthe Jacobian.By excluding parameters from optimization by the LMalgorithm on the basis of the sensitivity matrix, we as-sume that the sensitivity to a specific parameter not onlyvanishes locally at the current expansion point but ev-erywhere, i.e., ∂/∂p k M i ( (cid:126)p ) = 0 for all predictions M i and any parameter vector (cid:126)p in the admissible parameterdomain. The validity of this assumption can be assessedafter the LM optimization by checking whether suppos-edly insensitive parameters are also insensitive accordingto the sensitivity matrix evaluated at the optimized pa-rameter set. Given this assumption holds, the matrix U et al. in eq. (C29) can be written in the partitioned form U = (cid:18) U
00 0 (cid:19) (C32)due to the vanishing partial derivatives ∂ k ∂ t M i if p k or p t is an insensitive parameter. The matrix J T QJ is of thesame structure due to the structure of J . Consequently,the Hessian matrix reads H = − (cid:18) U + R + J T QJ R R R (cid:19) (C33)and the full posterior covariance matrix can be computedas P = − H − . Noteworthy, the blocks of the inverseprior parameter covariance matrix are different from theinverses of the blocks in the the prior parameter covariancematrix, i.e., R ij (cid:54) = ( P ij ) − .The full computation of the matrix U may be too costlyfor a large number of parameters in combination with anexpensive nuclear model. A seemingly possible solution isto calculate only the diagonal elements of U and assumeall off-diagonal elements being zero. In a simulation study,we found that this construction may lead to H not beingpositive-definite. A solution that still takes into accountsecond-order derivatives of the model is to approximate U by a diagonal matrix ˜ U whose elements are given by ˜ U ii = (cid:40) U ii if U ii > otherwise (C34)Considering that the resulting approximative posteriorcovariance matrix is of the form, ˜ P = (cid:16) ˜ U + J T QJ + R (cid:17) − , (C35)we can see that ˜ U enters in the same way as the inverseprior covariance matrix R . Therefore ˜ U can be interpretedas an additional prior information that enforces an upperbound on the diagonal elements of the posterior covariancematrix. Even though this approximation may be regardedas less conservative than the GLS formula ( J T QJ + R ) − for the posterior covariance matrix, we think that it is abetter approximation to the full Hessian matrix becauseit takes more information of the full Hessian matrix intoaccount. Preliminary simulation studies seem to supportthis view. Moreover, as can be understood by eq. (C38)of appendix C 5, an overestimation of posterior uncertain-ties as can be the case in the GLS approach leads to atoo large adjustment of the parameters.
5. Computation of the posterior expectation ofinsensitive paramaters
In appendix C 4 we explained the computation of thefull posterior covariance matrix including parameters opti-mized by the LM algorithm and those considered insensi-tive to the data and therefore excluded from optimization. However, prior correlations not only change the posteriorcovariances associated with insensitive parameters butalso the posterior expectation values.To find the update prescription for insensitive param-eters, we first consider the update for all parameters. Tothis end, we need to locate the maximum of the posteriorpdf in eq. (C31). As a necessary condition, the gradientof this pdf for the parameter vector associated with themaximum must vanish, G + H ( (cid:126)p − (cid:126)p ref ) = (cid:126) . (C36)Thus the parameter vector maximizing eq. (C31) reads (cid:126)p = (cid:126)p ref − H − G . (C37)If the full Hessian is not available due to the infeasibil-ity to compute all second-order derivatives of the nuclearmodel, we can adopt the approximation suggested in ap-pendix C 4. Recalling that P = − H − and employingthe approximation of the posterior covariance matrix ineq. (C35), the update formula can be restated as (cid:126)p = (cid:126)p ref + ˜ P G (C38)We already computed G in eq. (C27) for a posteriorpdf with a multivariate normal prior and likelihood. Werestate it here in matrix notation, G = J T Q ( σ − σ ref ) − R ( (cid:126)p ref − (cid:126)s ) , (C39)with J being the Jacobian matrix of the nuclear model.This gradient is evaluated at (cid:126)p ref , where the values of sensi-tive parameters are from the solution of the LM algorithmand the values of insensitive parameters are given by theprior expectations.All ingredients of the update formula in eq. (C38) arefully specified now. The final parameter vector (cid:126)p containsupdated values of both sensitive and insensitive param-eters. We recommend to discard the new values of thesensitive parameters and keep the ones obtained by theLM algorithm. The reason being that the LM algorithmlocates the maximum of the posterior pdf exactly whereaseq. (C38) employs an approximation of the posterior co-variance matrix. As the insensitive parameters were ex-cluded from optimization, we have to rely on eq. (C38)using the approximative posterior covariance matrix toupdate the insensitive parameters. Appendix D: DOWNLOADING THE PIPELINEAND PACKAGES
The internet is a fast moving medium. Domain namesappear and vanish after a while, and with them all theircontent. Also content evolves and the content accessibleunder a uniform resource locator (URL) or colloquiallyweb address may change within months or years. More-over, in many contexts, such as legal documents and sci-entific publications and data, it is important to have some42..Nuclear Data Evaluation Pipeline G. Schnabel et al. kind of assurance that the documents have not been tam-pered with. An interesting and promising approach to pro-tect against the accidental loss of information and camou-flaged modifications of content is the
InterPlanetary FileSystem (IPFS) [109]. It is a distributed file system buildon top of computers accessible over the internet. Files arenot accessed by a name chosen by someone but a uniqueidentifier directly derived from the file content by using acryptographic hash function.We anticipate that the pipeline and employed packageswill evolve over time. We want to make sure that thepipeline in the state associated with this version of thepaper can be easily located. Inspired by the idea of contentaddressable storage, the underlying concept of the IPFS,we follow a similar approach and provide identificationstrings that allow the localisation of the pipeline and sup-port packages on the internet. First, we elaborate on theidentification strings and how they can be used to locatethe pipeline, packages and files on the internet. Then weprovide a catalogue of components including the pipelineitself with a brief description and the identification stringsto locate the components.
1. About identification strings
We use two types of identification strings: (1) git com-mit hashes and (2) SHA-256 hashes:1. We used
Git as version control system during the de-velopment of the pipeline and all the packages. Dur-ing the development process after a certain numberof modifications, these modifications are committed.Each commit records the exact state of the files inthe project at a given time. Moreoever, each com-mit is associated with a unique identifier, which isnot only unique within a development project butunique among all the development projects on theplanet tracked by git. At least, the chance of obtain-ing identical identifiers of different projects with dif-fering file contents by chance is extremely unlikely .The Git commit hash is therefore a reasonable iden-tifier for files and projects at a given time if securityis not the highest priority.2. In the future Git may get upgraded or we maychange to another version control system. Therefore,there is a need to be able to generate an identifier,which is independent of Git. For the files referencedin this paper, we obtain this identifier in the follow-ing way: Researchers at CWI Amsterdam and Google have demonstrated afeasible approach to generate two documents with identical SHA-1 hash [110]. The SHA-1 hash function is also used by Git. Theamount of compute power required to generate two documentswith identical SHA-1 hashes has still to be considered enormousat the time of writing. c a t <( f i n d . − type d \ − not − path ' ∗ / \ . ∗ ' \ − p r i n t 0 | \LC_ALL=C s o r t − z ) \<( f i n d . − type f \ − not − path ' ∗ / \ . ∗ ' \ − p r i n t 0 | \LC_ALL=C s o r t − z | \x a r g s − − We refer to this identifier as the tar hash . If securityis a major concern, then the SHA-256 hash functionto obtain the tar hash offers far greater security thanthe SHA-1 hash function. Please note that hiddenfiles are excluded to avoid the consideration of the .git folder and the .gitignore file.At the time of writing, there are several ways to locatea resource by its identifier:• Using either the Git commit hash or the tar hash,the resource can be located via http://nugget.link/
It is possible tosupply only a part of the full hash instead as long asthe sequence is unique among the registered hashes.• Using a Git commit hash, it can be located onGitHub via the link: https://github.com/search?q=
For example, the Dockerfile for setting up thepipeline can be located by following https://github.com/search?q=d16d70579707&type=Commits . In this ex-ample we only used the initial part of the identifier,which also works.• The tar hash and Git commit hash may have alsobeen added to the webpage associated with a projector package. A search query with the tar hash orGit commit hash on Google, Bing, DuckDuckGo orany other popular search engine may therefore alsoreturn the relevant resource. At the time of writing,however, we have only tested this with the Googlesearch.Links to the Dockerfiles of the pipeline and relatedpackages, such as the MongoDB JSON database can alsobe found on .
2. Catalogue of pipeline modules
In this section we present the list of packages in table IXemployed in the pipeline along with their identificationstrings and a brief description. Also the repository of thepipeline with the sequence of retrieval, preprocessing andfitting scripts is listed.We decided to assign the version number 0.1.0 to thepackages in the state at the time of writing of this article.43..Nuclear Data Evaluation Pipeline G. Schnabel et al.
Most packages come with a README file with a short ex-ample of usage. Function documentation in most packagesis available, created with the help of
Roxygen2 [111].Having said that, the documentation of many compo-nents does probably not meet the standards expected insoftware engineering regarding detail and comprehensive-ness at present. As another note of caution, as the pack-ages have been developed along with the pipeline, thepackages have not been thoroughly tested in a contextoutside of the current implementation of the pipeline. 44..Nuclear Data Evaluation Pipeline G. Schnabel et al. eval-fe56-docker
A Dockerfile to install the pipeline and dependencies as one bundletar hash 9ac4ab657ecf4887f302d835c45f5e690c998a8a1b0d9c769f19c66d7a41cf3egit commit id d16d70579707892b64d5cd972fed02492cc9742e eval-fe56
The sequence of script files constituting the pipelinetar hash 4aa395ec1a7a121cae9fab09e687450bc19034c11c9daba81bcc525acd4887ddgit commit id 58a51e018d402362ecdc8d0851d597d9f7890c9b nucdataBaynet
Management of experimental and model data in the Bayesian contexttar hash 9930f956d4da44bca735672a58464ce844c84f456551c6a7985509c46dbd60a0git commit id 284500f7d09b142ab5a468687612a93d8c51d85b
TALYSeval
Prepare TALYS calculations and retrieve resultstar hash 963bc700027d16f33c06bb9194bf95b1e1512b6eeddfc05190d1a3ce32ae4f88git commit id 857583ba7e259bf5a8bad9fbd558adc77bfe3579 tasmanPatch
Patch to extend functionality of the TASMAN codetar hash cfd0a09ecc0b9b135d0f1cd03c1d370bf8ccb55e620cc44282546e4b8a164ed2git commit id 05d243c36708356142dc1828590c8d1f62a9bde4 interactiveSSH
Executing bash commands remotely via SSH from within Rtar hash 4b8af95035326a42542cefb186ab0237c827a31f5286bbbe40c0ad933f0d4e3dgit commit id 08adb733de05eb9b4acb76f76e64b0d9cd5a032f rsyncFacility
Wrapper package around the rsync command line tooltar hash 5cb3d8e1e9d367eb9d50aea963f9958841eed0987707892314b8549cd899e5fcgit commit id fe3825864ee93c0b9fcfd6cf68400051836982bb remoteFunctionSSH
Execute R functions remotely in the same way as local functionstar hash 8a8e49572ad0e99e8c56149c63329a951766b7a0c232452f3006b4c39870a42egit commit id 681c4db9f85a712a406b3d49ae47fd235bd9e6b3 clusterSSH
Execute R functions in parallel on a computer clustertar hash 78d7ff6205cb96b25335dfb0a96dfc3bcbefed22b10c48809285e7f7ce750f79git commit id c8cc3434aec82f46c003d88c35861757f39dc5b5 clusterTALYS
Launch TALYS calculations in parallel on a computer clustertar hash c8fe3f0abbde33b60017034f980ad0b843afcb200ed86e04c9d51aa9cd4df1eagit commit id 144cecbbcac73f21086020ff3a69958c837b2d74
MongoEXFOR
Convenience interface to the EXFOR JSON MongoDB databasetar hash cf4c020ea89d3cc2fa10c1298c98c8c4260bdd0708e88501a0d7bcf94b907ed7git commit id 2e2dea55aea90ddca7851f264580bb3ad158c45f exforUncertainty
Tentative rules for rule-based correction of experimental uncertaintiestar hash c917f70ac4520a0aef2cca8479e8e28644c8ab3c1062d439a7af9fad2d25a38agit commit id 7b48fa675396e4f2290e6a9361fefc77d66d5ffa talysExforMapping
Mapping of TALYS predictions to observables recorded in EXFOR subentriestar hash 5ac8f35314ff5c411cf4409ac48339b35784ac8de250be1de246ed90c0362cefgit commit id 2ef3b109b66e0fbe103ec6a6cd2a42b828252926 jsonExforUtils
Utility functions to deal with information in the EXFOR entriestar hash e856b7734bd99b81ba02cd002236590b31094ef0d91eb811200c0dc22b2885d7git commit id b464876d312a7316549ed0f3f48dea13217be83d exforParser
Read EXFOR entries in R and conversion to JSONtar hash 09104fd60025c4d655d7fad1306a2afc2f049c6293136ed4e4a48652f80eba63git commit id c732c4e05824ac68f9905d603bce10777b2ce9b9 createExforDb
Creation of the EXFOR JSON MongoDB database from EXFOR masterfilestar hash 10c618b46ef402be1fc82a0e6583ea09fccae23003c23f99b1448fedc027c29fgit commit id af6d2af843cae280f37292858fb2a7c4fb50aeb4 compEXFOR-docker
Dockerfile for a stand-alone installation of the MongoDB EXFOR JSON databasetar hash e1bc125d69477777a328cf7a20b255cbee29748b187a4c24676e58aa47b6e8ddgit commit id 6bd8232157a82e5dce38150e58f257554937d386 exfor-couchdb-docker