Multiplicity and transverse density as ordering variables of the $Λ/K^0_S$ ratio in oxygen-oxygen collisions at the LHC

main.pdf

Introduction

In collisions of heavy nuclei the baryon-to-meson ratios, such as p/$\pi$ and $\Lambda/\mathrm{K}^{0}_{\mathrm{S}}$, are enhanced at intermediate transverse momentum, $2 \lesssim p_{\mathrm{T}} \lesssim5$ GeV, relative to proton–proton collisions [1]. The enhancement is commonly attributed to collective radial flow, which pushes heavier hadrons to higher $p_{\mathrm{T}}$, and to hadronization by quark recombination [2, 3, 4]. Similar features, together with an enhancement of strange-hadron yields, are observed in high-multiplicity pp and p–Pb collisions [5, 6, 7, 8, 9], where they evolve smoothly with charged-particle multiplicity [10, 11]. Multiplicity has thus become the most widely used variable for comparing small and large systems.

For observables driven by collective expansion, however, the system’s transverse size is expected to matter as well. Elliptic flow has long been discussed in terms of the multiplicity per unit transverse area, $(1/S)\,\mathrm{d}N/\mathrm{d}y$ [12, 13, 14], and in hydrodynamics the mean transverse momentum is governed by the entropy density rather than the total entropy [15]. High-multiplicity p–Pb collisions, being more compact than Pb–Pb collisions of the same multiplicity, were predicted to develop a stronger radial flow [16, 17]; the mean transverse momentum [18] and the blast-wave velocities [6] measured in p–Pb collisions are indeed larger than in Pb–Pb collisions of equal multiplicity. How compact p–Pb collisions are depends on the initial-state geometry, in particular on sub-nucleonic fluctuations [19, 20, 21, 22]. In core–corona approaches, the fraction of hadrons produced in a dense core likewise depends on the local density rather than on the total multiplicity [23, 24, 25, 26].

Collisions of light nuclei at the LHC [27, 28, 29] make it possible to separate the two variables with data. Minimum-bias O–O collisions reach the multiplicities of high-multiplicity p–Pb collisions, but with a transverse overlap region similar to that of peripheral Pb–Pb collisions. CMS has measured a pronounced $\Lambda/\mathrm{K}^{0}_{\mathrm{S}}$ enhancement in minimum-bias O–O collisions at $\sqrt{s_{\mathrm{NN}}}=5.36$ TeV [29], and preliminary ALICE results indicate that the mean transverse momenta of strange hadrons in O–O collisions follow the Pb–Pb rather than the p–Pb trend [30]. The Ne–Ne collisions recorded in the same run, with the strongly deformed ${}^{20}\mathrm{Ne}$ nucleus, add a second light-ion system [31].

In this Letter we ask quantitatively, using only published data, whether the $p_{\mathrm{T}}$ shape of the minimum-bias O–O $\Lambda/\mathrm{K}^{0}_{\mathrm{S}}$ ratio can be obtained from p–Pb or Pb–Pb collisions of equal multiplicity (Sec. 3) or of equal transverse density computed with a Glauber model (Sec. 4). The same procedure is applied to EPOS4 and PYTHIA Angantyr events (Sec. 5), and the Ne– Ne ratio is predicted (§ 1.6). Further details are given in the Supplementary Material, cited below as Sec. S1–S6.

Data and method

All inputs are public HEPData records (list in Sec. S1): the CMS minimum-bias O–O K$^{0}_{\mathrm{S}}$ and $\Lambda$ spectra at $|y| < 1$ and the pp reference at the same energy [29], the O–O $\mathrm{d}N_{\mathrm{ch}}/\mathrm{d}\eta$ in 5% centrality classes [28], and the CMS K$^{0}_{\mathrm{S}}$ and $\Lambda+\bar{\Lambda}$ spectra and blast-wave parameters in classes of the number of reconstructed tracks $N_{\mathrm{trk}}^{\mathrm{offline}}$ ($|\eta| < 2.4$, $p_{\mathrm{T}} > 0.4$ GeV) for p–Pb at 5.02 TeV, pp at 7 TeV and Pb–Pb at 2.76 TeV [6]. All $\Lambda$ yields are corrected for feed-down from $\Xi$ and $\Omega$ decays. The central analysis uses only CMS data, since comparisons across experiments are sensitive to conventions (Sec. S1).

Minimum-bias O–O events are a mixture of multiplicities. If the spectra depended on multiplicity only, the minimum-bias ratio would be

where $c$ runs over twenty 5% O–O centrality classes, $N_c=\kappa(\mathrm{d}N_{\mathrm{ch}}/\mathrm{d}\eta)_c$ is the corresponding $N_{\mathrm{trk}}^{\mathrm{offline}}$ with $\kappa=2.75\pm0.35$, and the reference per-event spectra are interpolated between the measured classes linearly in $(\log N,\log\mathrm{yield})$. The 70–100% classes are extrapolated exponentially; beyond the highest reference class, either a log-linear extrapolation or a saturated shape is used.

