Deep learning predicted elliptic flow of identified particles in heavy-ion collisions at the RHIC and LHC energies

Recent developments on a deep learning feed-forward network for estimating elliptic flow ($v_2$) coefficients in heavy-ion collisions have shown us the prediction power of this technique. The success of the model is mainly the estimation of $v_2$ from final state particle kinematic information and learning the centrality and the transverse momentum ($p_{\rm T}$) dependence of $v_2$. The deep learning model is trained with Pb-Pb collisions at $\sqrt{s_{\rm NN}} = 5.02$ TeV minimum bias events simulated with a multiphase transport model (AMPT). We extend this work to estimate $v_2$ for light-flavor identified particles such as $\pi^{\pm}$, $\rm K^{\pm}$, and $\rm p+\bar{p}$ in heavy-ion collisions at RHIC and LHC energies. The number of constituent quark (NCQ) scaling is also shown. The evolution of $p_{\rm T}$-crossing point of $v_2(p_{\rm T})$, depicting a change in meson-baryon elliptic flow at intermediate-$p_{\rm T}$, is studied for various collision systems and energies. The model is further evaluated by training it for different $p_{\rm T}$ regions. These results are compared with the available experimental data wherever possible.


I. INTRODUCTION
Relativistic heavy-ion collisions have been studied extensively for decades in experiments at the Large Hadron Collider (LHC), CERN, Switzerland, and Relativistic Heavy Ion collider (RHIC), BNL, USA. The formation of a deconfined thermalized medium of quarks and gluons has already been established in such collisions [1]. This medium of hot and dense state of strongly interacting matter is known as the quark-gluon plasma (QGP). Direct probes for studying the properties of QGP are not available owing to the behavior of the strongly interacting matter. However, the signatures of the formation of QGP could be studied using various indirect effects such as jet quenching, strangeness enhancement, and quarkonia suppression, just to name a few. Another crucial observable, which is widely studied to investigate the properties of QGP in heavy-ion collisions, is the transverse collective flow [2]. This transverse collective flow is anisotropic and depends on the equation of state and transport coefficients of the system. Due to the almond-shaped, elliptical nuclear overlap region in noncentral heavy-ion collisions, the initial state has finite spatial anisotropy, further developing the final state momentum anisotropy of the emitted particles. This momentum anisotropy could be expressed as the coefficients of the Fourier expansion of the azimuthal momentum distribution of the produced particles. Finite azimuthal anisotropy has been well observed in heavy-ion collision experiments so far [3][4][5][6]. The second-order flow coefficient, also known as the elliptic flow (v 2 ), which is believed to be driven mainly by the geometry of the distributed nucleons in the nuclear overlap region during a noncentral heavy-ion collision, has the dominant contribution to the overall azimuthal anisotropy.
Different hydrodynamic model calculations suggest that elliptic flow builds up in the early partonic phases of the QGP and evolves through hadronic rescatterings in the hadronic phase. Theoretically, azimuthal anisotropy is thus understood to have consequences from (ideal) fluid dynamics applied to the QGP phase and kinetic descriptions applied to the microscopic hadron cascade phase [7][8][9][10][11]. Therefore, to understand the interplay of partonic and hadronic phases in the evolution of collective flow, it becomes necessary to study the elliptic flow for different identified particles. The competing effects of radial (symmetric) flow and hadronic rescatterings could also be studied by estimating the elliptic flow of identified particles.
Another exciting thing that has been observed in heavy-ion collisions is that v 2 of baryons is larger than that of mesons in the intermediate-p T range (2.0 p T 5.0 GeV/c). This is usually attributed to the particle production by constituent quark coalescence [12]. According to this theory, the development of hydrodynamic flow happens in the early deconfined partonic phase and is subsequently transferred to the hadrons by hadronization through the quark coalescence mechanism [13][14][15]. In an ideal case, this behavior should lead to the number of constituent-quark (NCQ) scaling of the observed flow for baryons and mesons and, thus, strongly hints at the appearance of collectivity at the early deconfined par-tonic level. In Au-Au collisions at √ s NN = 200 GeV, the scaling holds only approximately [16]. However, in Pb-Pb collisions at √ s NN = 2.76 and 5.02 TeV, the scaling holds only up to a level of ±20% [4,17,18].
There are several methods to estimate the flow coefficients, such as the complex reaction plane identification method [19], the cumulant method [20], the Lee-Yang zeroes method [21], and the principal component analysis method [22][23][24][25]. For the first time, it has been shown that deep learning models, based on machine learning (ML) framework, could also be used to estimate the flow coefficients, mainly the results of elliptic flow have been emphasized [26].
The motivation of the present study is to prepare a deep learning framework for the estimation of elliptic flow from final state particle kinematic information. This is an attempt to see whether a machine learning model can learn from final state particle correlations to predict any physical observable of interest. Since the elliptic flow is influenced by several factors such as the collision centrality, energy, system size, particle mass, particle species, and transverse momentum etc., it will be interesting to see how good the deep learning model is to preserve these dependencies. The trained model is found to be faster and more efficient in computing the flow coefficient than some of the traditional methods such as the multi-particle cumulant method.
The proposed deep learning model successfully learns and preserves the centrality and transverse momentum (p T ) dependence of v 2 . It is also shown to perform steadily and accurately when confronted with a data set with additional noise. Here, in this work, we estimate the elliptic flow coefficient for identified stable light hadrons like π ± , K ± , and p+p in Au-Au collisions at √ s NN = 200 GeV, Pb-Pb collisions at √ s NN = 2.76, and 5.02 TeV, and Xe-Xe collisions at √ s NN = 5.44 TeV for (0-10)%, (40)(41)(42)(43)(44)(45)(46)(47)(48)(49)(50)%, and (60-70)% centrality classes by using the deep learning estimator mentioned in Ref. [26]. The NCQ scaling is also studied for these collision systems. The evolution of the p T -crossing point of v 2 (p T ), depicting a change in baryon-meson elliptic flow at intermediatep T , is studied for the above-mentioned collision systems at different energies. The prediction capability of the model is further evaluated by training it with particles in different p T regions.
The paper is organized as follows. We briefly introduce the event generation using a multiphase transport model (AMPT) and the target observable, elliptic flow in Section II. The proposed deep learning estimator is described in Section III along with some performance plots of the model. In Section IV, we describe the results and conclude with a summary in Section V.

