Bayesian approach for uncertainty quantification of hybrid spectral unmixing in $\gamma$-ray spectrometry

This paper proposes and evaluates two Bayesian methods, Laplace approximation and Markov Chain Monte Carlo, for quantifying the uncertainty of hybrid spectral unmixing estimators in $\gamma$-ray spectrometry, demonstrating that while both perform well under ideal conditions, Markov Chain Monte Carlo remains robust when spectral constraints or dominant backgrounds create non-Gaussian posterior distributions where Laplace approximation fails.

The Gist

The Big Picture: Solving a "Noisy" Puzzle

Imagine you are trying to identify different types of fruit in a giant, blurry smoothie. You know the smoothie contains apples, bananas, and oranges, but the blender has mixed them up, and the container is slightly foggy.

In the world of Gamma-ray spectrometry, scientists face a similar problem. They have a detector that "sees" radiation from radioactive materials (like Cobalt-60 or Cesium-137). However, the radiation doesn't just travel in a straight line; it bounces off walls, gets absorbed by steel shielding, or scatters. This is like the "fog" or the "blur" in our smoothie.

Because of this, the "signature" (the unique pattern) of each radioactive material changes depending on what it's passing through. A few years ago, the authors developed a smart computer program called SEMSUN (a hybrid of machine learning and math) that can look at this blurry, mixed-up data and guess:

1. How much of each radioactive material is there? (The Counting)
2. How much has the signal been distorted by the environment? (The Variable $\lambda$)

The Problem: The SEMSUN program is great at making a guess, but it doesn't tell you how sure it is. If you are making safety decisions (like "Is this area safe to enter?"), you need to know the margin of error. This paper is about building a "confidence meter" for that program.

---

The Two Methods: The "Quick Sketch" vs. The "Deep Dive"

To figure out how confident we should be in the results, the authors tried two different approaches. Think of them as two ways to predict the weather.

1. The Laplace Approximation (LA) – "The Quick Sketch"

Imagine you are a meteorologist who wants to predict tomorrow's temperature. You look at the data and draw a nice, perfect Bell Curve (a smooth, symmetrical hill). You assume the temperature will likely be right in the middle, with fewer chances of it being very hot or very cold.

• How it works: This method assumes the uncertainty follows a perfect, smooth bell curve. It's very fast (less than a tenth of a second) and easy to calculate.
• The Catch: Real life isn't always a perfect bell curve. If the data is pushed against a wall (like a rule that says "you can't have negative radiation"), the curve gets squashed and looks weird. The "Quick Sketch" method doesn't handle these weird shapes well. It might say, "I'm 95% sure," when it's actually only 60% sure.

2. The Markov Chain Monte Carlo (MCMC) – "The Deep Dive"

Now, imagine a different meteorologist. Instead of drawing a curve, they run a supercomputer simulation 1,000 times. They say, "Okay, let's pretend it's a windy day, then a rainy day, then a sunny day..." and they generate thousands of possible outcomes based on the rules of physics.

• How it works: This method doesn't assume a shape. It literally samples the data thousands of times to see what the "real" distribution looks like. It builds a map of all possibilities.
• The Catch: It takes a long time (several minutes) and requires a lot of computing power. It's like running a marathon instead of a sprint.

---

The Experiment: Who Got It Right?

The authors tested both methods using a "Long-Run Success Rate" (LRSR). This is a fancy way of saying: "If we run this test 100 times, how often does our '95% confidence interval' actually catch the true answer?"

They wanted the answer to be 95.4 times out of 100.

The Results:

1. When things are easy (The "Open Field"):
   If the radiation levels are high and the distortion isn't hitting any "walls" (mathematical limits), both methods work great. The "Quick Sketch" (LA) and the "Deep Dive" (MCMC) give almost the same result.

   • Verdict: Use the Quick Sketch because it's fast.

2. When things get tricky (The "Cornered" or "Noisy" Field):

   • Scenario A: The radiation is very low (like a whisper in a noisy room).
   • Scenario B: The background noise is huge compared to the signal.
   • Scenario C: The distortion is at its maximum or minimum limit (hitting the "walls").

   In these cases, the "Quick Sketch" (LA) fails. Because it forces the data into a perfect bell curve, it gives a false sense of security. Its success rate drops way below 95%.
   However, the "Deep Dive" (MCMC) keeps its cool. Because it actually simulates the messy reality, it still hits the 95% target.

   • Verdict: Use the Deep Dive. It's slower, but it won't lie to you.

The "User Guide" Conclusion

The paper concludes with a practical rule for scientists:

• Step 1: Run the fast "Quick Sketch" (LA) method first.
• Step 2: Check if the data is "stuck" against a wall (e.g., is the estimated distortion at the very edge of possible values? Is the background noise drowning out the signal?).
• Step 3:
  • If the data looks "free" and smooth? Great! Trust the fast result.
  • If the data looks "stuck" or messy? Stop! Switch to the slow "Deep Dive" (MCMC) method to get a trustworthy answer.

Summary in One Sentence

This paper teaches us how to know when a fast, simple math trick is good enough to measure radiation uncertainty, and when we need to slow down and run a heavy-duty computer simulation to avoid making dangerous mistakes.

Technical Summary

1. Problem Statement

In $\gamma$-ray spectrometry, accurately identifying and quantifying radionuclides is challenging when spectral signatures are deformed by physical phenomena such as attenuation and Compton scattering in the source surroundings (e.g., shielding materials).