Since the reference spectra come from an earlier CMS analysis, each prediction is calibrated with the double ratio

whose denominator is the minimum-bias pp ratio of the earlier analysis, built like Eq. (1) from its six 7 TeV pp classes with their fractions $f_c$ of the minimum-bias sample [6] and multiplicities scaled by $r=0.93$, the measured energy dependence of the charged-particle density [32]. The factor $s_{\mathrm{cal}}$ is 0.94 at 1 GeV and 0.81–0.87 above 1.3 GeV; it changes the tilt defined below by −0.05. The uncalibrated predictions and the calibration with $r=1$ are included in the variations.

Data $D$ and prediction $P$ are compared with $\chi^2(a)=(D-aP)^{\mathrm{T}}(C_D+a^2C_P)^{-1}(D-aP)$, profiled over the free scale $a$, so that only the $p_{\mathrm{T}}$ shape is tested. Since HEPData provides only total systematic uncertainties, the covariances are built by decomposing each systematic uncertainty into normalization-like, linear-in-log $p_{\mathrm{T}}$, low-$p_{\mathrm{T}}$, and uncorrelated components, from fully uncorrelated to almost fully correlated; a “balanced” decomposition is the nominal one (Sec. S3). As a single measure of a shape difference, we fit

by generalized least squares with the same covariances. The tilt $s$ is insensitive to normalization; $s>0$ means the prediction is harder than the data, as expected for a stronger radial boost.

Is multiplicity sufficient?

Figure 1 compares the measured O–O ratio with the predictions. The Pb–Pb-based prediction follows the data over the full $p_{\mathrm{T}}$ range. The prediction based on energy-matched p–Pb data lies below the data for $p_{\mathrm{T}}<2$ GeV and above it for $p_{\mathrm{T}}>2.2$ GeV, increasingly so at high $p_{\mathrm{T}}$, the pattern expected from a stronger radial boost. In the shape test it is disfavored, with $\chi^2/\mathrm{ndf}=29-64/14$ ($p\leq0.01$) for uncorrelated uncertainties and 69–125/14 ($p<3\times10^{-9}$) for the balanced correlation model, whereas the Pb–Pb prediction gives $\chi^2/\mathrm{ndf}=1.0-19.5/14$ over all variations. In terms of the tilt, the p–Pb prediction is harder than the data by $s=+0.35\pm0.04$ ($+0.29$ to $+0.43$ over all variations, $8.7\sigma$), and the Pb–Pb prediction gives $s=+0.01\pm0.04$.

The Pb–Pb reference is corrected from 2.76 to 5.02 TeV with the change of the $\Lambda/K^{0}_{\mathrm{S}}$ shape at matched $\mathrm{d}N_{\mathrm{ch}}/\mathrm{d}\eta$ measured by ALICE [1] [33]. This factor deviates from unity by at most 6.5% (rms 3.4%) and shifts the Pb–Pb tilt by +0.03. Scaled in $\ln\sqrt{s_{\mathrm{NN}}}$ to the remaining step from 5.02 to 5.36 TeV, which also separates the p–Pb reference from the O–O data, it changes the tilts by less than 0.005 (Sec. S1).

Three of the twenty O–O classes (0–15%, 42% of the charged-particle yield) lie beyond the highest p–Pb class ($N_{\mathrm{trk}}^{\mathrm{offline}}\approx230$). Built from the seventeen interpolated classes alone, the p–Pb prediction is still harder than the data, $s=+0.24\pm0.04$, although it omits the most central classes, whose spectra are the hardest; the three extrapolated classes alone give $s=+0.53$ (+0.47 with the saturated shape). The extrapolation therefore increases the tension but does not create it. A prediction from the 7 TeV pp classes is not used in the tests, because most of its weight comes from classes beyond the measured pp multiplicity range (Sec. S2).

Transverse density as an alternative ordering variable

We use a Monte Carlo Glauber model [34, 35] in which each participant and each binary collision is a source of Gamma-distributed weight, with the fraction $f = 0.5$ fixed by the O–O centrality dependence of $\mathrm{d}N_{\mathrm{ch}}/\mathrm{d}\eta$ [28], which is reproduced to 2–3% (rms). The transverse area $S_{\perp} = 4\pi\sqrt{\sigma_x^2\sigma_y^2 - \sigma_{xy}^2}$ is computed from the source weights with a Gaussian nucleon width of 0.5 fm, and $\rho= \langle(\mathrm{d}N_{\mathrm{ch}}/\mathrm{d}\eta)/S_{\perp}\rangle$. Variations of $f$, of the fluctuations and of the nucleon width, and a version with three constituents per nucleon [36] with Gaussian spread $\sigma_{\mathrm{c}} = 0.4$ or 0.55 fm, are used as systematic uncertainties (§ 4.4). At $N_{\mathrm{trk}}^{\mathrm{offline}} \approx130$ we find $\rho= 6.9\ \mathrm{fm}^{-2}$ in p–Pb and $3.9\ \mathrm{fm}^{-2}$ in both peripheral Pb–Pb and O–O collisions. The absolute scale of $\rho$ depends on the source width, but the ratio between p–Pb and Pb–Pb at equal multiplicity is robust (1.8 without and 2.1 with sub-nucleonic structure). With constituents, the density of peripheral Pb–Pb classes becomes nearly independent of multiplicity, so density matching to Pb–Pb is ill-conditioned, and the p–Pb reference provides the more sensitive test.