II. EVENT GENERATION AND TARGET OBSERVABLE
In order to train the deep learning-based ML algorithm to estimate the elliptic flow, we have used a multiphase transport (AMPT) model to simulate the dataset, a short description of which is discussed in this section, along with the elliptic flow.

A. AMPT model
AMPT is a Monte Carlo-based event simulator that is used to generate ultra-relativistic p-A and A-A collisions at RHIC and LHC energies [27]. AMPT has four components, namely, initialization of collisions by the heavyion jet interaction generator model (HIJING) [28], parton transport by Zhang's Parton Cascade model (ZPC) [29], hadronization of the partons performed by spatial coalescence mechanism in string melting version of AMPT and Lund string fragmentation model in the default version of AMPT [30,31], and finally, the hadron transport using a relativistic transport model (ART) [32,33]. AMPT model is extensively used to study the properties of hot and dense medium formed in relativistic heavyion collisions, and the elliptic flow coefficient is one of them. String meting mode of AMPT describes the elliptic flow well in the intermediate transverse momentum region with quark coalescence mechanism for hadronisation [34][35][36]. Hence, in this study, we have incorporated the string melting mode with AMPT version 2.26t9b. AMPT settings used in this work are the same as reported in Refs. [37][38][39][40]. We have used the impact parameter slicing definition of the centrality, and the corresponding impact parameter values have been taken from Ref. [41].

B. Elliptic flow
Anisotropic flow is one of the key observables in the evolution of QGP medium in noncentral relativistic heavy-ion collisions. The pressure gradient formed in the hot and dense medium due to the initial state spatial anisotropy can transform into final state momentum azimuthal anisotropy. The azimuthal anisotropy depends upon the equation of state as well as the transport coefficients of the medium formed. It can be quantified by the coefficients of Fourier expansion of the azimuthal momentum distribution of particles, given by [42]: Here, v n = cos[n(φ − ψ n )] denotes n th order anisotropic flow coefficient, φ is the azimuthal angle, and ψ n is the n th harmonic symmetry plane angle. The second-order anisotropic flow coefficient (v 2 ) is called the elliptic flow and has a major contribution from the initial eccentricity. In order to calculate the elliptic flow event-by-event, we have used the event plane method [43]. Here, we have set ψ R = 0, which makes v 2 = cos(2φ) . The average is taken over all charged particles in an event.

