key: cord-0619537-szlq62gw
authors: Laloy, Eric; Rogiers, Bart; Bielen, An; Boden, Sven
title: Bayesian inference of 1D activity profiles from segmented gamma scanning of a heterogeneous radioactive waste drum
date: 2021-01-06
journal: nan
DOI: nan
sha: f9cac80a7b3f2ca781e6aa05e2308f65264a342c
doc_id: 619537
cord_uid: szlq62gw

We present a Bayesian approach to probabilistically infer vertical activity profiles within a radioactive waste drum from segmented gamma scanning (SGS) measurements. Our approach resorts to Markov chain Monte Carlo (MCMC) sampling using the state-of-the-art Hamiltonian Monte Carlo (HMC) technique and accounts for two important sources of uncertainty: the measurement uncertainty and the uncertainty in the source distribution within the drum. In addition, our efficiency model simulates the contributions of all considered segments to each count measurement. Our approach is first demonstrated with a synthetic example, after which it is used to resolve the vertical activity distribution of 5 nuclides in a real waste package.

Segmented gamma scanning (SGS) is a longstanding but still actively researched technique to assess the activity of radioactive waste drums in a non-destructive way [e.g., 2, 15, 20, 8 , and references therein]. Classical SGS interpretation methods either sum up all gamma spectra to derive representative activities for the whole drum or treat the measured counts associated with each segment individually. Efficiency calibration can be performed in a classical way using reference point sources that are placed in a waste drum matrix in a special way to simulate source homogeneity. Nowadays, Monte Carlo (MC) modeling techniques such as implemented in the commercially available ISOCS software [ISOCS: In Situ Object Counting System, 16] are commonly used to determine detector's efficiency. Though offering significant geometrical flexibility, the ISOCS variants cannot always fully account for the spatial distribution of activities within the drum. Furthermore, the uncertainty in the estimated activities by radionuclide quantification techniques is often calculated in a relatively simple way, whether based on first-order Taylor expansion or MC error propagation, ignoring for instance drawbacks related to non-physical values, non-Gaussian distributions or non-linear expressions, or the omission of, e.g., systematic uncertainties like that of the efficiency [e.g., 14] . Bayesian inference, also called Bayesian data inversion, is well suited to assess the uncertainty of quantities that result from the application of a mathematical model [see the book by 10, for a comprehensive description of Bayesian data analysis]. In our opinion, Bayesian inference is particularly attractive for radiological characterization for the following reasons.

• We are not simplifying the mathematical model, as is the case for first-order Taylor expansion, and can work with any type of distribution (this is the same for MC error propagation)

• We can account for prior information and expert opinion. Also MC error propagation allows this. The difference is however, that if there is information in the data, the Bayesian inference will update those prior beliefs into posterior beliefs, while in the MC error propagation, there is no way to do that. Also, for the Bayesian inference, we can express our prior beliefs on the quantities of interest (e.g. activities), while for the MC error propagation there is no way to do that.

• We do in principle not need to give special attention to performance characteristics like decision threshold, detection limit and minimum detectable activities (MDA). The Bayesian inference will, based on our prior beliefs, always provide the possible range of outcomes, even for cases where a classic workflow would have resulted in a reported MDA.

• With this approach it is straightforward to deal with situations for which multiple observations inform about the same quantities and/or processes. For instance, when one or more inferred radionuclide activities are associated with multiple energy peaks in gamma spectrometry. Here for a given nuclide, Bayesian analysis will provide a unique posterior activity distribution that is consistent with all the energy peaks. In contrast, this is something neither conventional first-order Taylor expansion nor MC uncertainty propagation can cope with. Indeed, in this case the uncertainty propagation would deliver separate activity distributions for each peak, with no straightforward way to combine them. Other situations of multiple observations informing about the same quantity arise when different radiological measurement techniques sensing the same radionuclide activities and/or relying on the same efficiency model are combined.

