Bayesian Model Selection and Uncertainty Propagation for Beam Energy Scan Heavy-Ion Collisions

This paper employs Bayesian model selection to optimize the energy dependence of parameters in a (3+1)-dimensional hybrid framework for relativistic heavy-ion collisions, subsequently using the resulting posterior distributions to predict various flow observables and quantify systematic uncertainties across different collision systems.

The Gist

Imagine you are a detective trying to figure out how a massive, invisible storm behaves. You can't see the storm itself (because it only lasts for a tiny fraction of a second), but you can see the debris it leaves behind: the particles flying out in all directions.

This paper is about a team of physicists acting as those detectives. They are studying Heavy-Ion Collisions—basically smashing gold atoms together at nearly the speed of light to create a super-hot, super-dense soup of energy called the Quark-Gluon Plasma (QGP). This soup is the state of matter that existed just microseconds after the Big Bang.

Here is the breakdown of their investigation, explained simply:

1. The Problem: Too Many Guesses

To understand this "soup," the scientists use a giant computer simulation. Think of this simulation like a complex recipe for a cake. But here's the catch: they don't know the exact measurements for 20 different ingredients (like temperature, viscosity, how much the ingredients mix, etc.).

In the past, they tried to guess these numbers by looking at data from one specific speed of collision. But they realized that the "recipe" might change depending on how hard they smash the atoms together. If they smash them gently (low energy), the rules might be different than if they smash them violently (high energy).

2. The Tool: The Bayesian "Taste Test"

The authors use a statistical method called Bayesian Model Selection.

• The Analogy: Imagine you are trying to find the perfect cake recipe. You have a "Base Recipe" (Model A) and a "Modified Recipe" that changes an ingredient based on the oven temperature (Model B).
• The Bayes Factor: This is like a super-smart taste test. It doesn't just ask, "Does Model B taste good?" It asks, "Does Model B taste significantly better than Model A, considering that adding a new ingredient makes the recipe more complicated?"
• The Result: The scientists used this test to see if they really needed to change their "ingredients" (model parameters) based on the collision energy. They found that for most ingredients, the simple recipe worked fine. However, for two specific ingredients—related to the size and shape of the initial "hot spots" where the collision starts—they did need to change the recipe depending on the energy. The data proved that the "hot spots" look different at low energies compared to high energies.

3. The Investigation: Adding More Clues

Once they fixed their recipe, they didn't stop there. They decided to test it against more "clues" (experimental data).

• Old Clues: They previously only looked at how many particles were produced.
• New Clues: They added data about the types of particles (pions, protons, etc.), how fast they were moving, and how they flowed in different directions.
• The Surprise: When they added these new clues, the "recipe" had to change again! Specifically, they found that the "switching energy" (the point where the fluid soup turns back into solid particles) needed to be lower than they thought.
• The Domino Effect: Because all the ingredients in the recipe are connected, changing the "switching energy" forced them to adjust other ingredients, like the "stickiness" (viscosity) of the soup. They found that at lower energies, the soup is actually "thicker" (more viscous) than previously thought to match the new data.

4. The Crystal Ball: Predicting the Future

Now that they have a refined, high-precision recipe, they used it to make predictions for experiments that haven't happened yet or haven't been fully analyzed.

• Longitudinal Flow Decorrelation: Imagine a long, twisting ribbon of smoke. If the ribbon twists differently at the top than at the bottom, it's "decorrelated." They predicted how this twisting happens in their soup.
• Small Systems: They predicted what happens when you smash smaller things together (like Oxygen or Deuterium nuclei). They found that even in these tiny collisions, the "soup" behaves like a fluid, which is a big deal for physics.
• Particle Fluctuations: They predicted a new way to measure how the "radial flow" (the outward push of the explosion) fluctuates.

5. The Conclusion

The main takeaway is that simplicity is good, but flexibility is necessary.

• They proved that their computer model is robust.
• They showed that you don't need to overcomplicate the model for everything, but you do need to let the initial size of the collision zone change with energy.
• By using this smart statistical "taste test," they narrowed down the uncertainty of their predictions. Now, when the experimentalists at the RHIC collider (Relativistic Heavy Ion Collider) run their next experiments, the physicists have a much sharper, more accurate map to compare their results against.