The blast-wave parameters published by CMS separate p–Pb and Pb–Pb strongly at equal multiplicity (0.58 versus 0.43 at $N_{\mathrm{trk}}^{\mathrm{offline}} \approx130$), but lie close to a common curve as a function of $\rho$ [Fig. 2(a),(b)].

Replacing $N_{\mathrm{c}}$ in Eq. (1) by the O–O class density $\rho_{\mathrm{c}}$, taking the reference spectra per unit $\mathrm{d}N_{\mathrm{ch}}/\mathrm{d}\eta$ at the same $\rho$ and rescaling them to the O–O class yield, gives the density-matched prediction [Fig. 2(c)]. For p–Pb, the shape $\chi^2$ decreases by 25–85 units. With nucleons as sources, the prediction becomes compatible with the data for uncorrelated uncertainties ($p = 0.15$–0.99), but a tension remains for the balanced correlation model ($p = 9 \times10^{-4}$–0.08). With sub-nucleonic structure, which makes high-multiplicity p–Pb collisions more compact, it agrees with the data: $\chi^2/\mathrm{ndf} = 2.9$–4.2/14 for uncorrelated and 11.4–15.9/14 ($p = 0.3$–0.65) for balanced correlations. Only the almost fully correlated decomposition still disfavors it, as it also disfavors the uncalibrated Pb–Pb prediction.

The tilt summarizes this progression (Table 1): density matching removes about 60% of the tilt with nucleons as sources and all of it with sub-nucleonic structure, $s = -0.02 \pm0.04$. In all variations of the procedure, the multiplicity-matched tilt stays above +0.29 and the density-matched tilt with constituents within $\pm0.04$ of zero; the dominant uncertainty of the density picture is the modeling of the p–Pb transverse size (§ 4.3, Table S1).

Prediction	$s$	range
p–Pb, $N_{\mathrm{trk}}^{\mathrm{offline}}$	$+0.35 \pm0.04$	$+0.29 \ldots+0.43$
p–Pb, $\rho$, nucleons	$+0.14 \pm0.04$	$+0.12 \ldots+0.19$
p–Pb, $\rho$, constituents	$-0.02 \pm0.04$	$-0.04 \ldots+0.03$
Pb–Pb, $N_{\mathrm{trk}}^{\mathrm{offline}}$ (energy corr.)	$+0.01 \pm0.04$	$-0.05 \ldots+0.09$

Model comparison

To test whether the density ordering emerges from known dynamics, we apply the same procedure to generator events: EPOS 4.0.3 [37], with core–corona separation, 3+1D viscous hydrodynamics for the core and the UrQMD afterburner [38, 39], and PYTHIA 8.317 Angantyr [40, 41], without hydrodynamic expansion, with the default settings, with rope hadronization and string shoving [42, 43], and with QCD-based color reconnection. Event classes, acceptance, and multiplicity estimators follow the data. Within each model, the minimum-bias O–O ratio is predicted from the model’s own p–Pb classes, matched in $N_{\mathrm{trk}}^{\mathrm{offline}}$ or in $\rho$ with the same Glauber tables as for the data, and the tilt is computed with bootstrap uncertainties (§ 4.5).

EPOS4 reproduces the pattern found in the data (Figure 3): matched in multiplicity, the p–Pb-based prediction is too hard ($s = +0.274 \pm0.011$; data $+0.35 \pm0.04$); matched in density with nucleon sources the tilt is reduced to $+0.128 \pm0.007$ (data $+0.14 \pm0.04$), and with constituent sources it vanishes, $s = -0.020 \pm0.009$ (data $-0.02 \pm0.04$). In EPOS4, the effect is thus related to the geometry of the hydrodynamic core: at equal multiplicity, the compact p–Pb core expands more strongly than the O–O core. Without the hadronic cascade, the tilts are $+0.312 \pm0.012$ and $+0.049 \pm0.009$, so the ordering arises before the hadronic phase. Angantyr, in contrast, shows no sensitivity to the transverse geometry: the tilts are small for both matchings ($-0.064 \pm0.008$ and $-0.022 \pm0.004$ for the default settings; similar for the other variants), and are opposite in sign to the data for multiplicity matching. In absolute terms, neither model describes the measured ratio in detail (§ 4.5, Fig. S2): EPOS4 reproduces its normalization and the position of the maximum but underestimates it by up to 30% below 1.5 GeV, whereas Angantyr underestimates the enhancement by about 40% with the default settings.

Prediction for Ne–Ne collisions