III. DEEP LEARNING ESTIMATOR
ML is a branch of artificial intelligence that learns correlations from the data to map the input and output observables. Boosted decision trees, deep neural networks, convolutional neural networks, etc., can learn the suitable correlations from the data through model training and help us map complex non-linear observables whose mapping function can not be trivially written in a sequential algorithm. Heavy-ion collisions at LHC energies produce thousands of final state hadrons, each carrying specific kinematic information and evolving through the different physical processes, which make the system complex, non-linear, probing the underlying dynamics challenging. In order to extract the desired physics information from such complex systems, ML has been used by the collider physics community for more than a decade now [44][45][46][47][48][49][50]. In this study, we apply a deep neural network (DNN) based ML algorithm on reckoning the elliptic flow coefficient, which takes the p T , mass and log (s NN /s 0 ) weighted (η − φ) binned values as different features of the input to the model as described in detail in Ref. [26]. The charged particle tracks in |η| < 0.8, with p T ≥ 0.15 GeV/c and φ ∈ [0, 2π] are considered for the training. Also, this is an inverse problem of supervised regression kind, thus the output labels of the DNN model, which are the elliptic flow coefficients, are estimated using the event plane method from the AMPT simulation on an event-by-event basis. The model is trained with minimum bias Pb-Pb collisions at √ s NN = 5.02 TeV. The input data set is divided by the number of events in 8:1:1 for making the training, validation, and testing data set, respectively. The DNN architecture used in this paper for the regression problem consists of one input layer, four dense hidden layers, and one output layer, as shown in Fig. 1. The input layer with 3072 features are mapped to the output layer via four hidden layers, one after another, each having 128, 256, 256, and 256 nodes, respectively. All the four hidden dense layers use the rectified linear unit as the activation function [51], and the output layer has a single node with a linear activation function. The DNN model uses the adam algorithm as the optimizer [52] with mean squared error as the loss function, the details of which can be found in Ref. [26]. We do not implement any L2 regularization or dropout as these seem to hamper the performance of the DNN, as shown in Fig. 2. The mean absolute error (MAE) of elliptic flow (∆v 2 ), defined in Eq. 2, decreases for an increase in λ and decrease in P values; however, ∆v 2 is found to be the smallest for the case with no L2 regularization and dropout. λ and P are the hyperparameters in L2 regularization and dropout, respectively.   Further quality assurance checks on the DNN model can be found in Ref. [26], which suggests that the model can successfully be used not only for different centralities but also in different energy data and gives reasonable results for p T dependent v 2 while having trained with minimum bias events in Pb-Pb collisions at √ s NN = 5.02 TeV data, event-by-event.