In a nutshell: They took a blurry, guess-heavy map of the early universe's "soup," used a smart statistical filter to sharpen it, and now they have a high-definition guide to help future explorers understand how matter behaves under the most extreme conditions in the universe.

Technical Summary

1. Problem Statement

The primary challenge in relativistic heavy-ion collision physics is extracting the properties of the Quark-Gluon Plasma (QGP), such as transport coefficients (shear and bulk viscosities) and the equation of state, from experimental data. This is an inverse problem requiring computationally intensive event-by-event simulations.

• Context: The Relativistic Heavy-Ion Collider (RHIC) Beam Energy Scan (BES) program explores the QGP at finite net-baryon densities ($\mu_B$) by varying the collision energy ($\sqrt{s_{NN}}$).
• Gap: Previous Bayesian analyses (e.g., Ref. [38]) treated model parameters as constant across all collision energies. However, theoretical expectations and preliminary data suggest that certain initial-state and medium properties may depend on collision energy. Introducing energy dependence increases model complexity, raising the risk of overfitting.
• Goal: The authors aim to determine statistically whether introducing collision energy dependence ($\sqrt{s_{NN}}$-dependence) for specific phenomenological parameters is justified by the data, and to quantify the impact of adding new experimental observables on model constraints.

2. Methodology

The study employs a rigorous Bayesian inference framework combined with Bayesian Model Selection.

A. Theoretical Framework

• Model: A (3+1)-dimensional hybrid model (iEBE-MUSIC) consisting of:
  • Initial State: 3D-Glauber model (dynamical initialization accounting for nuclear thickness).
  • Hydrodynamics: Second-order relativistic viscous hydrodynamics (MUSIC) with a lattice-QCD based equation of state (NEOS-BQS).
  • Particlization: Cooper-Frye freeze-out (iSS) with non-equilibrium corrections ($\delta f$).
  • Hadronic Afterburner: UrQMD transport model.
• Parameters: The model utilizes 20 free parameters (Table I), including initial state geometry ($\sigma_x, \sigma_\eta$), energy loss ($y_{loss}$), and temperature/chemical-potential dependent viscosities ($\eta, \zeta$).
• Emulation: A PCSK Gaussian Process (GP) emulator (from the surmise package) is trained on a 20-dimensional parameter space to act as a fast surrogate for the expensive event-by-event simulations.

B. Bayesian Model Selection

To decide if a parameter should depend on $\sqrt{s_{NN}}$, the authors use the Bayes Factor ($B_{A/B}$):
$$B_{A/B} = \frac{P(y_{exp}|A)}{P(y_{exp}|B)}$$

• Model A (Baseline): Parameters are independent of $\sqrt{s_{NN}}$.
• Model B (Extended): Specific parameters are allowed to vary with $\sqrt{s_{NN}}$.
• Criterion: Using Jeffreys' scale, a negative $\ln(B_{A/B})$ favors the extended model. This method penalizes unnecessary complexity, ensuring that added parameters are only included if they significantly improve the description of the data.

C. Data and Observables

The analysis incorporates a comprehensive dataset from RHIC BES (Au+Au collisions at 7.7, 19.6, and 200 GeV):

• Yields & Mean $p_T$: Six identified hadron species ($\pi, K, p, \bar{p}$).
• Fluctuations: Normalized variance of charged hadron mean $p_T$.
• Flow: Anisotropic flow coefficients ($v_n\{2\}$) and scalar-product differential flow $v_2\{SP\}(p_T)$.
• Progressive Calibration: Three posterior distributions are generated by progressively adding these observables to test constraining power.

D. Uncertainty Propagation

To estimate theoretical uncertainties without running thousands of full simulations, the authors use a Cluster Sampling method:

1. Select high-likelihood samples from the posterior.
2. Apply K-Means clustering to identify distinct regions in parameter space.
3. Run high-statistics simulations on the cluster centers (5 sets) to compute the mean and standard deviation of observables.

3. Key Contributions and Results