The pioneer conference paper by Clément et al. [6] introduced the use of Bayesian analysis in the radioactive waste characterization community, using a standard random walk Metropolis (RWM) Markov chain Monte Carlo (MCMC) algorithm [see 10, for details about MCMC]. In this work, we propose a fully Bayesian approach to infer the vertical distribution of the activities of the measured nuclides and to quantify the corresponding uncertainty. Our approach relies on state-of-the-art MCMC sampling using the Hamiltonian Monte Carlo (HMC) technique [19, 4] and accounts for the net and background count uncertainties (i.e., Poisson statistics) together with the uncertainty in the detector's efficiencies caused by the uncertainty in the source distribution. Furthermore, our efficiency model accounts for the contributions of all segments to each count measurement given the assumption of a constant matrix density and corresponding gamma ray attenuation. To the best of our knowledge, this is the first application of Bayesian inference to estimation and uncertainty quantification of spatially-distributed activity from gamma spectrometry.

The remainder of this paper is organized as follows. Section 2 describes the considered waste drum, the used SGS apparatus and the measured data, together with the considered model and the MCMC-based Bayesian inversion. Section 3 then details the results of our high-dimensional MCMC sampling for spatially-distributed radionuclides' quantification for both a synthetic and a real case, before section 4 provides some discussion and outlines some future developments. Finally, section 5 summarizes our main findings and provides a short conclusion.

The non-destructive technique used in this study to identify and quantify the gamma-emitters present in a radioactive waste drum is segmented gamma scanning (SGS), which scans the drum slice by slice (segment), while the waste drum is rotating. We have used the so-called 3AX SGS [5] . The 3AX device can scan in the horizontal plane while the drum is rotating. Here, we only used the rotation capability and varied the height of the 3AX detector along the vertical axis of the drum. A 3AX SGS is made of three parts: the mechanics, the detector system and the data acquisition and processing unit. The mechanical system consists of an evaluation unit for the turntable and two trolleys positioned respectively on the left and right hand side of the turntable opposite to one another. The right trolley holds the detector-collimator system while the left one can be equipped with a transmission source (which has not been used here). Both trolleys can move horizontally towards and from the drum. The SGS is equipped with a high purity germanium detector (HPGe) connected to a cooling unit using liquid nitrogen. This detector is fixed in a lead collimator, enabling it to "see" only a fraction of the waste package. The slice of a drum that can be seen by the detector is defined by the aperture of the collimator and the distance between detector and waste package. The output of the detector is amplified by a preamplifier and a spectroscopy linear amplifier. This pulse is converted into a digital number by the analogue-to-digital converter (ADC) for further computer-based processing, and is in our case combined with the multichannel analyzer (MCA). The used detector was calibrated by applying a custom correction function to another similar detector that is ISOCS/LabSOCS-calibrated [17] . Data acquisition is performed and controlled by the main program which communicates with the Genie-2000 [12] software for conventional activity estimation.

Note that for nuclear fuel measurements, better activity estimates could be obtained by using a dedicated isotopic composition software [e.g., FRAM, MGA 23, 18] instead of Genie-2000. However, here estimation accuracy is not as critical as usual because these values only serve to define a weakly informative prior distribution for the MCMC inference (see section 2.5).

We used the 3AX scanning results of a real non-conditioned radioactive waste package. The latter, a galvanized 200-liter drum (height: 890 mm; diameter: 580 mm; wall thickness: 1.25 mm), contains waste originating from the decommissioning of a MOX (Mixed Oxide Fuel) glove box. The matrix of the waste mainly consists of metals such as carbon steel (∼70 wt%) and lead (∼8 wt%). Additionally, a limited amount of mainly halogen-containing materials mixed with non-combustible waste (∼22 wt%) is present. The drum was filled up completely with a net weight of 54 kg. This roughly corresponds to an average matrix density of 0.27 g/cm 3 . The drum dimensions, filling degree, matrix composition and matrix density distribution are considered fixed for this study and the associated uncertainties are therefore not taken into account. Instead, we focus on one of the main sources of uncertainties for this kind of measurements: the source distribution uncertainty (see section 2.3 for details).