IV. RESULTS AND DISCUSSIONS
In this section, we begin with the comparison of the elliptic flow of π ± , K ± , and p +p predicted from the DNN estimator with the true values from AMPT simu-lation and ALICE results. We proceed to show the centrality, energy, and transverse momentum dependence of v 2 , and also study whether the scaling properties of v 2 are encoded in the neural network. The evolution of the crossing point in p T for baryon-meson v 2 separation in intermediate-p T and the effect of transverse momentumdependent training are also also described in this section.
A. Comparing DNN predictions to AMPT and experimental data Figure 4 shows the elliptic flow, v 2 (p T ) for identified π ± , K ± , and p +p in (40-50)% central Pb-Pb collisions at √ s NN = 5.02 TeV. To compare the AMPT and DNN outcomes with ALICE results [18], tracks with p T > 0.5 GeV/c in midrapidity, |y| < 0.5 are considered for this case. The magnitude of v 2 (p T ) increases with increasing p T for all three particle species till it attains a maximum value around p T ≈ 2.0 GeV/c, and then it starts to decrease beyond this point. The values of v 2 (p T ) from AMPT obtained in this region (i.e. p T 2.0−3.0 GeV/c) is comparable in magnitude with ALICE results for the individual particle cases. However, beyond this transverse momentum value, AMPT fails to describe the data as v 2 (p T ) falls faster with increasing p T since fragmentation becomes the preferred mode of hadronization at high-p T . It is interesting to note that DNN predictions agree with AMPT values quite nicely up to p T 4.0−6.0 GeV/c. Beyond this p T , the values from DNN start to differ from the AMPT values. The DNN model is trained with Pb-Pb collisions at √ s NN = 5.02 TeV minimum bias events with the selection of all charged particles having p T > 0.15 GeV/c in pseudorapidity, |η| < 0.8 simulated with AMPT. It is seen that the statistics in each p Tbin for all three particle types decreases to about a few 100 counts for p T 6.0 GeV/c, which presents a few instances to the DNN model during the training process. For this reason, the mismatch between DNN and AMPT comes into the picture for p T 6.0 GeV/c. For the above-mentioned reasons, this analysis is kept in the range p T < 3.0 GeV/c for the rest of the plots where the Monte Carlo and data are comparable.
B. Centrality and energy dependence of elliptic flow Figure 5 represents the centrality dependence of v 2 (p T ) for π ± , K ± , and p +p, and all charged hadrons (h ± ) in Au-Au collisions at then it starts to decrease. This is well understood from the fact that initial geometrical anisotropy in the nuclear overlap region keeps increasing for peripheral collisions. However, for more peripheral collisions, i.e. (60-70)%, v 2 (p T ) gets reduced as the smaller size and shorter lifetime of the system do not allow the generation of large v 2 . For p T 1.5 GeV/c, there is a mass ordering of v 2 (p T ) for different particles, meaning lighter particles have more v 2 than the heavier ones, here This is understood as an effect of strong radial flow, which imposes an azimuthally symmetric velocity boost to all particles along with the anisotropic flow in the medium. In the intermediate- The baryon-meson v 2 separation could arise from the constituent quark coalescence mechanism embedded in the AMPT string melting mode [57]. The DNN estimator, trained at the LHC energy, in Pb-Pb collisions at √ s NN = 5.02 TeV minimum bias data set, can predict the identified particle v 2 for different collision systems at different collision energies. The mass ordering at low-p T and the grouping of particles based on their constituent quarks in the intermediate-p T is learned and preserved by the DNN estimator. Figure 6 shows the ratio of v 2 (p T ) estimated with the DNN model over the AMPT values for π ± , K ± , and p+p, and all charged hadrons (h ± ). For the (40-50)% centrality class, the DNN predictions for all particle types are in better agreement with the AMPT values as compared to (0-10)% and (60-70)% centrality classes. In Ref. [26], while dealing with the centrality dependence of v 2 for all charged hadrons, one can observe that the predictions for the most central and peripheral cases have larger statistical uncertainty than the mid-central cases. Here, Fig. 6 also suggests similar behavior. We argue that this is entirely statistics driven as the number of events in the extreme centrality bins is less in a minimum bias data set; hence, with the minimum bias training, the DNN model does not get enough examples to learn from these extreme domains. A similar argument could also be given for the slight mismatch of proton v 2 at lowp T and pion v 2 at intermediate-p T , where the respective statistics is comparably less. The solution to this problem could come from various collaborative learning models [58]. The same DNN model could be trained separately with centrality-wise data sets, unlike the minimum bias training to obtain a common set of model parameters. This could possibly take out the training sample bias due to the different number of events in different centrality bins, leading to an unequal level of information passing to the DNN training in different centrality bins.