A. Identification of $\sqrt{s_{NN}}$-Dependent Parameters

• Initial Geometry: The Bayes factor analysis provides moderate evidence ($\ln(B_{A/B}) \approx -2$) that the initial hotspot transverse size ($\sigma_x$) and longitudinal size ($\sigma_\eta$) depend on collision energy.
  • Physical Interpretation: At lower energies (7.7 GeV), hotspots are larger, indicating a transition from partonic to nucleonic degrees of freedom.
• Viscosity: The data favors a non-trivial dependence of shear viscosity on baryon chemical potential ($\mu_B$). However, a collision-energy-dependent bulk viscosity peak ($\zeta_{max}$) is not strongly supported; its effects are largely compensated by the energy dependence of $\sigma_x$.
• Remnants: The parameter controlling energy loss of collision remnants ($\alpha_{rem}$) shows strong energy dependence, dropping near zero at low energies.

B. Impact of Additional Observables

The study compares three posterior distributions (Posterior 1, 2, and 3) with increasing data constraints:

• Switching Energy Density ($e_{sw}$): Including identified particle yields (Posterior 2) shifts the preferred switching energy density from $e_{sw} \approx 0.31$ GeV/fm$^3$ (Posterior 1) to $e_{sw} \approx 0.16$ GeV/fm$^3$. This lower value is crucial for correctly describing hadron chemistry (particle ratios) at low energies.
• Shear Viscosity ($\eta_0$): The shift to a lower $e_{sw}$ implies a longer hydrodynamic lifetime, requiring a larger specific shear viscosity at $\mu_B=0$ to reproduce the observed anisotropic flow.
• Tension: Including $p_T$-differential flow ($v_2\{SP\}(p_T)$) in Posterior 3 reveals a tension: the model struggles to simultaneously describe the collision energy dependence and the pseudorapidity dependence of flow, suggesting the need for temperature-dependent $\eta/s$ in future models.

C. Model Predictions and Uncertainties

Using the optimized model and cluster sampling, the authors provide predictions for upcoming STAR measurements:

1. Longitudinal Flow Decorrelation ($r_n(\eta)$): The model predicts that $r_2(\eta)$ scales with beam rapidity ($\eta/y_{beam}$) at higher energies but shows deviations at lower energies. The results are sensitive to $e_{sw}$; the lower $e_{sw}$ in Posterior 2 improves agreement with preliminary data.
2. Small Systems (O+O and d+Au):
   • Predictions for $v_n\{2\}$ and $v_n\{4\}$ show that elliptic flow in central d+Au collisions is 30-50% larger than in O+O due to more eccentric initial geometry.
   • The ratio $v_2\{4\}/v_2\{2\}$ indicates significantly larger flow fluctuations in O+O compared to d+Au.
3. Identified Particle $v_0(p_T)$: A new observable probing radial flow fluctuations. The model predicts that $v_0(p_T)$ crosses zero at $p_T = \langle p_T \rangle$. As collision energy decreases, the zero-crossing shifts to lower $p_T$ due to reduced radial flow.

4. Significance

• Methodological Advancement: This work demonstrates the power of Bayesian Model Selection in high-energy nuclear physics. It moves beyond simple parameter fitting to objectively determine the necessary complexity of the physical model, preventing overfitting while capturing essential energy dependencies.
• Physical Insights: It establishes that the initial state geometry (hotspot size) is a critical energy-dependent parameter in the BES regime. It also highlights the interplay between the switching energy density ($e_{sw}$), hadronic chemistry, and shear viscosity.
• Future Guidance: The identified tensions (e.g., between $v_2(\eta)$ and $v_2(p_T)$) suggest that future models must incorporate temperature and $\mu_B$ dependence in the shear viscosity. The provided predictions for $v_0(p_T)$ and small-system flow serve as benchmarks for upcoming RHIC BES data, helping to constrain the QGP phase diagram and transport properties.

In summary, the paper refines the (3+1)D hybrid model for the RHIC BES program, proving that while the model is robust, specific initial-state parameters must be energy-dependent to accurately describe the transition from high-energy partonic matter to lower-energy hadronic matter.