Considering 20 individual measurements starting from the bottom of the drum and raising the detector 46.8 mm for each next measurement, the drum can be discretized into 20 horizontal segments ( Figure 1 ): two segments of 23.4 mm height for bottom and top measurement (S1 & S20) and 18 segments of 46.8 mm height (S2 − S19). Measurement setup details are as follows:

• The detector is surrounded by a rectangular lead slit collimator with a 1 mm copper coating (5 • opening; 18 mm x 110 mm; 100 mm depth); • The 200-liter drum is placed on a turning table and is continuously rotated during the measurements. • The measured 20 segments result in 20 individual spectra (M1 up to M20) and measurement time is 300 s. For the first measurement, the bottom of the drum is placed in front of the middle of the detector (M1). For each following measurement the detector is raised by 46.8 mm. In this way, the top of the drum is placed in front of the middle of the detector for the 20 th measurement (M20, see Figure 1 ). Table 1 lists the identified nuclides in the drum and their associated energy peaks. 

To account for the effect of the uncertainty associated with the source distribution within the drum on the detector efficiencies, we modeled the whole measurement system using the complex cylinder model of the Geometry Composer V4.3 library of the ISOCS/LabSOCS software by ISOCS/LabSOCS -Mirion Technologies (Canberra), Inc. [17] . Taking the drum rotation during the measurement into account, we modeled two extreme source distributions ( Figure 1 ):

• Maximum efficiency: a homogeneous source distribution (hypothesis h) within each segment ( Figure 1a ).

• Minimum efficiency: a point source (hypothesis p) placed on the drum axes on the border of a segment (on the bottom of the segment for S1 up to S10 and on the top of the segment for S11 up to S20, see Figure 1b )

Since there is symmetry between measurements M1-M10 and measurements M11-M20 (Figure 2 ), the efficiency calculations for the lower part of the drum could be duplicated to the upper part of the drum. We designed an individual model for every individual measurement at the lower part of the drum (M1 up to M10), placing the source stepwise in segments S1 up to S20 for the homogeneous (hypothesis h) and point source (hypothesis p) scenarios, respectively. This means that we designed a total of 400 models (10 detector locations × 20 source locations × 2 source distributions). Overall, for both the h and p assumptions, the detector's efficiency was thus calculated for each configuration of source location (S, the source is located in a given segment while the other 19 segments do not contain any source), detector location (M), and photon energy (among the 12 energies listed in Table 1 ). These efficiencies were then encapsulated into 2 n source × n peaks × n seg 3D arrays, E h and E p , for the h and p assumptions, respectively. Here n source = n seg = 20 and n peaks = 12.

By using a λ ∈ [0, 1] coefficient, one can balance the efficiency (e) associated with a given detection at energy e between the h and p assumptions for a given segment :

. This gives rise to our following count simulation model for a given count rate, c rate

where the subscript i indexes the detector location, i = 1, · · · , 20, the subscript j indexes the considered nuclide, j = 1, · · · , 5, the subscript k indexes the considered energy peak, k = 1, · · · , 12, a j,k=1···n k ,l denotes the activity of nuclide j at segment location l that is associated with the k = 1 · · · n k enery peaks, and P γ j,k is the emission probability of nuclide j at energy k.

The simulated gross count at detector location i for the energy k is then calculated as c gross

where t m is the measurement time and b i,k is the background continuum count at detector location i for the energy k.

A common representation of the forward problem is

is a deterministic forward model with parameters θ and the noise term e lumps all sources of errors.

In the Bayesian paradigm, parameters in θ are viewed as random variables with a posterior probability density function (pdf), p (θ|d), given by