C. Constituent quark number scaling
The study of centrality and p T dependence of v 2 for various identified particles seems incomplete without exploring its constituent quark number (NCQ) scaling behavior. By recalling our understanding of the specific mass ordering of v 2 (p T ) at lower p T regime, if the mass ordering is driven by a hydrodynamic pressure gradient, then v 2 (KE T ) at lower-p T should scale to a common set of values for all the particles, which is indeed observed in experiments [59]. Here, KE T (= m T − m 0 ) is the transverse kinetic energy, and m T (= p 2 T + m 2 0 ) and m 0 are the transverse mass and rest mass of the particle, respectively. At the intermediate-p T regime, v 2 (KE T ) splits into two branches separately grouped by baryons and mesons. In this domain, the constituent quarks of the particles play a more substantial role in generating v 2 in hadrons than its mass. When the NCQ-scaling is applied to v 2 (KE T ), it shows a better scaling behaviour than v 2 (p T ) at lower-p T [16,59]. Figure 7 shows the centrality dependence of v 2 /n q as a function of KE T /n q . Here, the number of constituent quarks, n q = 2 for mesons and n q = 3 for baryons. One can see that, at lower-KE T /n q , the scaling is valid for all the particle types. However, the proton seems to break the scaling at higher-KE T /n q . The violation of KE T scaling is observed in experiments, so the results from AMPT are in line [16,60]. Here, it is interesting to note that the DNN can learn particle species-dependent scaling behavior as it closely agrees with the AMPT values for all energy and centrality classes.

D. Evolution of the crossing point
The separation of baryon-meson elliptic flow at the intermediate-p T is usually attributed to the quark coalescence/recombination mechanism of hadronization. It suggests an apparent dependence of hadron flow on the number of constituent quarks of the hadron. This leads to a higher flow of baryons than mesons and a relative enhancement of baryon yield over meson yield in the intermediate-p T regime. Now, by looking into the crossing point in p T to know precisely where the separation of baryon-meson flow occurs, one can confirm the coalescence picture if the crossing occurs at somewhat higher transverse momenta for higher centralities. Figure 8 shows the centrality dependence of the p T -crossing point of pion and proton elliptic flow for various collision systems at various energies. The vertical error bars represent the difference between the pion-proton and kaon-proton p T -crossing points. One can observe that the crossing point gradually shifts to higher-p T as one move from peripheral to central collisions for all collision energies. The shift of the curve to higher-p T with increasing energy could hint towards producing a denser partonic medium and the effect of increased radial flow with increasing collision energy. Here, the AMPT model qualitatively reproduces the trend of p T -crossing curve as seen in the ALICE results. A similar observation is made from the DNN predictions, which also closely match the AMPT curves for the respective collision systems and energies. The DNN model not only learns the splitting of the elliptic flow into two branches grouped separately by baryons and mesons but also preserves the same p T -crossing point of the separation, thus adding a more quantitative hold on the estimation of elliptic flow for identified charged hadrons.