The 2025 LHC light-ion run also delivered Ne–Ne collisions at $\sqrt{s_{\mathrm{NN}}} = 5.36$ TeV. The ${}^{20}\mathrm{Ne}$ nucleus is larger than ${}^{16}\mathrm{O}$ and strongly deformed, so it is natural to ask whether it probes a different transverse density at the same multiplicity. No strange-hadron measurement exists yet for Ne–Ne; we therefore predict its minimum-bias $\Lambda/K^0_{\mathrm{s}}$ with the methods tested above. ${}^{20}\mathrm{Ne}$ is described as a deformed Woods–Saxon distribution with the parameters of TGlauberMC v3.3 [44], with a more strongly deformed and a spherical nucleus as variations. The per-source yield is taken from the O–O calibration, so that the Ne–Ne centrality dependence of $\mathrm{d}N_{\mathrm{ch}}/\mathrm{d}\eta$ is a prediction; it agrees within a few percent with the Ne–Ne/O–O ratio measured by ATLAS [45], and the two deformed parametrizations bracket the data (§ 4.4). At equal multiplicity, the two systems have almost the same transverse density [Fig. 4(a)]: the Ne–Ne/O–O ratio of $\rho$ is 0.97–1.01 for $20 \leq\mathrm{d}N_{\mathrm{ch}}/\mathrm{d}\eta\leq100$ with the deformed nucleus and 0.89–0.96 with the spherical one, because a given multiplicity is reached in more peripheral Ne–Ne than O–O collisions, with overlap regions of similar size.

Each method of Secs. 3 and § 1.4 is applied to the Ne–Ne classes and anchored to the measured O–O ratio,

where $P$ is the prediction of the same method for each system, so that method-specific calibrations cancel. All methods give double ratios between 0.97 and 1.04 over $0.8<p_{\mathrm{T}}<5\ \mathrm{GeV}$ (up to 1.06 for the Pb–Pb reference), with tilts between $+0.036$ and $+0.046$ for the nominal settings and between $-0.03$ and $+0.05$ over all variations. The slight hardening reflects the higher mean multiplicity of Ne–Ne collisions; multiplicity and density matching differ by less than 0.01 in the tilt, because the two systems have the same density at equal multiplicity. The minimum-bias Ne–Ne $\Lambda/K^0_{\mathrm{s}}$ is therefore predicted to agree with O–O within a few percent Figure 4(b). This is a consistency test of the picture rather than a discriminating one: Ne–Ne should show the same tension with p–Pb at equal multiplicity and the same agreement at equal $\rho$ as O–O.

We generated $8.99\times10^5$ minimum-bias Ne–Ne events with EPOS4, using generator input for the 20+20 system provided by the EPOS authors (K. Werner, private communication), and $5\times10^6$ events for each Angantyr variant (Sec. S6). In EPOS4 the Ne–Ne/O–O ratio of $\Lambda/K^0_{\mathrm{s}}$ rises from 0.97 at 0.9 GeV to 1.03–1.07 above 2.5 GeV Figure 4(c), with a tilt of $+0.050\pm0.006$, consistent with the data-driven prediction ($+0.040\pm0.006$ without the hadronic cascade); in Angantyr the ratio is flat (tilts $+0.002$ to $+0.005$, $\pm0.003$). Within EPOS4, the p–Pb-based predictions of the Ne–Ne ratio reproduce the O–O pattern: $s=+0.283\pm0.017$ matched in $N_{\mathrm{trk}}^{\mathrm{offline}}$, $+0.107\pm0.007$ in $\rho$ with nucleon sources, and $-0.040\pm0.008$ with constituents (O–O: $-0.020\pm0.009$), whereas Angantyr again shows no ordering (Table S2).

Model	Ne/O $\mathrm{d}N_{\mathrm{ch}}/\mathrm{d}\eta$	Ne/O tilt	$s$, $N_{\mathrm{trk}}^{\mathrm{offline}}$	$s$, $\rho$
EPOS4	1.268	+0.050(6)	+0.283(17) (+0.274)	−0.040(8) (−0.020)
EPOS4, no UrQMD	1.268	+0.040(6)	+0.326(18) (+0.312)	+0.034(8) (+0.049)
Angantyr default	1.272	+0.005(3)	−0.087(11) (−0.064)	−0.024(4) (−0.022)
Angantyr rope	1.269	+0.003(3)	−0.086(8) (−0.070)	−0.032(4) (−0.032)
Angantyr QCD CR	1.272	+0.002(3)	+0.013(10) (+0.012)	−0.007(4) (−0.012)

Discussion and conclusions

A possible reason why the transverse density orders the $\Lambda/K^0_{\mathrm{s}}$ shape better than multiplicity is the build-up of radial flow. In hydrodynamics, the radial boost is generated by pressure gradients, and the mean transverse momentum is, to good approximation, a function of the entropy density $s\propto(\mathrm{d}N/\mathrm{d}y)/S_{\perp}$ [15], so that the compact p–Pb system at high multiplicity expands more violently than the larger O–O or peripheral Pb–Pb system [16, 17]. Because the heavier A is pushed further in $p_{\mathrm{T}}$ than the $K^0_{\mathrm{s}}$, the $\Lambda/K^0_{\mathrm{s}}$ ratio is a species-differential probe of this difference, complementary to the mean transverse momenta [18] [30]. A similar ordering is expected in core–corona approaches [23, 24, 25, 26, 46]. Our Glauber density is a proxy for both pictures, which the present data cannot distinguish; EPOS4, which combines both, reproduces the ordering.

