[Paper Review] pyerrors: a python framework for error analysis of Monte Carlo data

Fabian Joswig, Simon Kuberski+2 more|arXiv (Cornell University)|Sep 28, 2022

Forecasting Techniques and Applications4 citations

TL;DR

pyerrors is a Python framework that enables precise error analysis of Monte Carlo data using linear error propagation via automatic differentiation and the $Γ$-method for autocorrelation time estimation. It integrates seamlessly into the scientific Python ecosystem, allowing reliable uncertainty propagation across observables from different simulation sources while handling irregular or gapped Monte Carlo chains.

ABSTRACT

We present the pyerrors python package for statistical error analysis of Monte Carlo data. Linear error propagation using automatic differentiation in an object oriented framework is combined with the $Γ$-method for a reliable estimation of autocorrelation times. Data from different sources can easily be combined, keeping the information on the origin of error components intact throughout the analysis. pyerrors can be smoothly integrated into the existing scientific python ecosystem which allows for efficient and compact analyses.

Motivation & Objective

• To address the challenge of accurate statistical error estimation in Markov Chain Monte Carlo (MCMC) simulations, particularly under critical slowing down and increasing autocorrelation in lattice QCD.
• To provide a user-friendly, extensible Python package that integrates naturally with existing scientific Python tools like NumPy for efficient data analysis.
• To enable reliable uncertainty propagation across observables derived from simulations with different parameters or configurations, preserving error component provenance.
• To support the analysis of irregular or gapped Monte Carlo histories, where data is missing or sampled on subsets of configurations.
• To offer a robust, machine-precision error propagation system using automatic differentiation, reducing underestimation of statistical errors.

Proposed method

• Introduces a new data type, `Obs`, which behaves like a NumPy array and supports automatic differentiation for linear error propagation to machine precision.
• Employs the $Γ$-method for estimating autocorrelation times, which outperforms traditional methods like Jackknife and Bootstrap in the presence of critical slowing down.
• Uses a modified Fourier-based approach via `rfft` to compute normalized autocorrelation functions efficiently, even for irregular chains.
• Adapts the $Γ$-method to handle gapped or irregular Monte Carlo histories by introducing indicator functions $n_{eta}^{i}$ and rescaling fluctuations to preserve error consistency.
• Replaces standard normalization $N - t$ with $\mathcal{N}_{\alpha}(t) = \sum_{i=1}^{N-t} n_f^i n_f^{i+t}$ to correctly account for missing configurations in autocorrelation estimation.
• Enables seamless integration with scientific Python workflows by allowing `Obs` objects to be used directly in NumPy operations and supporting advanced features like effective mass computation and GEVP solutions.

Experimental results

Research questions

• RQ1 How can statistical errors in MCMC simulations be reliably estimated when autocorrelations increase due to critical slowing down?
• RQ2 Can the $Γ$-method be effectively implemented in a high-level language like Python while maintaining computational efficiency and numerical accuracy?
• RQ3 How can error propagation be consistently tracked across observables derived from simulations with different parameters or configurations?
• RQ4 What is the impact of missing or irregularly sampled configurations on error estimation, and how can it be corrected?
• RQ5 To what extent can a Python-based framework achieve performance and reliability comparable to low-level language implementations (e.g., Fortran, Julia) in lattice field theory error analysis?

Key findings

• The pyerrors framework successfully implements the $Γ$-method for autocorrelation time estimation with high numerical accuracy, outperforming traditional resampling techniques like Jackknife and Bootstrap in the presence of critical slowing down.
• Error propagation via automatic differentiation in pyerrors achieves machine-precision accuracy, minimizing numerical errors in uncertainty estimation.
• The framework supports the combination of observables from different simulation sources while preserving the provenance of error components, enabling traceable uncertainty analysis.
• The implementation correctly handles irregular Monte Carlo chains by introducing rescaled fluctuations and configuration indicators, ensuring consistent error propagation even with missing data.
• The package has been successfully applied in multiple scientific publications, demonstrating its utility and robustness in real lattice QCD analyses.
• The use of `rfft`-based computation enables efficient autocorrelation estimation, making the framework suitable for large-scale simulations.