E. Effect of transverse momentum dependent training on DNN
In Figure 9, we present the predictions for elliptic flow as a function of transverse momentum from three differ-ent DNN models, each trained with a different set of inputs based on specific p T -cuts. These models are trained on the minimum bias Pb-Pb collisions at √ s NN = 5.02 TeV. For this investigation, the same model architecture with similar hyperparameter settings is used for training. The only difference would be the obtained weights and biases of the network for different p T -cuts.
The first row of Figure 9 refers to the predictions from the DNN model trained with particle kinematic information for tracks in |η| < 0.8 with p T > 0.15 GeV/c. We name it the unbiased case. It shows the evolution of v 2 (p T ) in (40)(41)(42)(43)(44)(45)(46)(47)(48)(49)(50) including the transition behavior around p T | v max 2 . Thus the ratio of DNN to AMPT remains fairly close to unity for the entire p T -range (shown up to 7 GeV/c) for all the particle species. So far, all the above-discussed results are obtained from the DNN model with this unbiased training. Now, we study the effect of training the DNN model in different transverse momentum regions. The aim is to obtain two different DNN models trained with and without the transition point of p T | v max 2 and then test their applicability to the estimation of v 2 (p T ) curve for the full range of p T . We can obtain two groups of particles in an event based on the transverse momentum cuts applied to the tracks and we name them as low-p T cut group (p Train T ≤ 3.0 GeV/c) and high-p T cut group (p Train T > 3.0 GeV/c). The estimation of elliptic flow is also performed using only these selected particles for the respective groups. The low-p T cut group includes the p T | v max 2 , thus, the transition behaviour of v 2 (p T ) around p T | v max 2 . However, the high-p T cut group does not include this transition behavior directly from particle kinematics. After the training of these two separate DNN models with low-p T and high-p T particles from Pb-Pb collisions at √ s NN = 5.02 TeV minimum bias events, we apply each model to predict the entire v 2 (p T ) curve for the full range of p T (shown up to 7 GeV/c) for the (40-50)% central Pb-Pb collisions at √ s NN = 5.02 TeV.
From the second and third row of Figure 9, it is quite interesting to see that both the models, although being trained with two disjoint sets of particle groups, can successfully reproduce the entire v 2 (p T ) curve up to a reasonable extent in their respective cases. Somehow the models learn and complement the missing information of the v 2 (p T ) curve. For the low-p T cut group, the missing factor in training is the high-p T particles, and for the high-p T cut group, it is the low-p T particles, and also the absence of the p T | v max 2 . As the information encoded in the network behaves like a generalized global curve of v 2 applicable to the entire p T domain, the DNN models can extrapolate any regions of interest regardless of their training. This global curve learned by the DNN models represents the dependence of v 2 on transverse momentum and also its transition behavior around the p T | v max 2 .
Although the mathematical function for such a curve is not available in the literature, the DNN can map this function from the existing input-output correlations irrespective of its training domain.
Another observation is the behavior of elliptic flow for different particle species. When the DNN model trained with low-p T particles is applied to predict v 2 for pion and kaon, a slight deviation is observed between DNN and AMPT starting from intermediate to high-p T . This is visible in the bottom ratio plots in the second row of Fig. 9. It hints that the elliptic flow for pion and kaon primarily originates from low-to intermediate-p T . However, in the same figure, the prediction for proton deviates even at very low-p T , although the model includes this low-p T region in training. It suggests that the contribution to elliptic flow for proton mainly originates from high-p T . Thus by removing this region in the low-p T cut DNN model, the accuracy is reduced for the entire curve for proton v 2 .
This hypothesis can be verified with the DNN model trained with high-p T particles. One can see from the third row of Fig. 9 that the prediction of elliptic flow for proton is now more accurate for the entire p T range as this DNN model includes high-p T particles. This observation supports our earlier statement of DNN captures a global curve of v 2 (p T ) and proton v 2 has a dominant contribution from high-p T . However, in the same figure, pion and kaon elliptic flow deviate more in the low-p T as the necessary contribution to elliptic flow is absent in this region for the high-p T cut model, and same is visible in the ratio plots.

V. SUMMARY
To sum everything up that has been stated so far, this work demonstrates the applicability and accuracy of a DNN-based machine learning model to evaluate the second-order anisotropic flow coefficient (v 2 ) for identified particles from final state particle kinematic information in heavy-ion collisions. The DNN model can eminently estimate v 2 for light-flavor identified particles such as π ± , K ± , and p +p in heavy-ion collisions at RHIC and LHC energies. The prediction accuracy of the DNN model is not only limited to Pb-Pb collisions at √ s NN = 5.02 TeV at which the model is trained but also applicable to other systems such as Pb-Pb collisions at √ s NN = 2.76 TeV, Xe-Xe collisions at √ s NN = 5.44 TeV, and Au-Au collisions at √ s NN = 200 GeV. From the minimum bias training, the model could successfully learn and preserve the centrality, energy, and transverse momentum dependence of elliptic flow for other collision systems at various energies. This becomes a striking feature of DNN-based machine learning, where, while dealing with multihadron production dynamics, the machine is trained with the control parameters like the collision energy, impact parameter (here, the number of participants or the final state charged particle multiplicity), particle mass, and the transverse momentum. In addition, the DNN model retains the constituent quarks number scaling behavior for elliptic flow as a function of transverse kinetic energy for the respective collision systems and energies. The DNN model also quantitatively preserves the p T -crossing point of baryon-meson v 2 curves at intermediate-p T , which scales with centrality and centerof-mass energy of collisions for various systems. Finally, by training the DNN models with different kinematic regions based on specific p T -cuts, and then applying it to the full p T range, we arrive at the conclusion that there exists a global functional dependency of elliptic flow with transverse momentum, as the DNN model trained in a certain kinematic region could describe the overall curve for the full kinematic domain. We also learn that the elliptic flow of pion and kaon originates mainly from the low-p T region, whereas the elliptic flow of the proton is dominated due to the high-p T region.