The conclusions are subject to several limitations. (i) The transverse density is computed with a Glauber model. With nucleons as sources, a tilt of $+0.14 \pm0.04$ remains, and only with sub-nucleonic structure does it vanish; the full agreement is therefore conditional on a compact p–Pb geometry, which is favored by other studies [20, 21, 22] but not determined by the present analysis. (ii) The $p_{\mathrm{T}}$ correlations of the systematic uncertainties had to be modeled; the multiplicity-matched p–Pb prediction is disfavored for all correlation models, but the significance of the residual deviations of the density-matched predictions depends on the model. (iii) A residual energy dependence at high multiplicity cannot be excluded, but even one as large as the full measured 2.76–5.02 TeV change in Pb–Pb would shift the tilts by about 0.03, below their uncertainties of $\pm0.04$ and an order of magnitude below the multiplicity-matched p–Pb tilt; ALICE finds the multiplicity evolution of the $p_{\mathrm{T}}$-differential ratios in pp collisions to be similar at 7 and 13 TeV [9] [47]. (iv) Only the minimum-bias O–O ratio has been measured, which averages over a wide range of multiplicities and relies on the interpolation and, in part, extrapolation of the reference classes. (v) Only the $p_{\mathrm{T}}$ shape is tested.

In summary, using only published data, we find that the $p_{\mathrm{T}}$ shape of the $\Lambda/K^0_{\mathrm{s}}$ ratio in minimum-bias O–O collisions at $\sqrt{s_{\mathrm{NN}}}=5.36$ TeV is not reproduced by p–Pb collisions of the same multiplicity (tilt $+0.35 \pm0.04$), but is consistent with that of peripheral Pb–Pb collisions. In a Glauber model calibrated on O–O data, O–O and peripheral Pb–Pb collisions of equal multiplicity have similar transverse densities, about half of that in p–Pb. Matching the p–Pb reference in density reduces the discrepancy by about 60% with nucleons as sources; it is fully removed only when sub-nucleonic structure is included, so that this result is conditional on a compact geometry of high-multiplicity p–Pb collisions. The published radial-flow velocities of p–Pb and Pb–Pb collisions are close to each other at equal density. Within the uncertainties of the data and of the Glauber modeling, the transverse density therefore orders the $p_{\mathrm{T}}$ shape of $\Lambda/K^0_{\mathrm{s}}$ in these systems better than multiplicity does. EPOS4 shows the same pattern, already before the hadronic cascade, whereas PYTHIA Angantyr shows no dependence on the transverse geometry. For Ne–Ne collisions, the Glauber density at equal multiplicity is that of O–O, and the minimum-bias $\Lambda/K^0_{\mathrm{s}}$ ratio is predicted to agree with the O–O ratio within a few percent, slightly harder at high $p_{\mathrm{T}}$, as EPOS4 also finds. The centrality-differential measurements in O–O and Ne–Ne collisions will test the density ordering class by class and constrain the initial-state geometry of high-multiplicity p–Pb collisions. More generally, comparisons of small and large systems “when compared at the same multiplicity” of $p_{\mathrm{T}}$-differential observables that are sensitive to radial flow may need to account for their different transverse sizes; multiplicity remains a good ordering variable for $p_{\mathrm{T}}$-integrated strangeness enhancement [8].

Data availability

All experimental inputs are public HEPData records, except the pp class fractions taken from Table 1 of Ref. [6]. The analysis code is available from the corresponding author upon reasonable request. Supplementary material (Secs. S1–S6, Figs. S1–S2, Tables S1–S2) accompanies this article.

Funding

This research did not receive any specific grant from funding agencies in the public, commercial, or not-for-profit sectors.

Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Acknowledgments

We thank the CMS and ALICE Collaborations for making their results available on HEPData, and the EPOS authors, in particular K. Werner, for providing the EPOS4 input for the Ne–Ne system and for advice.

supplementary.pdf

Section, equation, figure and table numbers without the prefix “S” refer to the Letter.

Inputs, conventions, calibration and energy dependence

All inputs are taken from HEPData:

• O–O at $\sqrt{s_{\mathrm{NN}}}=5.36$ TeV, minimum bias: $K^0_{\mathrm{S}}$ and $\Lambda$ spectra and their ratio at $|y|<1$, and the pp reference at the same energy [1]; $\mathrm{d}N_{\mathrm{ch}}/\mathrm{d}\eta$ at $|\eta|<0.5$ in 5% centrality classes from 0 to 70% [2].

• p–Pb at 5.02 TeV, pp at 7 TeV and Pb–Pb at 2.76 TeV: $K^0_{\mathrm{S}}$ and $\Lambda+\bar{\Lambda}$ spectra at $|y_{\mathrm{cm}}|<1$ in classes of the number of reconstructed tracks $N_{\mathrm{trk}}^{\mathrm{offline}}$ ($|\eta|<2.4$, $p_{\mathrm{T}}>0.4$ GeV), and the blast-wave parameters obtained by CMS [3].