where L (θ|d) ≡ p (d|θ) signifies the likelihood function of θ. The normalization factor p (d) = p (θ) p (d|θ) dθ is not required for parameter inference when the parameter dimensionality is fixed. In the remainder of this paper, we will thus focus on the unnormalized density p (θ|d) ∝ p (θ) L (θ|d).

If we assume d to follow a Poisson process, which is the norm for count data, L (θ|d) can be written as

Here the vector of inferred variables, θ, consists of the 100 activities, a (5 nuclides × 20 segments), 237 background continuum counts b (12 energy peaks × 20 segments less the 3 measurements for which neither a background continuum count nor a net count could be measured) and 20 λ coefficients. As of the d vector, it comprises N d = 237 gross counts obtained from the sum of the 237 net counts and 237 background continuum counts. In addition, our derived posterior distribution is not only conditioned to d but also to the chosen E h and E p arrays. This results in assessment of the

The marginal posterior distribution of a given quantity, say a 1 , is obtained by integrating the posterior distribution over all other inferred variables

where the a ∼1 vector contain all elements of a but a 1 .

We assume the prior distributions for a, b and λ to be independent, which gives

For the activities, a, we have some prior knowledge in the form of estimates,â, derived by application of the Genie-2000 [12] procedure to each segment separately. Assuming that these estimates have an accuracy of ± 1 order of magnitude gives an accuracy in base 10 logarithmic scale of: log 10 (â) ± 1. Further assuming that log 10 (a) is normally distributed and that log 10 (â) ± 1 provides a 99% uncertainty interval induces a standard deviation, σ a , for log (a) of σ a = 1 3 log (10) where the log (10) term accounts for the change of basis between log 10 (·) and log (·). This forms the rationale for using a lognormal prior for a: p (log (a)) = N log (â) , σ 2 a I . As demonstrated in section 3 this is a weakly informative prior from which the posterior distribution will easily depart if needed to appropriately fit the count data. Note that the values inâ are likely to be overestimations as interpreting each segment separately implies that all the measured counts at a given detector location are assigned to the activity of the corresponding segment, which is incorrect.

Cases for which the detected net count(s) for a given (set of) measurement(s) {i, k = 1, · · · , n k }, are too small compared to the derived background continuum count(s) for the Genie-2000 procedure to estimate an activity value, are easily included in the analysis as for non-influential activity value(s), the Bayesian approach will simply return the prior distribution. Overall, these "non-detects" concerned 43 out of the 100 inferred activities, including the 20 Pu-240 activities. In these cases the Genie-2000 procedure returns a calculated minimum detectable amount (MDA) which means that in absence of any other information, the "true" value could be anything between zero and something somewhat close to the calculated MDA. For those activities, a mda , we decided to also use a lognormal prior but with mean parameter set to MDA/2: p (log (a mda )) = N log (0.5MDA) , σ 2 a I . Indeed, some practitioners use half the MDAs as working values. This is arguably subjective and another choice could be to use an uniform prior between zero and some upper bound related to the MDA, say: p (a mda ) = U (0, 2MDA). Note however that no matter the specified upper bound, the inference should point out automatically that activities above the MDA are inconsistent with the observed counts. This is further discussed in section 4.

Regarding the 237-dimensional background continuum count vector, b, the obvious choice is to set p (b) as the product of independent Poisson prior distributions with shape parameter equal to the (rounded) measured background continuum counts: P ois (b). However, for technical reasons linked to the used software (see next section) we cannot use a Poisson prior (in short, Poisson priors cannot be assigned to real variables). Instead we use for p (b) an uncorrelated and independent normal prior distribution with mean and variance vectors both equal to the measured count values:

provides an increasingly accurate approximation to P ois (x) as x increases [see, e.g., 3]. The approximation is commonly deemed excellent for x > 1000 and reasonably accurate for x > 10 [24] . Lastly, the prior distribution for the 20-dimensional λ vector is taken as a bounded uniform distribution: U (0, 1).

