Determination of Nuclear PDFs using Markov Chain Monte Carlo Methods

1. Problem Statement

Global Quantum Chromodynamics (QCD) analyses of Nuclear Parton Distribution Functions (nPDFs) have traditionally relied on the Hessian method for uncertainty estimation. This approach assumes a Gaussian approximation of the likelihood function around a single minimum, expanding the $\chi^2$ function to quadratic terms. However, this methodology faces significant limitations in the context of nPDFs:

• Non-Gaussian Likelihoods: Nuclear data often suffer from limited constraints and tensions between datasets, leading to likelihood surfaces that are skewed, non-linear, or possess multiple competing minima.
• Inadequacy of Local Approximations: The Hessian method fails to capture complex parameter-space structures (e.g., multimodal distributions, asymmetric tails) and relies on an empirically chosen tolerance parameter ($\Delta\chi^2$) to inflate uncertainties, which weakens the probabilistic interpretation of confidence intervals.
• Complex Correlations: The presence of both protons and neutrons in nuclei, linked by isospin symmetry, introduces intricate correlations between parameters that are difficult to model with simple error ellipsoids.

The paper addresses the need for a statistically rigorous framework that can directly sample the full posterior distribution to accurately quantify uncertainties in nPDFs, particularly for the lead (Pb) nucleus.

2. Methodology

The authors implement a Bayesian framework using Markov Chain Monte Carlo (MCMC) techniques within the nCTEQ framework.

• Formalism:

  • Nuclear PDFs: The nPDFs are defined via effective bound-proton and bound-neutron PDFs, assuming isospin symmetry. The $x$-dependence is parametrized at an input scale ($Q_0 = 1.3$ GeV) using a functional form based on the CJ15 proton baseline.
  • Nuclear Modifications: Parameters are made dependent on the nuclear mass number $A$ via a logarithmic expansion ($c_j(A) = p_j + a_j \ln A + b_j \ln^2 A$).
  • Likelihood: The posterior distribution is proportional to the likelihood, constructed from the standard $\chi^2$ function incorporating statistical, correlated systematic, and normalization uncertainties (using the T-method).

• Algorithm:

  • The study employs an Adaptive Metropolis–Hastings (aMH) algorithm. Unlike standard MCMC, aMH updates the proposal covariance matrix recursively during the sampling process to adapt to the local structure of the parameter space.
  • Diagnostics: The authors rigorously apply burn-in removal (thermalization), convergence checks (trace plots, cumulative means), and thinning (to reduce autocorrelation) to ensure the final ensemble represents independent samples from the stationary distribution.

• Analysis Strategy:
  Two complementary analyses were performed to disentangle methodological effects from modeling assumptions:

  1. Pb-only Fit: Uses exclusively data from the Lead ($A=208$) nucleus (LHC vector boson production, heavy-flavor production, and neutrino DIS). This isolates the parameter-space structure of a single nuclear system.
  2. Multi-Nuclei Fit: Incorporates data from a wide range of nuclei (Deuterium to Lead) using the standard analytic $A$-dependence to constrain the lead PDFs indirectly.

3. Key Contributions

• First MCMC-based nPDF Determination: This is the first comprehensive global analysis of nPDFs using MCMC techniques with real-world experimental data.
• Discovery of Complex Parameter Structures: The MCMC approach revealed a highly non-trivial likelihood landscape for nPDFs, including multiple minima and pronounced non-Gaussian behavior, particularly for valence quark distributions.
• Comparison of Uncertainty Methods: The paper provides a direct, systematic comparison between MCMC-derived uncertainties and traditional Hessian uncertainties, highlighting the limitations of the latter in non-Gaussian regimes.
• Single-Nucleus Extraction: The Pb-only analysis represents the first extraction of nuclear PDFs for a single nucleus (Lead) using only Lead data, made possible by the wealth of recent LHC proton-lead collision data.

4. Key Results

A. Posterior Distribution and Multiple Minima

• Valence Quarks: The posterior distributions for valence quark parameters (specifically $a_{uv}^3$ and $a_{dv}^3$) exhibit significant deviations from Gaussianity.
  • Multiple Modes: The analysis identified a secondary minimum in the parameter space. For the $dv$ parameters, a distinct region with $a_{dv}^3 > 0.8$ was found, separated by a barrier from the main minimum. This secondary minimum corresponds to a different shape of the valence PDFs (enhancement at low-$x$, depletion around $x \sim 10^{-2}$).
  • Barriers: For $uv$ parameters, a clear barrier exists around $a_{uv}^3 \approx 0$, splitting the posterior into positive and negative regions.
• Sea Quarks and Gluons: In contrast, the sea quark and gluon parameters exhibit nearly Gaussian distributions with simple correlations, as they are tightly constrained by the data.

B. Uncertainty Quantification

• MCMC vs. Hessian:
  • For gluons and sea quarks, MCMC and Hessian uncertainties are consistent.
  • For valence quarks, the Hessian method significantly underestimates uncertainties in regions of non-Gaussian behavior. The MCMC method (specifically the cumulative-$\chi^2$ approach) yields broader, asymmetric uncertainty bands that correctly capture the probability mass in the tails and secondary minima.
  • The Hessian method fails to detect the existence of the secondary minimum, effectively ignoring a valid region of the parameter space that contributes to the true uncertainty.

C. Impact of Multi-Nuclei Data

• Uncertainty Reduction: Including data from lighter nuclei reduces uncertainties for quark distributions (valence and sea) across the $x$-range, as expected from the added constraints.
• Shape Modification: The inclusion of lighter nuclei alters the shape of the extracted Lead PDFs. Specifically, the nuclear modification effects (shadowing/anti-shadowing) appear "inverted" or weakened compared to the Pb-only fit. This suggests that the assumed analytic $A$-dependence introduces systematic biases when extrapolating constraints from light nuclei to heavy ones.
• Gluon Stability: The gluon distribution and its uncertainty remain largely unchanged between the Pb-only and multi-nuclei fits, confirming that gluon constraints are dominated by Lead-specific data (LHC heavy-flavor data).

5. Significance and Conclusion

This work establishes that MCMC methods provide a more reliable and rigorous framework for nPDF uncertainty quantification compared to the traditional Hessian approach.

• Robustness: MCMC naturally handles non-Gaussian likelihoods, multimodal structures, and complex parameter correlations without relying on arbitrary tolerance parameters.
• Methodological Insight: The study demonstrates that the "true" uncertainty in nPDFs, particularly for valence quarks, is significantly larger and more complex than previously estimated by Hessian methods.
• Future Directions: The ability to extract PDFs for a single nucleus (Lead) using only Lead data marks a paradigm shift, reducing reliance on phenomenological $A$-dependence assumptions. The paper concludes that while MCMC is computationally more intensive, it is essential for precision QCD phenomenology where data constraints are limited and likelihoods are complex.

The authors emphasize that the complex landscape of nPDF parameter spaces is not unique to this problem but is likely common in other areas of particle physics with limited data, making MCMC a critical tool for future global analyses.