• pp at 5.02 TeV [4]; ALICE Pb–Pb at 2.76 and 5.02 TeV [5, 6], used for the collision-energy correction; ALICE pp $\mathrm{d}N_{\mathrm{ch}}/\mathrm{d}\eta$ at 5.02 and 13 TeV [7], used for the energy dependence of the pp multiplicity in the calibration.

The only inputs not taken from HEPData are the fractions of the minimum-bias sample in the CMS pp multiplicity classes (0.93, 0.06, $6\times10^{-3}$, …), given in Table 1 of Ref. [3]. The CMS O–O measurement reconstructs $\Lambda$ only (no $\bar{\Lambda}$); the $\Lambda+\bar{\Lambda}$ spectra of the other analyses are divided by two. All $\Lambda$ yields are corrected for feed-down from $\Xi$ and $\Omega$ decays. Comparisons across experiments are sensitive to conventions (rapidity versus pseudorapidity normalization changes $\Lambda/K^0_{\mathrm{S}}$ by up to 60% at 0.7 GeV), so the central analysis uses only CMS data.

Calibration. The reference spectra come from a different CMS analysis (2016) than the O–O data (2026), and the pp reference classes are at 7 TeV. Both analyses measured minimum-bias pp collisions, which is the basis of the calibration of Eq. (2). Because a minimum-bias ratio is a yield-weighted average, it is not equivalent to the reference taken at the mean multiplicity; the denominator of Eq. (2) is therefore built from the six pp classes with their event fractions. The class multiplicities are scaled by $r=N(5.36\ \mathrm{TeV})/N(7\ \mathrm{TeV})=0.93$, so that the calibration does not absorb the genuine multiplicity difference between the two energies. Omitting this scaling ($r=1$) changes the multiplicity-matched p–Pb tilt from +0.352 to +0.334; without any calibration it is +0.405.

Pb–Pb energy correction. The change of the $\Lambda/K^0_{\mathrm{S}}$ shape between 2.76 and 5.02 TeV at matched $\mathrm{d}N_{\mathrm{ch}}/\mathrm{d}\eta$ is taken from ALICE [5, 6] and averaged over $\mathrm{d}N_{\mathrm{ch}}/\mathrm{d}\eta\approx55$ and 205. Normalized to unit mean, this factor deviates from unity by at most 6.5% (rms 3.4% over the $p_{\mathrm{T}}$ range); for $\mathrm{d}N_{\mathrm{ch}}/\mathrm{d}\eta\approx205$ alone it reaches 9% (rms 3.9%), for $\mathrm{d}N_{\mathrm{ch}}/\mathrm{d}\eta\approx55$ 6.1% (rms 3.2%). Applied to the multiplicity-matched Pb–Pb prediction it changes the tilt from −0.014 to +0.014.

5.02 to 5.36 TeV. The remaining step from 5.02 to 5.36 TeV, which also separates the p–Pb reference from the O–O data, is 0.11 of the 2.76–5.02 TeV step in $\ln\sqrt{s_{\mathrm{NN}}}$. Scaling the measured factor accordingly changes the tilts by less than 0.005: p–Pb matched in $N_{\mathrm{trk}}^{\mathrm{offline}}$ +0.352 → +0.357, in $\rho$ with nucleons +0.138 → +0.142 and with constituents −0.021 → −0.016; Pb–Pb +0.014 → +0.017. At matched multiplicity, the multiplicity evolution of the $p_{\mathrm{T}}$-differential ratios in pp collisions at 13 TeV is similar to that at 7 TeV, and the $p_{\mathrm{T}}$-integrated ratios agree at similar multiplicities [8, 9], which suggests that the energy dependence at fixed multiplicity is small.

Extrapolated classes. With the saturated shape beyond the highest p–Pb class the full multiplicity-matched tilt is +0.33; for $\kappa= 2.40$ (3.10), with two (four) extrapolated classes carrying 30% (52%) of the charged-particle yield, it is +0.32 (+0.38).

The pp-based prediction

Figure S1 adds the prediction based on the 7 TeV pp classes, which is shown for completeness but not used in the tests. The CMS pp classes extend only to $N_{\mathrm{trk}}^{\mathrm{offline}} \approx137$, whereas the most central O–O classes correspond to $N_{\mathrm{trk}}^{\mathrm{offline}} \approx370$: six of the twenty O–O classes, which carry two thirds of the charged-particle yield, lie beyond the measured pp range. There the prediction is an extrapolation in multiplicity by up to a factor of 2.7, which amplifies the 2–4% statistical uncertainties of the highest pp classes to 6–17% on the minimum-bias prediction (for p–Pb and Pb–Pb, whose classes reach $N_{\mathrm{trk}}^{\mathrm{offline}} \approx230$, it is 0.6–2% and 1–5%). In addition, the calibration to the 5.36 TeV pp measurement corrects the energy difference at minimum-bias pp multiplicity, but not a possible energy dependence at high multiplicity.