As no analytical solution of the 357-dimensional p a, b, λ|d, E h , E p distribution is available, we sample from p a, b, λ|d, E h , E p by MCMC simulation [see 10] using a state-of-the-art implementation of the HMC sampler [see 19, 4, for an extensive description]. Convergence of the MCMC to the posterior target is monitored by means of the potential scale reduction factor,R [9, 10] , using four independent Markov chains evolved in parallel (see next section). TheR statistic compares for each parameter of interest the average within-chain variance to the variance of all the Markov chains mixed together. The closer the values of these two variances, the closer to one the value of theR diagnostic. Values ofR jointly smaller than 1.02 for all sampled variables are deemed to indicate convergence to a limiting distribution.

We used the open-source greta package [13] to perform the HMC-based MCMC sampling.The greta package is an R [22] interface to some of the MCMC sampling algorithms implemented in the Tensorflow-probability package [TFP, 7] which itself relies on the Tensorflow (TF) machine learning platform [1] . The standard language to interact with TFP is Python [21] and the TFP-Python syntax is somewhat complicated. The main advantage of greta is to make it straightfoward to build probabilistic models by constructing a directed acyclic graph and to perform MCMC sampling with TFP, using a seamless R-like syntax. The most useful MCMC sampler available through greta and used herein is HMC [19] . The HMC sampler has proven to be quite efficient when the gradient of the likelihood or objective function can be calculated by applying the chain rule, using an automatic differentiation (AD) technique such as implemented in TFP. Not all numerical models are suited to AD (for instance, numerical solvers of partial differential equations are not) but numerical models used for routine interpretation of radiological characterization measurements typically are. The TFP MCMC algorithms can be run in parallel on both CPUs and GPUs which makes greta to be quite computationally efficient. In this study, we ran four separate HMC trials in parallel over 4 CPUs.

The four parallel HMC runs were performed for a total of 30,000 warmup iterations each, after which 10,000 posterior samples were collected per run. TheR convergence criterion computed using these 4 Markov chains was below 1.02 for all of the 357 sampled variables, indicating official convergence to the posterior target by the four runs. The collected 40,000 samples were thus used to approximate p a, b, λ|d, E h , E p . On the used 4-core notebook (equipped with Intel i7-6820HQ CPU @ 2.70GHz), achieving these 40,000 iterations per/run incurred a computational time of about 2 hours in total.

