Error Estimation for Moments Analysis in Heavy-Ion Collision Experiments

Fluctuations of conserved quantities are predicted to be sensitive to the correlation length and connected to the thermodynamic susceptibility. Thus, moments of net-baryon, net-charge and net-strangeness have been extensively studied theoretically and experimentally to explore phase structure and bulk properties of QCD matters created in heavy ion collision experiment. As the moments analysis is statistics hungry study, the error estimation is crucial to extract physics information from the limited experimental data. In this paper, we will derive the limit distributions and error formula based on Delta theorem in statistics for various order moments used in the experimental data analysis. The Monte Carlo simulation is also applied to test the error formula.

💡 Research Summary

The paper addresses a crucial methodological issue in the analysis of event‑by‑event fluctuations of conserved charges (net‑baryon, net‑charge, net‑strangeness) in relativistic heavy‑ion collisions: the reliable estimation of statistical uncertainties for higher‑order moments and their products. Because the moments (variance σ², skewness S, kurtosis κ, and various cumulant ratios) are highly sensitive to the correlation length of the system and thus to the QCD phase structure, precise error quantification is indispensable for drawing physics conclusions from the limited data sets typical of collider experiments.

The authors begin by defining central moments μ̂_r = ⟨(N−μ̂)^r⟩ and cumulants C_n in terms of the event‑wise net‑proton number N. They then express the commonly used observables—σ = √μ₂, S = μ₃/μ₂^{3/2}, κ = μ₄/μ₂²−3—as functions of these cumulants, and introduce the products σ·S, κσ², κσ/S that directly correspond to ratios of thermodynamic susceptibilities (e.g., χ₃/χ₂, χ₄/χ₂).

The statistical core of the work is the application of the Delta theorem. First, Theorem A is presented, stating that the vector of sample central moments (μ̂₂, μ̂₃, …, μ̂_k) converges in distribution to a multivariate normal with a covariance matrix Σ that can be expressed analytically in terms of the true moments of the underlying distribution. Using this result, the authors derive both the one‑dimensional and multivariate forms of the Delta theorem, which allow the propagation of uncertainties through arbitrary smooth transformations g(·).

For the basic set {σ, S, κ}, the Jacobian matrix D = ∂g/∂μ evaluated at the true moments is computed explicitly. The resulting covariance matrix Γ = D Σ Dᵀ / n (where n is the number of events) yields analytic expressions for Var(σ), Var(S), Var(κ) and the covariances among them. In the special case of a symmetric distribution (e.g., a normal distribution), the off‑diagonal elements vanish, confirming that σ, S, and κ become statistically independent. The authors provide compact formulas for the normal case: Var(σ) = σ²/(2n), Var(S) = 6/n, Var(κ) = 24/n.

The same methodology is extended to the products σ·S, κσ², and κσ/S. Here the transformation vector g includes ratios of cumulants, leading to a more involved Jacobian. The derived covariance matrix Π = D Σ Dᵀ / n gives the variances of these composite observables. Again, for a Gaussian underlying distribution the covariances disappear, and the variances reduce to simple multiples of σ² or σ⁴ divided by n.

Higher‑order cumulant ratios, specifically C₆/C₂ and C₈/C₂, are treated in Section 3.3. The authors construct a four‑dimensional vector (μ₂, μ₃, μ₄, μ₆) and its covariance Ω, then apply the Delta theorem to the function g(μ₂, μ₃, μ₄, μ₆) = (μ₆−15μ₂μ₄−10μ₂³+30μ₃²)/μ₂ (and similarly for the eighth‑order ratio). The resulting variance formulas are lengthy polynomials in the normalized central moments m_r = μ_r/σ^r, but they simplify dramatically for a normal distribution: Var(C₆/C₂) = 720 σ⁸ / n and Var(C₈/C₂) = 40320 σ¹² / n.

To validate the analytic error formulas, the authors perform a Monte‑Carlo study using the Skellam distribution, which models the difference of two independent Poisson processes and is a realistic approximation for net‑proton fluctuations. They set the Poisson means to μ₁ = 4.11 and μ₂ = 2.99, values typical of central Au+Au collisions at √s_NN = 200 GeV measured by the STAR experiment. By generating 30 million Skellam‑distributed events, they compute the sample moments and cumulant ratios, and compare the empirical relative errors with the predictions from the Delta‑theorem‑based formulas. The simulation confirms that the analytic expressions accurately describe the scaling of uncertainties with event count, and it highlights the rapid growth of relative errors for higher‑order ratios (C₆/C₂, C₈/C₂), underscoring the need for very large data samples when probing these observables.

In summary, the paper delivers a rigorous, general framework for error estimation in higher‑order moment analyses of heavy‑ion collision data. By grounding the derivations in the Delta theorem, the authors provide closed‑form variance and covariance expressions for all commonly used observables, including individual moments, their products, and high‑order cumulant ratios. The Monte‑Carlo validation demonstrates the practical applicability of the formulas to realistic experimental conditions. This work equips the heavy‑ion community with the statistical tools necessary to quantify uncertainties reliably, thereby strengthening the interpretation of fluctuation measurements in the search for the QCD critical point and the characterization of the phase diagram.