Covariance models and systematic variations

HEPData provides only total systematic uncertainties. Each relative systematic uncertainty is therefore decomposed into a normalization-like component (fully correlated in $p_{\mathrm{T}}$), a component linear in $\log p_{\mathrm{T}}$ (tilt), a component restricted to low $p_{\mathrm{T}}$, and an uncorrelated remainder. Results are given for several decompositions, from fully uncorrelated to almost fully correlated; each component is fully correlated in $p_{\mathrm{T}}$ except the last, and statistical uncertainties are added on the diagonal. The nominal “balanced” decomposition assigns variance fractions of 0.3, 0.3, 0.2 and 0.2 to the normalization, tilt, low-$p_{\mathrm{T}}$ and uncorrelated components.

Table S1 lists the variations considered for the two key tilts. The experimental uncertainties of the O–O and reference ratios enter through the covariance model; the remaining entries are variations of the procedure. The multiplicity-matched tilt stays above +0.29 in every variation, and the density-matched tilt with constituent sources stays within $\pm0.04$ of zero. The density-matched tilt depends mainly on how the transverse size of high-multiplicity p–Pb collisions is modeled, which is the dominant uncertainty of the density picture.

Variation	$N_{\mathrm{trk}}^{\mathrm{offline}\text{-}\mathrm{matched}}$	$\rho$-matched
Nominal (balanced covariance)	$+0.352 \pm0.040$	$-0.021 \pm0.040$ (constituents)
Uncorrelated uncertainties	$+0.329 \pm0.052$	$-0.035 \pm0.052$ (constituents)
$\kappa$, extrapolation, calibration	$+0.29 \ldots+0.43$	$-0.04 \ldots+0.03$ (constituents)
Constituent spread $\sigma_c = 0.55$ fm	—	$-0.036 \pm0.040$
Nucleon sources, $w = 0.5$ fm, $f = 0.5$	—	$+0.138 \pm0.040$
Nucleon sources, $w = 0.4$ fm	—	$+0.097 \pm0.040$
Nucleon sources, $f = 0.4$	—	$+0.132 \pm0.040$

The covariance model uses the uncertainties of the reference ratio in one class. To include the statistical fluctuations of all reference classes as they propagate through the interpolation, and beyond the highest class through the extrapolation, we refit the tilts on 400 bootstrap replicas in which every reference spectrum is fluctuated within its statistical uncertainty. The spread of the refitted tilts is 0.007–0.009 for the p–Pb-based predictions and 0.008–0.015 for the Pb–Pb-based one, with no shift of the central values; added in quadrature, it changes the uncertainties of Table 1 by at most 0.003.

Glauber model

The Monte Carlo Glauber model [10, 11] uses Woods–Saxon ${}^{208}\mathrm{Pb}$, harmonic-oscillator ${}^{16}\mathrm{O}$ and $\sigma_{\mathrm{NN}} = 61.8$, 67.6 and 68.2 mb at 2.76, 5.02 and 5.36 TeV. Event classes are defined with estimators that mimic the experimental ones; the model reproduces the mean $N_{\mathrm{trk}}^{\mathrm{offline}}$ of the CMS Pb–Pb classes and those of p–Pb to within 8%. In the constituent version [12] each nucleon is replaced by three constituents with Gaussian spread $\sigma_c = 0.4$ or 0.55 fm, the constituent–constituent cross section (23.4 or 15.6 mb at 5.02 TeV) is tuned to reproduce $\sigma_{\mathrm{NN}}$, and the O–O centrality dependence of $\mathrm{d}N_{\mathrm{ch}}/\mathrm{d}\eta$ is then best described with $f = 0.8$ (1.5–2.2% rms). It reproduces $\sigma_{\mathrm{inel}} = 2.09$ b for p–Pb and 1.35 b for O–O. Our own blast-wave implementation [13] reproduces the CMS parameters of Figure 2(a) when fitted with statistical uncertainties only, as in Ref. [3].

Ne–Ne. ${}^{20}\mathrm{Ne}$ is implemented as a deformed Woods–Saxon distribution with random orientation, with the parameters of TGlauberMC v3.3 [14] ($R = 2.724$ fm, $a = 0.498$ fm, $\beta_2 = 0.490$, $\beta_3 = 0.216$, $\beta_4 = 0.306$); a more strongly deformed nucleus ($\beta_2 = 0.721$) and a spherical two-parameter Fermi distribution are used as variations. Since Ne–Ne and O–O are measured at the same energy, the per-source yield is taken from the O–O calibration, so that the Ne–Ne centrality dependence of $\mathrm{d}N_{\mathrm{ch}}/\mathrm{d}\eta$ is a prediction. We test it against the Ne–Ne/O–O ratio of $\mathrm{d}N_{\mathrm{ch}}/\mathrm{d}\eta$ measured by ATLAS in 5% centrality classes from 0 to 80% [15], in which most of the centrality-determination and normalization uncertainties cancel. For the nominal nucleus the model exceeds the measured ratio by 3% on average (rms deviation 3.4%, 5.7% in the 0–5% class) with nucleon sources and by 4.5% with three constituents; the more deformed nucleus lies 3–4% below the data, and the spherical nucleus agrees best (rms 2.5% and 1.7%). The two deformed parameterizations thus bracket the data. The ATLAS and CMS O–O measurements themselves differ by $-5\%$ (0–5%) to $+27\%$ (65–70%), so absolute values from the two experiments are not combined: the Ne–Ne class multiplicities are the CMS O–O values multiplied by the ATLAS Ne–Ne/O–O ratio. The mean minimum-bias $\mathrm{d}N_{\mathrm{ch}}/\mathrm{d}\eta$ is then 47.0 for Ne–Ne and 39.8 for O–O with the same construction.