Before proceeding with the inversion of the real data, we study in this section whether our proposed approach is able to correctly recover the activities for a synthetic problem for which the true values are known. We used the same drum, efficiency data and count simulation model as for the real case, but considered only one nuclide and energy peak: Pu-239 at 413.71 keV. The "true" λ was set to 1 in every segment except for segment 16 where it was set to 0. The Pu-239 activity was set to 1×10 8 Bq in every segment except for segments 5 and 16 were it was set to 1×10 11 Bq. The resulting i = 1, · · · , 20 simulated net counts, c net true,i , were used as mean parameters of 20 Poisson distributions from which the 20 observed nets counts were drawn: c net obs,i ∝ P ois c net true,i . The 20 observed background continuum counts, b obs,i , were also sampled from a Poisson distribution with a mean parameter inspired from the real measurements: c back obs,i ∝ P ois (100). Note that due to the added errors in the observed net counts, no perfect inverse solution exists which is required for Bayesian inference. The same normal p (b) prior distribution and uniform p (λ) prior distribution as for the real case were used. We defined p (a) = N log (â synth ) , 3σ 2 a I with the 20 elements of theâ synth vector being 1 × 10 9 Bq. This way, the prior mean of the Pu-239 log-activity is (i) one order of magnitude too large for the 18 segments where the "true" Pu-239 activity of the waste (1 × 10 8 Bq) is small compared to the background activity, and (ii) two order of magnitude too large for the 2 segments where the "true" Pu-239 activity is large (1 × 10 11 Bq). In addition, setting the standard deviation of the log-prior activity to 3σ 2 a (that is, 1 log (10)) induces a wide and thus quite uncertain prior. (Figure 2a) , it is seen that the MCMC sampling recovers the 2 large true activities in segments 5 and 16 very accurately, with a high degree of certainty. This even though the mode of the lognormal prior is way off. The posterior log-activity distributions in the other segments (all but 5 and 16) often peak near the true values and are somewhat asymetric with a larger probability mass at the right side. This is because here the MCMC sampling finds a balance between simulating sufficiently small gross counts and honoring the prior. This nicely illustrates how Bayesian analysis handles MDAs (non-detects). If the range of values spanned by the prior distribution is sufficiently low, then it simply returns the prior distribution. Otherwise, the posterior distribution forms an intermediate distribution between the prior and the likelihood (equation (4)). The posterior background continuum count distributions globally peak around the true values, with slight deviations from their associated priors ( Figure 2b ). As expected, the posterior λ distributions are the same as their prior counterparts for the low activities (Figure 2c ). For the 2 segments where the true activity is large, the inference put most of the probability mass towards the true boundary values of λ, λ = 1 for segment 5 and λ = 0 for segment 16. Overall these results demonstrate that our approach can correctly infer the joint posterior distribution of a, b and λ. The black dots denote the true values used to generate the gross count data. Note that the x-axis of subplot (a) has a base 10 logarithmic scale. For visual convenience, the background continuum counts have been standardized using the normal prior parameters. For the λ variable, we show a posterior histogram instead of a kernel density estimate. This is because when applied to a bounded data sample, kernel density smoothing tends to create artifacts near the bounds. (Table 1) , are globally well resolved. This is indicated by a posterior distribution that is narrower than the prior distribution. Strikingly, for all 100 activities but that of Pu-239 in segment 16, the posterior mode is lower than the prior mean of the log-activity (log (â) or log (0.5MDA), see section 2.5). As written above, this is because processing each segment individually, as done to define p (a), is likely to overestimate the actual activities as it overlooks the fact that the detected counts at a given detector location are not only caused by the activity in the considered segment but also by the activities in the neighboring segments. In contrast, for a given nuclide our net count rate simulation model (see equation (1) and associated text) accounts for the contributions from all activities to each measurement. For Pu-239 at segment 16 only (Figure 3) , the posterior mode of the log-activity is above the selected prior mean: about 3.9 × 10 9 Bq against 1.13 × 10 9 Bq. Note that since this activity is among those for which the Genie-2000 procedure returned a MDA value, here the posterior mode is about two times the MDA (prior mean of the log-activity was taken as log (M DA/2)).

The "single-energy" Pu-240 and Pu-241 nuclides are relatively well resolved though, perhaps not surprisingly, some are more uncertain than the multi-energy nuclides. This is especially the case for Pu-240 at segment 1 and Pu-241 at segments 13 − 16 and 20. Overall, at these locations the 95% uncertainty intervals for Pu-240 and Pu-241 cover about one order of magnitude. Table 2 lists the corresponding 95% posterior uncertainty (or credible) intervals of the total activity of each nuclide over the whole drum. For this package, the Pu-241 nuclide has by far the largest activity, one to two orders of magnitude higher than for the other nuclides. As written earlier, the measured net count for the 20 Pu-240 activities is zero and, in this case, the prior distribution was set to N log (0.5MDA) , σ 2 a I . Since the range spanned by the posterior Pu-240 log-activities is systematically lower than the range spanned by the prior Pu-240 log-activities, our results show that for the considered case study setting the Pu-240 non-detects to half the MDA is overly conservative. (20) . Note that the x-axis of subplot (a) has a base 10 logarithmic scale. Note that for the λ variable, we show a posterior histogram instead of a kernel density estimate. This is because when applied to a bounded data sample, kernel density smoothing tends to create artifacts near the bounds.