• The Core Issue: Traditional methods assume fixed spectral signatures. However, in in situ measurements, signatures vary, requiring a "semi-blind" unmixing approach where both the radionuclide counts ($a$) and the deformation parameters ($\lambda$) must be estimated simultaneously.
• The Gap: While a hybrid machine learning algorithm (SEMSUN) was previously developed to solve this inverse problem, it lacked a robust framework for Uncertainty Quantification (UQ). Specifically, there was no method to calculate reliable Coverage Intervals (CI) (equivalent to Credible Intervals in Bayesian statistics) for the estimated counts and deformation parameters, which is essential for metrological decision-making.

2. Methodology

The authors propose a Bayesian framework to quantify the uncertainty of the estimators $\theta = (a, \lambda)$ obtained from the SEMSUN algorithm. They develop and compare two distinct Bayesian methods:

A. The Model

The observed spectrum $y$ is modeled as a Poisson distribution $y \sim P(Xa)$, where $X$ is the spectral signature matrix.

• Spectral Variability: The matrix $X$ is not fixed; it depends on a latent variable $\lambda$ (representing physical conditions like shield thickness) via a non-linear function $f(\lambda)$ learned by an Interpolating AutoEncoder (IAE).
• Optimization: The SEMSUN algorithm finds the Maximum A Posteriori (MAP) estimate, which coincides with the Maximum Likelihood Estimate (MLE) under uniform priors.

B. Uncertainty Quantification Methods

1. Laplace Approximation (LA):

   • Principle: Approximates the posterior distribution as a Gaussian centered at the MAP estimate.
   • Implementation: Uses the inverse of the Hessian of the log-posterior (equivalent to the Fisher Information Matrix) to estimate the covariance. Since the spectral deformation function $f(\lambda)$ is a neural network, analytical derivatives are unavailable; the authors use PyTorch's autograd to compute gradients numerically.
   • Assumption: The posterior is Gaussian. This holds only if constraints (e.g., $a \geq 0$, $0 \leq \lambda \leq 1$) are inactive.

2. Markov Chain Monte Carlo (MCMC):

   • Principle: Directly samples the posterior distribution without assuming a specific shape (e.g., Gaussian).
   • Implementation: Uses the No-U-Turn Sampler (NUTS), an extension of Hamiltonian Monte Carlo, implemented via the Pyro package.
   • Output: Generates a chain of samples to construct the Highest Posterior Density (HPD) interval, which is the narrowest interval containing the specified probability mass.

C. Validation Metric: Long-Run Success Rate (LRSR)

To validate the Coverage Intervals, the authors use the Long-Run Success Rate (LRSR).

• Procedure: They simulate thousands of $\gamma$-spectra based on known "true" values of $a$ and $\lambda$.
• Calculation: For each simulation, they calculate the CI using LA and MCMC. The LRSR is the percentage of times the calculated CI contains the true value.
• Target: A 95.4% CI should yield an LRSR close to 95.4%.

3. Key Contributions

• First UQ for Hybrid Unmixing: Provides the first uncertainty quantification framework specifically for the SEMSUN hybrid algorithm, addressing the complex inverse problem of spectral deformation.
• Comparative Analysis: Rigorously compares the computationally efficient Laplace Approximation against the robust but expensive MCMC method.
• Applicability Criterion: Establishes a practical criterion to determine when the LA method is valid (i.e., when constraints are inactive and the posterior is effectively Gaussian) versus when MCMC is required.
• User Guide: Introduces a workflow for applying these UQ methods in real-world $\gamma$-ray spectrometry.

4. Results

The study evaluated two scenarios: a simple mixture (Background + $^{60}$Co) and a complex mixture (Background + $^{60}$Co, $^{137}$Cs, $^{88}$Y, $^{99m}$Tc).

• Zone 1 (Inactive Constraints): When the estimated counts are high and the deformation parameter $\lambda$ is far from boundaries (0 or 1), the posterior distribution is Gaussian.

  • Result: Both LA and MCMC yield LRSR values close to the expected 95.4%.
  • Efficiency: In this zone, LA is highly advantageous, taking <0.1 seconds compared to several minutes for MCMC. The Total Variation Distance (TVD) between the two distributions is negligible.

• Zones 2, 3, & 4 (Active Constraints/Low Statistics):

  • Scenario: Occurs when $\lambda$ is near 0 or 1 (extreme shielding) or when background counts dominate the signal (low signal-to-noise ratio).
  • Result: The posterior distribution becomes non-Gaussian (truncated or skewed).
  • LA Performance: The LA method fails, producing LRSR values significantly lower than 95.4% (e.g., dropping to ~60% in extreme low-count cases) because the Gaussian assumption breaks down.
  • MCMC Performance: The MCMC method remains robust, maintaining LRSR close to 95.4% by accurately capturing the non-Gaussian shape of the posterior.

• High Statistics Anomaly: At very high counting rates ($\sim 10^5$), both methods showed slight deviations from the expected LRSR due to reconstruction errors inherent in the IAE model, though the bias remained low.

5. Significance

• Robust Decision Making: The study demonstrates that relying solely on Gaussian approximations (LA) for uncertainty in spectral unmixing can lead to false confidence in critical scenarios (e.g., low statistics or heavy shielding).
• Practical Workflow: The authors propose a hybrid strategy:
  1. Run the fast LA method first.
  2. Check the constraint activation criterion (e.g., is $\hat{a} - 3\sigma < 0$?).
  3. If constraints are inactive, use the LA result.
  4. If constraints are active or background dominates, switch to MCMC to ensure accurate uncertainty quantification.
• Metrological Impact: This work bridges the gap between advanced machine learning spectral unmixing and strict metrological standards (GUM), enabling the reliable use of automated algorithms in nuclear security, decommissioning, and environmental monitoring.