Generators

EPOS 4.0.3 [16] is run in its full configuration (core–corona separation, $3+1$D viscous hydrodynamics for the core, UrQMD [17, 18]), with ten particlization samplings per hydrodynamic evolution. PYTHIA 8.317 Angantyr [19, 20] is run with the default Monash settings, with rope hadronization and string shoving [21, 22] (untuned, as stated in the PYTHIA manual), and with QCD-based color reconnection. We generated $9 \times10^{5}$ minimum-bias O–O events at 5.36 TeV and $9 \times10^{5}$ p–Pb events at 5.02 TeV (half of them restricted to $b < 4$ fm to populate high multiplicities) with EPOS, and $5 \times10^{6}$ events per system and variant with Angantyr. Acceptance ($|y_{\mathrm{cm}}| < 1$, $\Lambda$ only, $\Xi$ and $\Omega$ stable so that no weak feed-down enters) and multiplicity estimators (forward transverse energy for O–O, the laboratory-frame $N_{\mathrm{trk}}^{\mathrm{offline}}$ for p–Pb) follow the data. Tilt uncertainties are obtained by bootstrap over independent hydrodynamic evolutions or event blocks. Reference classes are interpolated by a weighted quadratic fit in $\log N_{\mathrm{trk}}^{\mathrm{offline}}$ or $\log\rho$; linear interpolation changes the tilts by less than 0.02. With the larger constituent spread the EPOS4 density-matched tilt is $-0.040 \pm0.009$.

Three of the twenty O–O classes in EPOS (two in Angantyr) lie beyond the most central p–Pb class with sufficient statistics, so the multiplicity-matched prediction is extrapolated there and fluctuates above $p_{\mathrm{T}} \approx2.8$ GeV. Restricting the tilt fit to $p_{\mathrm{T}} < 2.8$ GeV changes the EPOS values to $+0.20$ (multiplicity) and $-0.04$ (density, constituent sources); the conclusions are unchanged. The Angantyr variant with QCD-based color reconnection gives tilts consistent with zero for all matchings ($+0.012 \pm0.009$ for multiplicity, $-0.012 \pm0.004$ for density matching).

Absolute comparison. In absolute terms neither model describes the measured O–O ratio in detail (Figure S2). EPOS4 reproduces its normalization (best-fit scale 1.02) and the position of the maximum, but underestimates the ratio below 1.5 GeV by up to 30% and overestimates it by about 10% at 3–3.6 GeV (shape $\chi^2/\mathrm{ndf}=67/14$ with uncorrelated data uncertainties). Without the hadronic cascade the low-$p_{\mathrm{T}}$ part is closer to the data (34/14). Default Angantyr underestimates the ratio by about 40% at the maximum and has almost no intermediate-$p_{\mathrm{T}}$ enhancement (54/14); rope hadronization raises it by only about 10% (65/14). The QCD color-reconnection variant gives the poorest description of the measured shape ($\chi^2/\mathrm{ndf}=143/14$): it overestimates the ratio below 1.5 GeV and underestimates it by about 20% above 2.5 GeV.

Ne–Ne generator results

The EPOS4 Ne–Ne sample ($8.99 \times10^{5}$ minimum-bias events at 5.36 TeV, configuration of § 4.5) uses generator input for the 20+20 system provided by the EPOS authors (K. Werner, private communication); two of 1800 jobs, which stopped inside the generator, were discarded. Angantyr: $5 \times10^{6}$ events for each of the three variants. Both generators describe ${}^{20}\mathrm{Ne}$ as a spherical Woods–Saxon nucleus. Both overestimate the Ne–Ne/O–O ratio of $\mathrm{d}N_{\mathrm{ch}}/\mathrm{d}\eta$ in the 0–5% class (1.27 for EPOS4 and all Angantyr variants, ATLAS 1.21), similar to our Glauber model with the deformed nucleus. In Angantyr the Ne–Ne/O–O ratio of $\Lambda/K^{0}_{\mathrm{S}}$ is flat, consistent with the weak multiplicity dependence of its $\Lambda/K^{0}_{\mathrm{S}}$ shape; the slight hardening in EPOS4 is consistent with the higher Ne–Ne multiplicity combined with a multiplicity-dependent radial boost, which also underlies the data-driven prediction. Table S2 summarizes the results.