The marginal posterior background continuum activities, p b|d, E h , E p , are depicted in Figure 6 . The posterior distributions are generally close to their associated prior distribution but some depart from it. Consistently with our Bayesian framework, these deviations are needed to maximize the posterior density. We would like to stress that the prior, p (b), is based on a measured background continuum count that is considered to be a realization of a Poisson distribution of which the shape parameter is unknown. By setting p (b) = P ois (b) (or p (b) = N (b, C b ) ≈ P ois (b) as done herein) one assumes that the measured background continuum count is equal to the mean of its underlying distribution, which is obviously not necessarily the case. Some deviations from p (b) should thus come at no surprise. If deemed necessary, an alternative solution to enforce a tight closeness between N (b, C b ) ≈ P ois (b) and p b|d, E h , E p is discussed in section 4. (1) to top (20) . Figure 5b presents the marginal posterior densities of the elements of the λ vector, which is our parameter that accounts for the uncertainty in the efficiencies. For segments 2 to 19 the posterior probability mass is mostly concentrated on the upper bound of 1, meaning that the MCMC is strongly favoring the hypothesis of a homogeneous radioactive source distribution in these locations. In contrast, for the bottom (1) and top (20) segments every intermediate state between a fully homogeneous source distribution and a fully point source distribution has non-zero probability. it is important to note that our efficiency model (see section 2.3) might be biased, especially for the top segment. Hence, a smaller filling height than the drum top might cause the large uncertainty in the source distribution for top segment 20 and the long-tail uncertainty distribution in the underlying segment 19.

The fit of the posterior gross count simulations, F a, b, λ, E h , E p to the measured gross counts, d is depicted in Figure 7 using a base 10 logarithmic scale to make the discrepancies between (very) small measured and simulated counts visible. The vertical bars denote the 95% posterior uncertainty intervals associated with the simulated counts. The observed counts are mostly well fitted and the 95% uncertainty intervals are relatively tight and most of the time include the 1:1 line. Furthermore, the most important discrepancies, in relative terms, are for some rather small gross counts, below 10 − 15. This is consistent with both our expectations and the Poisson count statistics. 

We have shown that vertical activity profiles can be probabilistically inferred from SGS data accounting for the effect of the source distribution uncertainty on the detector efficiency. The following points and/or potential improvements nevertheless deserve further attention.

As for any prior distribution, the choice of the priors for the activities flagged as MDAs by the standard analysis procedure [12] applied to each segment separately, a mda , is a relatively subjective expert-based decision. Note, however, that if the prior is sufficiently weakly informative or wide and the information content of the count measurement data is large enough (sufficiently large net counts compared to the corresponding background counts and not too many missing values), then the actual shape of the prior should not substantially influence the posterior results. An additional MCMC run with using p (a mda ) = U (0, 2MDA) instead of p (log (a mda )) = N log (0.5MDA) , σ 2 a I shows that, unsurprisingly, p (a mda ) = U (0, 2MDA) leads to posterior distributions for the elements of a mda that are close to their associated uniform prior (Figure 8 ). For the non-MDAs, the two posterior distributions are very similar (Figure 8 ).

In some relatively few cases, the marginal posterior background continuum activities, p b|d, E h , E p , differ largely from the used prior distribution for b, N (b, C b ) ≈ P ois (b). This is probably due to some small inconsistencies between our model and the real drum. Notwithstanding, if having p b|d, E h , E p closely approximating N (b, C b ) ≈ P ois (b) is judged necessary, then setting p (b) = N (b, s b C b ) with s b < 1 can do the job. Indeed, we found for the considered case study that using s b = 0.01 allows for p b|d, E h , E p to closely reproduce N (b, C b ) while only slightly changing p a|d, E h , E p . However, this comes at the expense of a less good fit of the simulated gross counts, F a, b, λ, E h , E p to the measured gross counts, d (not shown).

This study only considers the count measurement uncertainty and the uncertainty of the source distribution (albeit in a simplified way), the latter being one important source of uncertainty in the detector's efficiency. Nonetheless, the matrix density distribution of the package is assumed to be homogeneous. Since the matrix-related uncertainties can be relatively large and thereby impacting the detector's efficiency too, more work is needed to account for it within the inference. This could be approached by using a Monte Carlo particle transport code to simulate the impact of various matrix and source configurations on the efficiency of the considered detector. The so-derived dataset could then be used to construct a computationally cheap nonlinear regression model (also called proxy or metamodel) that would predict the simulated efficiencies from the sampled matrix and source parameters, for each considered nuclide and segment location. Both matrix-related and source distribution uncertainty could thus be accounted for within the Bayesian inference by sampling the surrogate model parameters along with the other variables. This warrants further investigations.

More generally, an extended Bayesian approach based on the one presented herein could also be used to infer spatial activity distributions in 1D, 2D and 3D within a drum waste package from combined measurements from different techniques. This will be investigated in future work, in the framework of the MICADO EU project (see section 6). (1) p (log (a mda )) = N log (0.5MDA) , σ 2 a I (Posterior 1) and (2) p (a mda ) = U (0, 2MDA). For visual convenience, the lower limit of the x-axis is thresolded to 1 × 10 7 Bq. Also, note the base 10 logarithmic scale of the x-axis.

We propose a Bayesian approach to probabilistically infer vertical activity distributions within a radioactive waste drum from segmented gamma scanning (SGS) measurements. Our approach relies on state-of-the-art Markov chain Monte Carlo (MCMC) sampling using the Hamiltonian Monte Carlo (HMC) technique and accounts for the net and background count measurement uncertainty together with the uncertainty in the detector efficiencies caused by the uncertainty in the source distribution within the drum. Furthermore, our efficiency model accounts for the contributions of all considered segments to each count measurement. For the considered case study, our approach is used to resolve the vertical activity distribution of 5 nuclides over 20 locations and produce sound uncertainty estimates. Main topics for future research include how to incorporate the impact of the uncertainty in matrix density distribution on the modeled efficiencies, to extend our setup from 1D to 2D and/or 3D, and to combine measurements from different techniques for better radionuclides' quantification.

This work received funding by the EU project MICADO Consortium. 2018: "Horizon 2020, Nfrp 2018-10, Proposal 847641: Measurement and Instrumentation for Cleaning and Decommissioning Operations".

Large-Scale Machine Learning on Heterogeneous Distributed Systems. ArXiv e-prints

An improved method for the nondestructive characterization of radioactive waste by gamma scanning

Poisson Approximation

A Conceptual Introduction to Hamiltonian Monte Carlo

Synopsis of gamma scanning systems -comparison of gamma determining systems and measuring procedures for radioactive waste packages. European Network of Testing Facilities for the Quality Checking of Radioactive Waste Packages

Stochastic approach for radionuclides quantification

A new gamma spectroscopy methodology based on probabilistic uncertainty estimation and conservative approach

Inference from iterative simulation using multiple sequences

Customisation Tools Manual 9230847F 3/01. Mirion Technologies (Canberra)

greta: simple and scalable statistical modelling in R

Minimum detectable activity, systematic uncertainties, and the ISO 11929 standard

Reconstruction of the activity of point sources for the accurate characterization of nuclear waste drums by segmented gamma scanning

Situ Object Counting System) Calibration Software

ISOCS/LabSOCS Detector Characterization Report. CSO #93483. 2020. Mirion Technologies

Group Analysis Software

Handbook of Markov Chain Monte Carlo

Segmented gamma-ray assay of large volume radioactive waste drums containing plutonium lumps

Python Software Foundation. 2020

R: A language and environment for statistical computing. R Foundation for Statistical Computing

Application guide to gamma-ray isotopic analysis using the FRAM software, LA-14018