11institutetext: State Key Laboratory of Particle Detection and Electronics,
University of Science and Technology of China, Hefei, China

# Rescattering effect on the measurement of $K^{\ast 0}$ spin alignment in
heavy-ion collisions

Ziyang Li  Wangmei<EMAIL_ADDRESS>and Zebo<EMAIL_ADDRESS>
###### Abstract

Spin alignment of vector mesons in noncentral relativistic heavy ion
collisions provides a novel probe of the global quark polarization and
hadronization mechanism. In this paper, we discuss the spin alignment of a
short-lived vector meson, namely $K^{\ast 0}$, arising from the hadronic
rescattering process using the UrQMD model. This spin alignment is not related
to global quark polarization but cannot be distinguished from those arising
from global quark polarization in experiments. The spin alignment parameter
$\rho_{00}$ is found to deviate from 1/3 by up to -0.008 (0.03) with respect
to the reaction (production) plane. These deviations are much larger than the
expected signal from all the theoretical models implementing the conventional
global quark polarization mechanism as well as the current experimental
precision, and should be considered seriously when comparing measurements with
theoretical calculations.

##  1 Introduction

Relativistic heavy-ion collisions create a hot, dense medium with partonic
degrees of freedom, called the quark-gluon plasma (QGP). It provides unique
opportunities to study quantum chromodynamics (QCD). In noncentral heavy-ion
collisions, a large orbital angular momentum ($\sim$ $10^{4}~{}\hbar$)
liang2005globally ; Becattini:2007sr  and a strong magnetic field ($\sim$
$10^{14}~{}T$) Kharzeev:2007jp  are also expected. While the magnetic field is
expected to be short-lived, part of the angular momentum is conserved and
could be felt throughout the evolution of the system formed in the collision
PhysRevLett.125.012301 . It is predicted that particles produced in heavy-ion
collisions with large global angular momentum could be polarized globally in
the direction of the orbital angular momentum due to spin-orbit coupling
liang2005globally . In 2017, the STAR Collaboration reported the discovery of
global polarization for $\Lambda$ and $\bar{\Lambda}$ hyperons in Au+Au
collisions adamczyk2017global . According to the flavor-spin wave function,
the polarization of the $\Lambda$($\bar{\Lambda}$) hyperon is carried solely
by the strange quark s($\bar{s}$), indicating the global polarization of the
s($\bar{s}$) quark.

The polarized quarks can also form polarized vector mesons such as
$\phi(1020)$ and $K^{\ast}(892)$. These mesons originate predominantly from
primordial production and are less affected by feed-down contribution compared
to $\Lambda$ and $\bar{\Lambda}$ hyperons. Furthermore, as spin-1 particles,
their daughter’s polar angle distribution is an even function; thus, there is
no local cancellation as associated with spin-1/2 hyperons when integrating
over time and phase space.

The spin state of a vector meson can be described by a 3$\times$3 Hermitian
spin density matrix $\rho$ with unit trace. The diagonal elements
$\rho_{-1-1}$ and $\rho_{11}$ cannot be measured separately. Consequently,
$\rho_{00}$ is the only independent element SCHILLING1970397 . It can be
determined from the angular distribution of the decay products,

$\frac{dN}{d\cos\theta^{*}}\propto(1-\rho_{00})+(3\rho_{00}-1)\cos^{2}\theta^{*},$
(1)

where $\theta^{*}$ is the angle between the polarization direction and the
momentum direction of one of the decay daughters in the rest frame of the
vector meson. $\rho_{00}$ is 1/3 in the absence of spin alignment, and a
deviation of $\rho_{00}$ from 1/3 signals net spin alignment.

The ALICE Collaboration recently published the measurements of vector mesons
$K^{\ast 0}$ and $\phi$ spin alignment in Pb+Pb collisions at 2.76 TeV
PhysRevLett.125.012301 . $\rho_{00}$ values are found to decrease with
decreasing transverse momentum ($p_{T}$), consistent with 1/3 at high $p_{T}$
but less than 1/3 at low $p_{T}$ ($p_{T}<2~{}\text{GeV}/c$) at a level of
3$\sigma$ (2$\sigma$) for $K^{\ast 0}$ ($\phi$), respectively. The values at
low $p_{T}$ are systematically smaller for $K^{\ast 0}$ than for $\phi$, and
with respect to the production plane than that with respect to the event
plane. The difference is on the order of a few percent, although comparable to
the uncertainties. The deviation of $\rho_{00}$ from 1/3 at low $p_{T}$ for
the vector mesons is qualitatively consistent with the expectation from
theoretical models, which attribute it to a polarization of quarks in the
presence of angular momentum in heavy-ion collisions and subsequent
hadronization by the process of recombination. However, the measured spin
alignment is unexpectedly large compared to the global polarization of
$\Lambda$ hyperons ALICE:2019onw . This is therefore puzzling.

The STAR Collaboration measured the spin alignment for $K^{\ast 0}$ and $\phi$
in Au+Au collisions at $\sqrt{s_{\mathrm{NN}}}$ from 11.5 to 200 GeV
STAR:2022spinalignment . $\rho_{00}$ for $K^{\ast 0}$ is largely consistent
with 1/3 ($0.3356\pm 0.0034~{}\textrm{(stat.)}\pm 0.0043~{}\textrm{(syst.)}$
at $\sqrt{s_{\mathrm{NN}}}\leq 54.4~{}\textrm{GeV}$), while $\rho_{00}$ for
$\phi$ is above 1/3 with a significance of $7.4\sigma$ ($0.3512\pm
0.0017~{}\textrm{(stat.)}\pm 0.0017~{}\textrm{(syst.)}$ at
$\sqrt{s_{\mathrm{NN}}}\leq 62~{}\textrm{GeV}$). The large difference between
these two vector mesons cannot be described by any theoretical models with
conventional global quark polarization mechanisms PhysRevD.105.L011901 ;
PhysRevC.98.014910 ; adamczyk2017global ; PhysRevD.104.076016 ;
PhysRevD.102.056013 ; XIA2021136325 ; PhysRevC.97.034917 ; liang2005spin ;
sheng2020can . A recent theoretical model based on $\phi$-meson vector field
coupling to strange quark can qualitatively describe the difference, but more
developments are needed for quantitative description STAR:2022spinalignment ;
sheng2020can . On the other hand, it is important to survey the other possible
sources contributing to the difference between $K^{\ast 0}$ and $\phi$ mesons
spin alignment.

The lifetime of $K^{\ast 0}$ is approximately 4 fm/$c$, which is approximately
10 times shorter than the $\phi$ lifetime and comparable to the time between
chemical freeze-out and kinetic freeze-out, making it sensitive to hadronic
interactions. In experiments, the vector mesons are usually reconstructed via
their strong decay to hadrons. If a $K^{\ast 0}$ decays before chemical
freeze-out, the decay daughters are subject to rescattering by hadrons via
elastic or inelastic collisions. If one or more daughter particles are
rescattered, either its momentum is altered or completely absorbed, and the
parent resonance cannot be reconstructed experimentally. The measurements of
centrality and $p_{T}$ dependence of $K^{\ast 0}$ yield in heavy-ion
collisions show strong evidence STAR:2004bgh ; PhysRevC.84.034909 ;
ALICE_Kstar_rescattering  of the significant rescattering effect on $K^{\ast
0}$ production.

The rescattering probability depends on the particle density profile, which is
related to the geometry of the medium. In noncentral heavy ion collisions, the
rescattering probability as a function of azimuthal angle may be an
anisotropic distribution due to the initial geometry of the collisions. Since
$\rho_{00}$ of vector mesons are measured via the angular distribution of the
daughter particle with respect to the event plane or production plane, it may
be distorted by the rescattering effect for reconstructable resonances.
Consequently, the rescattering effect could produce artificial spin alignment
for resonances, such as $K^{\ast 0}$.

In this paper, the spin alignment of $K^{\ast 0}$ arising from hadronic
rescattering in Au+Au collisions at
$\sqrt{s_{\mathrm{NN}}}=200~{}\textrm{GeV}$ is studied with the
Ultrarelativistic Quantum Molecular Dynamics (UrQMD) model. The analysis
method is introduced in Sec. 2. Sec. 3 includes the results and discussions.
The paper is summarized in Sec. 4.

##  2 Analysis Method

Figure 1: Invariant mass distribution of $K^{*0}\rightarrow K\pi$ for
reconstructable $K^{*0}$ (red solid line) and unreconstructable $K^{*0}$ (blue
dot-dashed line) in 0-20$\%$ centrality. The unreconstructable $K^{*0}$ is
those with one or two daughters experiencing at least one elastic scattering.
The distributions are scaled by their integral.

(a)

(b)

Figure 2: (a) The yield of all (red line) and reconstructable (black line)
$K^{\ast 0}$ as a function of the time when they decayed in Au+Au collisions
at $\sqrt{s_{\mathrm{NN}}}$ = 200 GeV for 0-80$\%$ centrality from the UrQMD.
(b) The ratio of reconstructable and all $K^{\ast}$ as a function of the time
for various centralities.  Figure 3: The yield ratio of reconstructable
$K^{\ast 0}$ over $K^{-}$ in Au+Au collisions from UrQMD divided by that
measured in $p+p$ collisions as a function of $\langle
N_{\textrm{part}}\rangle$, compared with measurements in Au + Au and Cu + Cu
collisions PhysRevC.84.034909 . The center-of-mass energy is
$\sqrt{s_{\mathrm{NN}}}$=200 GeV.

The spin alignment of $K^{\ast 0}$ arising from hadronic rescattering is
investigated with the UrQMD model. UrQMD is a microscopic transport approach
based on the covariant propagation of hadrons and strings. All cross-sections
are calculated by the principle of detailed balance or are fitted to available
data. One important feature of UrQMD is that it has the option to save all the
collision history of particles during the evolution of the system, making it
an ideal tool to study the hadronic rescattering effect.

As mentioned above, $K^{\ast 0}$ is usually reconstructed via the decay of
$K^{\ast 0}\rightarrow K\pi$. In a heavy-ion collision such as Au+Au collision
at $\sqrt{s_{\mathrm{NN}}}$ = 200 GeV, hundreds or even thousands of pions and
kaons are produced. In the experiment, to reconstruct $K^{\ast 0}\rightarrow
K\pi$, one needs to pair all the pion candidates with all the kaon candidates
in the same event and calculate their invariant masses. With this method, the
decay of $K^{\ast 0}\rightarrow K\pi$ can be reconstructed, together with a
huge random combinatorial background. Although the random combinatorial
background can be reproduced with the like-sign or mixed-event technique and
statistically subtracted STAR:2004bgh , the measurement of $K^{\ast 0}$ is
very challenging as the typical signal-to-background ratio is on the order of
1/100. In many heavy-ion collision event generators, only the final particles
are saved. One needs to reconstruct $K^{\ast 0}$ in the same way as in
experiments. Alternatively, one can turn off the decay of $K^{\ast 0}$ and
directly access the information of $K^{\ast 0}$ from the output file. However,
in this case, $K^{\ast 0}$ does not suffer from the rescattering effect, which
is the main focus of this paper.

Thanks to the collision history file of UrQMD, one can easily study the
rescattering effect of resonance. The history file contains each binary
interaction, resonance decay, and string excitation that occurred in the
course of the heavy-ion reaction. First, we find all the $K^{\ast
0}\rightarrow K^{\pm}\pi^{\mp}$ decays with $|y(K^{\ast 0})|<0.5$ from the
collision history file. We then follow the collision history of the decay
daughters. If any of the decay daughters participates in any (elastic or
inelastic) collision with other hadrons, the $K^{\ast 0}$ is categorized as
“unreconstructable”. If neither of the two decay daughters participates in any
collision with other hadrons, the $K^{\ast 0}$ is categorized as
“reconstructable”. Since it has been found that a finite $\eta$ acceptance
will lead to an increase in the observed $\rho_{00}$ LAN2018319 , we apply no
$\eta$ acceptance cut to decay daughters. Figure 1 shows the invariant mass
distribution of $K^{\ast 0}\rightarrow K^{\pm}\pi^{\mp}$ in 0-20$\%$
centrality. The solid red line depicts the invariant mass distribution of
$\pi$ and K from the decay of reconstructable $K^{\ast 0}$. It follows the
Breit-Wigner distribution. The dashed blue line depicts the invariant mass
distribution of $\pi$ and K from $K^{\ast 0}$ with at least one of them
participating in an elastic collision with other hadrons and the momentum
after the collision is used in the invariant mass calculation. The
distribution is significantly broadened. In reality, the distribution should
be even broader as the decay daughter measured by detectors may have
experienced more than one elastic collision. These $K^{\ast 0}$ cannot be
reconstructed in the presence of a large random combinatorial background and
significantly correlated background in heavy ion collisions. Not presented in
Fig. 1 is the $K^{\ast 0}$ with one or two decay daughters absorbed in the
following collisions. It is completely lost. Fig. 2(a) shows the number of all
and reconstructable $K^{\ast 0}$ as a function of the time when they decay to
K and $\pi$ in 0-80% Au+Au collisions at $\sqrt{s_{\mathrm{NN}}}$ = 200 GeV.
The mean of the time distribution for both all and reconstructable $K^{\ast}$
is between 10 and 30 fm/$c$, which means that most $K^{\ast 0}$ decays during
this time. Figure 2(b) shows the ratio of reconstructable over all $K^{\ast
0}$ as a function of time in various centrality classes. It increases with
increasing time and toward peripheral collisions, showing strong evidence of a
significant rescattering effect implemented in the UrQMD model, especially at
low $p_{T}$ and in central collisions. Figure 3 shows the $K^{\ast 0}/K^{-}$
ratio in UrQMD normalized by the ratio measured in $p+p$ collisions at
$\sqrt{s_{\mathrm{NN}}}$ = 200 GeV PhysRevC.84.034909 . The double ratio
decreases with increasing system size, consistent with the experimental
results in Au + Au and Cu + Cu collisions PhysRevC.84.034909 . This comparison
verifies that the rescattering processes implemented in the UrQMD model are
reasonably good.

Figure 4: The $\cos\theta^{*}$ distributions for (a, c) all $K^{\ast}$ and (b,
d) reconstructable $K^{\ast}$ for 0.8 ¡ $p_{T}$ ¡ 1.2 GeV/c in 20-60$\%$
centrality w.r.t. the reaction plane and production plane, respectively. The
solid lines represent the fit to Eq. 1 and the dot-dashed lines represent the
fit to a constant.

The spin alignment parameter $\rho_{00}$ of all and reconstructable $K^{\ast
0}$ is extracted with respect to the reaction plane and production plane,
respectively, according to Eq. 1. The reaction plane is defined by the beam
momentum and the impact parameter, while the production plane is defined by
the beam momentum and $K^{\ast 0}$ momentum. The polarization direction is
chosen to be the direction perpendicular to the planes. Fig. 4 shows an
example of the $\cos\theta^{*}$ distribution with respect to the reaction
plane and production plane for all (a, c) and reconstructable (b, d) $K^{\ast
0}$ at $0.8<p_{T}<1.2~{}\textrm{GeV}/c$ in 20-60% Au+Au collisions at
$\sqrt{s_{\mathrm{NN}}}$ = 200 GeV. The solid lines represent the fit to Eq.
1.

##  3 Results and Discussions

###  3.1 Spin Alignment with Respect to the Reaction Plane

(a)

(b)

Figure 5:  $\rho_{00}$ with respect to the reaction plane of all (open) and
reconstructable (filled) $K^{\ast 0}$ as a function of $p_{T}$ in 20-60$\%$
centrality (a) and that of reconstructable $K^{\ast 0}$ as a function of
centrality in two different $p_{T}$ intervals (b) in Au+Au collisions at
$\sqrt{s_{\mathrm{NN}}}$ = 200 GeV.

We first investigated the spin alignment of all $K^{\ast 0}$. As shown in
Figure 4(a), the $\cos\theta^{*}$ distribution is flat for all $K^{\ast 0}$ in
20-60% Au+Au collisions and the extracted $\rho_{00}=0.3325\pm 0.0009$ is
consistent with 1/3 within uncertainties. The $p_{T}$ dependence is shown in
Fig. 5(a) as open circles. No spin alignment for all $K^{\ast 0}$ is observed,
which proves that the initially produced $K^{\ast 0}$ has no spin alignment in
the UrQMD model as expected.

We then studied the rescattering effect on the $\rho_{00}$ measurement. Figure
4(b) shows that there is less reconstructable $K^{\ast 0}$ toward
$|\cos\theta^{*}|=1$ or $\theta^{*}=0~{}\textrm{or }\pi$. This is because
there are more particles out of the reaction plane than in the reaction plane
in noncentral heavy ion collisions. The decay daughters of $K^{\ast 0}$ have a
larger rescattering probability along the polarization direction than
perpendicular to the polarization direction. Consequently, there should be
less reconstructable $K^{\ast 0}$ at $|\cos\theta^{*}|\sim 1$ than
$|\cos\theta^{*}|\sim 0$. It is noted that the picture is oversimplified, the
angular distribution should be diluted by the momentum of $K^{\ast 0}$ and the
rescattering probability depends not only on the hadron density profile but
also on the momentum distribution. The UrQMD model takes all these factors
into account. The resulted in $\rho_{00}$ in the particular $p_{T}$, and the
centrality bin is $0.3275\pm 0.0011$, which is lower than 1/3 by 0.0058
($>5\sigma$). Fig. 5 shows its $p_{T}$ and centrality dependence. The
deviation of reconstructable $K^{\ast 0}$ $\rho_{00}$ from 1/3 is more
significant in intermediate $p_{T}$ region and (semi-)central collisions,
reaching as large as -0.008 at $0.8<p_{T}<1.2~{}\textrm{GeV}/c$ and 0-20%
centrality.

Such a magnitude of deviation is significant compared with various theoretical
sources. In previous theoretical calculations, the sources include the
vorticity, electromagnetic fields, vector meson fields, fragmentation of
polarized quarks, axial or helicity charge, and turbulent color fields. Since
the different kinds of meson fields do not have a large correlation in space,
the contributions from vector meson fields are absent for $K^{\ast 0}$. The
spin alignment of $K^{\ast 0}$ is dominated by the vorticity field. The most
obvious source would be the electric part of the vorticity tensor, which gives
a negative contribution on the order of $10^{-4}$ PhysRevD.102.056013 ;
sheng2020can . The deviation of $\rho_{00}$ from the rescattering effect we
measured is about an order of magnitude higher than this most obvious source.
Then the magnetic part of the vorticity tensor gives a negative contribution,
while the electric field gives a positive contribution, and the fragmentation
of polarized quarks can give either a positive or negative contribution. These
three sources are at the level of $10^{-5}$ sheng2020can ; liang2005spin .
Helicity polarization gives a negative contribution at all centralities
PhysRevD.104.076016 . The exception from turbulent color fields is
$\rho_{00}(K^{*0})<\rho_{00}(\phi)<1/3$ PhysRevD.105.L011901 .

Figure 6: Eccentricity of $\pi^{\pm}$ as a function of the evolution time in
20-60$\%$ centrality

(a)

(b)

Figure 7: $\pi^{\pm}$ distribution in the transverse plane at (a) t = 5 fm/c
and (b) 20 fm/c.

Although the deviation of $\rho_{00}$ of $K^{\ast 0}$ due to hadronic
rescattering is significantly large compared to various theoretical
calculations, it is much smaller compared to the naive estimation based on the
initial eccentricity of the collisions or the anisotropy ($v_{2}$) of the
final hadrons. Since the $\pi\pi$ cross section is much higher than that of
$K\pi$, the rescattering effect of $K^{\ast 0}$ is dominated by the scattering
of the $\pi$ daughter with $\pi$s in the system. To better understand the
behaviors of $K^{\ast 0}$ $\rho_{00}$, we calculate the eccentricity of
charged $\pi$ in coordinate space as a function of the evolution time in
20-60$\%$ centrality, as shown in Fig. 6. The eccentricity is defined as

$\frac{\Sigma y^{2}-\Sigma x^{2}}{\Sigma y^{2}+\Sigma x^{2}},$ (2)

where $x$ is the position in the direction of the impact parameter, and $y$ is
the position in the direction of the reaction plane. Eccentricity tends to be
zero or even negative with time, which means the anisotropy of the medium
decreases with the fireball evolution. Figure 7 shows charged $\pi$
distributions in the transverse plane. The distribution is still an ellipse at
5 fm/$c$ but close to a circle at 20 fm/$c$. According to Fig. 2, at the time
when most $K^{\ast 0}$ decays, the eccentricity of the medium is already
small, and the loss of reconstructable $K^{\ast 0}$ caused by rescattering is
not so significantly different along different directions.

###  3.2 Spin Alignment with Respect to the Production Plane

(a)

(b)

Figure 8: $\rho_{00}$ with respect to the production plane of all and
reconstructable $K^{\ast 0}$ as a function of $p_{T}$ in 20-60$\%$ centrality
(a) and that of reconstructable $K^{\ast 0}$ as a function of centrality in
two different $p_{T}$ intervals (b) in Au+Au collisions at
$\sqrt{s_{\mathrm{NN}}}$ = 200 GeV.

Since the initially produced $K^{\ast 0}$ has no spin alignment in the UrQMD
model, the measurements of all $K^{\ast 0}$ with respect to the production
plane are similar to those with respect to the reaction plane. As shown in
Fig. 4(c), the $\cos\theta^{*}$ distribution is flat for all $K^{\ast 0}$ in
20-60% Au+Au collisions, and the extracted $\rho_{00}=0.3329\pm 0.0009$ is
consistent with 1/3 within uncertainties. The $p_{T}$ dependence is shown in
Fig. 8(a) as open circles.

We then studied the $\rho_{00}$ value after the rescattering effect with
respect to the production plane. Figure 4(d) shows that there is less
reconstructable $K^{\ast 0}$ toward $|\cos\theta^{*}|=0$ or
$\theta^{*}=\pi/2$, which with the results with respect to the reaction plane.
The extracted $\rho_{00}$ value in this particular $p_{T}$ and centrality bin
is $0.3506\pm 0.0011$, which is larger than 1/3 by 0.017. Figure 8 shows its
$p_{T}$ and centrality dependence. The deviation of reconstructable $K^{\ast
0}$ $\rho_{00}$ from 1/3 is more significant in intermediate $p_{T}$ region
and (semi-)central collisions, reaching as large as 0.03 at
$0.8<p_{T}<1.2~{}\textrm{GeV}/c$ and 0-20% centrality. This is because the
decay daughters are uniformly distributed in the rest frame of $K^{\ast 0}$ as
we proved before, and daughters that move back to the $K^{\ast 0}$ momentum
direction in its rest frame will have relatively lower momentum when boosted
to the laboratory frame. Then, these daughters are expected to be rescattered
more easily since rescattering affects more for low momentum particles. Since
the polarization direction is always perpendicular to the $K^{\ast 0}$
momentum direction, there are fewer reconstructable $K^{\ast 0}$s at
$|\cos\theta^{*}|\sim 0$ than $|\cos\theta^{*}|\sim 1$. To look it in more
detail, we calculate the distribution of $\theta^{**}$, which is the angle
between the momentum direction of the parent $K^{\ast 0}$ in the laboratory
frame and the momentum direction of a daughter particle in the $K^{\ast 0}$
rest frame. Figure 9(a) shows the $\cos\theta^{**}$ distribution for all K or
$\pi$, and the distribution is flat. Figure 9(b) shows the $\cos\theta^{**}$
distribution for scattered K and $\pi$, and there are more scattered daughters
toward $\cos\theta^{**}=-1$ or in the opposite direction of $K^{\ast 0}$. For
the results with respect to the reaction plan, as the $K^{\ast 0}$ momentum
direction is more or less isotropic with respect to the reaction plane, this
effect is small.

(a)

(b)

Figure 9: The $\cos\theta^{**}$ distribution for (a) all decay daughters and
(b) scattered decay daughters for 0.8 ¡ $p_{T}(K^{\ast 0}$ ¡ 1.2 GeV/c in
20-60$\%$ centrality.

##  4 Summary

In summary, we have studied the rescattering effect on the spin alignment
parameter $\rho_{00}$ of $K^{\ast 0}$ in Au+Au collisions at
$\sqrt{s_{\mathrm{NN}}}$ = 200 GeV based on the UrQMD model. The results show
that the rescattering effect gives a negative contribution to $\rho_{00}$ with
respect to the reaction plane. The $p_{T}$ dependence and centrality
dependence are also measured. The maximum deviation reaches as large as
-0.008, which is significant compared with various theoretical sources and
comparable to the current experimental precision.The results with respect to
the production plane show that the rescattering effect gives a positive
contribution to $\rho_{00}$, and the deviation is much larger than that with
respect to the reaction plane. This is mainly due to the momentum dependence
of the rescattering probability of $K^{\ast 0}$ decay daughters. These results
indicate that the rescattering effect should be considered in $K^{\ast 0}$
spin alignment. When preparing this manuscript, we noticed a study of the same
effect using the AMPT model Shen:2021pds . The deviation they found is
different from what we find from this study. The difference may come from the
different implementations of hadronic rescattering processes in these two
models.

## Acknowledgment

The authors thank Drs. Bedangadas Mohanty, Jinhui Chen, Aihong Tang and
Zhangbu Xu for helpful discussions. This work is supported in part by the
Strategic Priority Research Program of Chinese Academy of Sciences with Grant
No. XDB34030000 and the National Natural Science Foundation of China with
Grant No. 11720101001 and Natural Science Foundation of Anhui Province with
Grant No. 2208085J23.

## References

  * (1) Zuo-Tang Liang and Xin-Nian Wang.  Globally polarized quark-gluon plasma in non-central A+A collisions.  Phys. Rev. Lett., 94:102301, 2005.  [Erratum: Phys.Rev.Lett. 96, 039901 (2006)].
  * (2) F. Becattini, F. Piccinini, and J. Rizzo.  Angular momentum conservation in heavy ion collisions at very high energy.  Phys. Rev. C, 77:024906, 2008.
  * (3) Dmitri E. Kharzeev, Larry D. McLerran, and Harmen J. Warringa.  The Effects of topological charge change in heavy ion collisions: ’Event by event P and CP violation’.  Nucl. Phys. A, 803:227–253, 2008.
  * (4) Zuo-Tang Liang and Xin-Nian Wang.  Globally polarized quark-gluon plasma in non-central A+A collisions.  Phys. Rev. Lett., 94:102301, 2005.  [Erratum: Phys.Rev.Lett. 96, 039901 (2006)].
  * (5) L. Adamczyk et al.  Global $\Lambda$ hyperon polarization in nuclear collisions: evidence for the most vortical fluid.  Nature, 548:62–65, 2017.
  * (6) K. Schilling, P. Seyboth, and Guenter E. Wolf.  On the Analysis of Vector Meson Production by Polarized Photons.  Nucl. Phys. B, 15:397–412, 1970.  [Erratum: Nucl.Phys.B 18, 332 (1970)].
  * (7) Shreyasi Acharya et al.  Global polarization of $\Lambda\bar{\Lambda}$ hyperons in Pb-Pb collisions at $\sqrt{s_{NN}}$ = 2.76 and 5.02 TeV.  Phys. Rev. C, 101(4):044611, 2020.
  * (8) MS Abdallah, BE Aboona, J Adam, L Adamczyk, JR Adams, JK Adkins, G Agakishiev, I Aggarwal, MM Aggarwal, Z Ahammed, et al.  Observation of global spin alignment of $\phi$ and $k^{*0}$ vector mesons in nuclear collisions.  arXiv preprint arXiv:2204.02302, 2022.
  * (9) Berndt Müller and Di-Lun Yang.  Anomalous spin polarization from turbulent color fields.  Phys. Rev. D, 105:L011901, Jan 2022.
  * (10) STAR Collaboration.  Global polarization of $\mathrm{\Lambda}$ hyperons in au \+ au collisions at $\sqrt{{s}_{NN}}=200$ gev.  Phys. Rev. C, 98:014910, Jul 2018.
  * (11) Jian-Hua Gao.  Helicity polarization in relativistic heavy ion collisions.  Phys. Rev. D, 104:076016, Oct 2021.
  * (12) Xin-Li Sheng, Qun Wang, and Xin-Nian Wang.  Improved quark coalescence model for spin alignment and polarization of hadrons.  Phys. Rev. D, 102:056013, Sep 2020.
  * (13) Xiao-Liang Xia, Hui Li, Xu-Guang Huang, and Huan Zhong Huang.  Local spin alignment of vector mesons in relativistic heavy-ion collisions.  Physics Letters B, 817:136325, 2021.
  * (14) Yang-guang Yang, Ren-hong Fang, Qun Wang, and Xin-nian Wang.  Quark coalescence model for polarized vector mesons and baryons.  Phys. Rev. C, 97:034917, Mar 2018.
  * (15) Zuo-Tang Liang and Xin-Nian Wang.  Spin alignment of vector mesons in non-central A+A collisions.  Phys. Lett. B, 629:20–26, 2005.
  * (16) Xin-Li Sheng, Lucia Oliva, and Qun Wang.  What can we learn from the global spin alignment of $\phi$ mesons in heavy-ion collisions?  Physical Review D, 101(9):096005, 2020.
  * (17) J. Adams et al.  K(892)* resonance production in Au+Au and p+p collisions at s(NN)**(1/2) = 200-GeV at STAR.  Phys. Rev. C, 71:064902, 2005.
  * (18) M. M. Aggarwal et al.  $K^{*0}$ production in Cu+Cu and Au+Au collisions at $\sqrt{s_{N}N}=62.4$ GeV and 200 GeV.  Phys. Rev. C, 84:034909, 2011.
  * (19) Shreyasi Acharya et al.  Evidence of rescattering effect in Pb-Pb collisions at the LHC through production of $\rm{K}^{*}(892)^{0}$ and $\phi(1020)$ mesons.  Phys. Lett. B, 802:135225, 2020.
  * (20) Shaowei Lan, Zi-Wei Lin, Shusu Shi, and Xu Sun.  Effects of finite coverage on global polarization observables in heavy ion collisions.  Physics Letters B, 780:319–324, 2018.
  * (21) Diyu Shen, Jinhui Chen, and Zi-Wei Lin.  The effect of hadronic scatterings on the measurement of vector meson spin alignments in heavy-ion collisions.  Chin. Phys. C, 45(5):054002, 2021.


EUROPEAN ORGANIZATION FOR NUCLEAR RESEARCH (CERN)

 February 28, 2022

Quantum Machine Learning for $b$-jet identification

Alessio Gianelle1, Patrick Koppenburg2, Donatella Lucchesi1,3, Davide
Nicotra3,4, Eduardo Rodrigues5, Lorenzo Sestini1, Jacco de Vries4, Davide
Zuliani1,3,6

1INFN Sezione di Padova, Padova, Italy
2Nikhef National Institute for Subatomic Physics, Amsterdam, Netherlands
3Università degli Studi di Padova, Padova, Italy
4Universiteit Maastricht, Maastricht, Netherlands
5University of Liverpool, Liverpool, United Kingdom
6European Organization for Nuclear Research (CERN), Geneva, Switzerland

Machine Learning algorithms have played an important role in hadronic jet
classification problems. The large variety of models applied to Large Hadron
Collider data has demonstrated that there is still room for improvement. In
this context Quantum Machine Learning is a new and almost unexplored
methodology, where the intrinsic properties of quantum computation could be
used to exploit particles correlations for improving the jet classification
performance. In this paper, we present a brand new approach to identify if a
jet contains a hadron formed by a $b$ or $\bar{b}$ quark at the moment of
production, based on a Variational Quantum Classifier applied to simulated
data of the LHCb experiment. Quantum models are trained and evaluated using
LHCb simulation. The jet identification performance is compared with a Deep
Neural Network model to assess which method gives the better performance.

Submitted to JHEP

© 2024 CERN for the benefit of the LHCb collaboration. CC BY 4.0 licence.

##  1 Introduction

Machine Learning (ML) methods are widely used in experimental particle physics
[1] data manipulation. One of the most successful applications is the
classification of hadronic jets at the Large Hadron Collider (LHC) experiments
[2]. Jets are streams of particles produced via fragmentation and
hadronization of quarks and gluons that emerge from particle collisions,
_e.g._ proton-proton collisions at LHC. They are complex objects formed by
many detectable particles, and it is possible to identify their properties by
exploiting particle content and correlations, commonly referred to as jet
substructure. Typical jet classification problems are the identification of
the heavy-flavour hadron produced in the jet hadronization (_e.g._ $b$ hadron
vs $c$ hadron) or the identification of the charge of the heavy-flavour quark
that constitutes this hadron (_e.g._ $b$ vs $\bar{b}$). State-of-the-art ML
methods, such as Deep Neural Networks (DNN) [3], Graph Neural Networks [4],
Tensor Networks [5] and Particle Clouds [6], have been applied to jets data
collected by the LHC experiments, with a clear improvement of the
classification performance with respect to classical non-ML methods.

Recently, Quantum Computing (QC) has set the scene for a revolution in ML. The
new approach consists of using quantum circuits to tackle classification
tasks, in the framework of Quantum Machine Learning (QML) [7]: data are
embedded into a quantum state, which is then passed to a variational quantum
circuit, and by varying the circuit parameters a training procedure is
performed by means of minimising a classical loss function. Probability
measurements of the final state are then used to perform the classification.
Given the intrinsic properties of quantum computation, namely superposition
and entanglement, the new approach could lead to new insights from the
classification point of view. Jets that originate from gluons, or quarks of a
certain charge and flavor, would have a characteristic particle content and
correlations between them, which could be exploited to aid the identification
of the original particle. It is evident that this will open new physics
researches possibilities not even considered nowadays.

QML techniques have recently been applied to solve High Energy Physics (HEP)
problems, such as signal versus background separation [8, 9, 10, 11], and
particle track reconstruction [12, 13]. A more detailed review of QML
applications to HEP can be found in Ref. [14]. This paper presents the first
application of QML to the task of jet identification. QML methods are
performed on simulators and applied to simulated LHCb samples, to identify the
charge of the $b$ quark that forms the $b$ hadron produced in the jet
hadronization. In the rest of the paper, this task is simply referred as
$b$-jet tagging.

The paper is structured in the following way: Sec. 2 provides a description of
the LHCb jet reconstruction and identification together with the used dataset.
In Sec. 3 the considered QML algorithms are presented while the analysis flow
is described in Sec. 4. The results are discussed in Sec. 5. The conclusions
and future developments are presented in Sec. 6.

##  2 Jet reconstruction and identification at LHCb

LHCb [15] is a single-arm spectrometer designed to study $b$ and $c$ hadrons
in the forward region of proton-proton collisions. The reference system used
at LHCb is defined by the $z$-axis (the beam axis) parallel to the proton
beam, the $x$-axis parallel to the gravity acceleration, and the $y$-axis
perpendicular to the other two. The direction of particle momentum is
identified by the angles $\theta$ and $\phi$, where $\theta$ is the angle
between the momentum and the $z$-axis, and $\phi$ is the azimuthal angle
between the projection of the momentum in the $xy$ plane and the $y$-axis. The
pseudorapidity $\eta$ is defined as
$\eta=-\mathrm{log}\left[\mathrm{tan}\left(\frac{\theta}{2}\right)\right]$.
The LHCb detector covers the region in the pseudorapidity range $2<\eta<5$,
and consists of a tracking system and a particle identification system [16].
The tracking system is formed by a vertex detector, several tracking stations
and a dipole magnet. The vertex detector efficiently reconstructs the decay
vertex of $b$\- and $c$-hadron decays, while the tracking stations measure the
trajectories (tracks) of charged particles and their momenta. The particle
identification system is formed by two Ring Imaging Cherenkov detectors, two
calorimeters (electromagnetic and hadronic) and a muon detector, that allows
to precisely determine the type of the particles produced in the collision.

Jets are reconstructed using inputs selected by the Particle Flow algorithm
[17]. These are charged particles detected by the tracking system, and neutral
particles reconstructed in the calorimeter system as energy clusters isolated
from tracks [18]. The selected particles represent the input of the
anti-$k_{t}$ clustering algorithm [19] for jet clusterization. The jet
momentum is defined as the sum of the momenta of the particles forming the
jet. The jet axis is defined as the direction of the jet momentum. Particles
inside a jet are approximately contained in a cone structure with a distance
from the jet axis $\Delta R=\sqrt{(\Delta\eta)^{2}+(\Delta\phi)^{2}}=0.5$,
where $\Delta\eta$ is the difference in pseudo-rapidity and $\Delta\phi$ is
the difference in the azimuthal angle with respect to the jet axis.

The goal is to distinguish between jets that contain a $b$ or $\bar{b}$ hadron
just after the hadronization, _i.e._ in the instant of $b$-hadron production
and not at decay, since neutral $B$-mesons can undergo flavour oscillation.
Therefore the analysis is restricted to a sample of jets that belong to these
two categories, labelled as $b$ jets and $\bar{b}$ jets.

This preliminary selection is performed in two steps:

  * •

reconstruction of a vertex (secondary vertex) significantly displaced from the
primary proton-proton interaction point, using tracks detected by the vertex
detector [20], representing the $b$-hadron decay point;

  * •

identification of the jet that contains the secondary vertex within its cone.

The $b$-jet tagging subsequently becomes a binary classification problem where
the jet can belong to one of the two exclusive categories: $b$ jets or a
$\bar{b}$ jets. The charge of the $b$ quark at production is correlated to the
charge of the $b$-hadron decay products. This correlation is not perfect,
since neutral $B$-mesons can oscillate, and the charge of the $b$ quark at
production may be different from the charge at decay. As an example, in semi-
leptonic decays, the $b$ hadron can produce a muon, whose charge is directly
related to the charge of the $b$ quark. However, due to possible $B$-meson
oscillations and to the large number of particles, including muons, in a jet
the information is diluted and that has to be taken into account.

Two types of $b$-jet tagging algorithms are currently used in LHCb:

  * •

_exclusive algorithms_ , based on information coming from particles inside the
jet strictly correlated with the $b$-hadron decay, such as the muon;

  * •

and _inclusive algorithms_ , which aim to exploit the jet sub-structure, i.e.
information coming from the jet constituents, as shown in Fig. 1.

Figure 1: Sketch representing possible jet tagging methods. In the exclusive
method the information comes from a particle, e.g. the muon, whose charge is
correlated to the $b$ hadron (lower jet); in the inclusive method, the
information is extracted from the jet constituents (upper jet). The magnitude
of the particle momentum transverse to the jet axis is labelled as
$p^{\mbox{rel}}_{\mathrm{T}}$.

The QML approach presented in this paper belongs to the category of _inclusive
algorithms_. However, for the sake of comparison, the results are compared to
the _muon tagging_ [21], which is an exclusive algorithm. This tagger selects
the muon with the highest momentum with respect to the $z$ axis,
$p_{\mathrm{T}}$, inside the jet. This simple requirement is sufficient to
identify the muon coming from a $b$ hadron and therefore to infer the quark
charge by measuring the muon charge. The efficiency of this algorithm is
limited by the probability that a $b$ hadron decays to a final state with a
muon, which is $\sim 10\%$ [22].

The performance of different $b$-jet tagging algorithms are compared using the
tagging power $\epsilon_{\text{tag}}$ [23, 24, 25], defined as

$\epsilon_{\text{tag}}=\epsilon_{\text{eff}}(2a-1)^{2}$ (1)

as the figure of merit, where $\epsilon_{\text{eff}}$ is the tagging
efficiency, _i.e._ the fraction of jets where the classifier takes a decision,
and $a$ is the accuracy, _i.e._ the fraction of correctly tagged jets with
respect to the tagged jets. The tagging power is the effective fraction of
events that contribute to the statistical uncertainty in a measurement where
the $b$-jet tagging is applied.

###  2.1 Data samples description

LHCb simulated samples are used in the studies presented in this paper. The
$b\bar{b}$ di-jets samples are produced within the LHCb simulation framework
[26], which uses Pythia8 [27] to generate proton-proton interactions and jet
fragmentation and hadronization at a center-of-mass energy of $13\,$TeV, and
EvtGen [28] to simulate $b$ hadron decays. The Geant4 software [29], embedded
in the LHCb framework, is used to simulate the detector response. Pairs of $b$
and $\bar{b}$ jets are selected by requiring for each jet a $p_{\mathrm{T}}$
greater than 20 GeV/c and a $\eta$ in the range $2.2<\eta<4.2$, to ensure that
they are well inside the instrumented part of the detector. After the pre-
selection, a fixed number of $16$ different features related to the jet
substructure are used as input to the classifiers. Among the reconstructed
particles inside a jet the muon, kaon, pion, electron and proton with the
highest $p_{\mathrm{T}}$ are selected. For each particle three physical
variables are considered: the magnitude of the transverse momentum to the jet
axis ($p^{\mbox{rel}}_{\mathrm{T}}$), the charge ($q$), and the distance,
measured in the ($\eta$,$\phi$) space, between the particle and the jet axis
($\Delta R$). If a particle type is missing, the relative features are set to
$0$. The last feature is the weighted jet charge $Q$, defined as the sum of
the charges of the particles inside the jet weighted with the particles
$p^{\text{rel}}_{\mathrm{T}}$

$Q=\frac{\sum_{i}(p^{\text{rel}}_{\mathrm{T}})_{i}q_{i}}{\sum_{i}(p^{\text{rel}}_{\mathrm{T}})_{i}}~{}.$
(2)

The analysis for the $b$-jet identification is performed with two datasets.
The complete dataset includes the events selected with the 16 features
described above. The muon dataset contains jets with at least one muon and
only four features: $p_{\mathrm{T}}^{\text{rel}}$, $\Delta R$, $q$ of the muon
and the weighted jet charge $Q$. Table 1 summarises the characteristics of the
data samples.

| Muon | Kaon | Pion | Electron | Proton |
---|---|---|---|---|---|---
Dataset | $p_{\mathrm{T}}^{\mathrm{rel}}$ | $q$ | $\Delta R$ | $p_{\mathrm{T}}^{\mathrm{rel}}$ | $q$ | $\Delta R$ | $p_{\mathrm{T}}^{\mathrm{rel}}$ | $q$ | $\Delta R$ | $p_{\mathrm{T}}^{\mathrm{rel}}$ | $q$ | $\Delta R$ | $p_{\mathrm{T}}^{\mathrm{rel}}$ | $q$ | $\Delta R$ |  $Q$
Complete |  ✓ |  ✓ | ✓ |  ✓ |  ✓ | ✓ |  ✓ |  ✓ | ✓ |  ✓ |  ✓ | ✓ |  ✓ |  ✓ | ✓ |  ✓
Muon |  ✓ |  ✓ | ✓ |  |  |  |  |  |  |  |  |  |  |  |  |  ✓

Table 1: Summary of the features contained in the two datasets.

##  3 Quantum Machine Learning models

A quantum algorithm is implemented by means of a quantum circuit, namely a
collection of linked quantum gates acting on a $n$-qubit quantum state: the
measurements on the final state represent the outcome of the quantum
algorithm. Parametrized Quantum Circuits (PQCs) [30] are a type of circuit
that contains adjustable gates with tunable parameters. The Variational
Quantum Classifier (VQC) [31] is a hybrid quantum-classical algorithm to
perform classification tasks using a Machine Learning model based on a PQC
with the following structure:

Data encoding

data $x$, the features representing the jet substructure in this application,
are pre-processed and encoded into a subset of the parameters of a PQC. The
stage produces a quantum state $|x\rangle$ representing the input jet.

Variational circuit

the state $|x\rangle$ is processed by a PQC, $U(\theta)$, featuring trainable
parameters $\theta$ to be optimised during the training phase. This stage
produces a final state $|\psi\rangle=U(\theta)|x\rangle$.

Prediction

expectation values computed on the final state $|\psi\rangle$ are mapped to
probabilities for the two labels, $P_{b}$ and $P_{\bar{b}}$. The training
process aims to match the label predictions with the true charge of the
$b$-hadron in the jet in the instant of production, available in the
simulations.

Two different PQC models are studied in this work: Amplitude Embedding and
Angle Embedding, described below.

###  3.1 Amplitude Embedding

The Amplitude Embedding model consists in a PQC made by an embedding circuit
followed by a variational circuit. The schematic representation of this model
is shown in Fig. 2. The embedding circuit consists of an Amplitude Encoder
that encodes up to $2^{n}$ features into the amplitude of a $n$-qubit quantum
state, or equivalently, a vector of $N$ features can be encoded using
$\lceil\log_{2}N\rceil$ qubits:

$|x\rangle=\sum_{i=1}^{2^{n}}x_{i}\,|n_{i}\rangle$ (3)

where $x_{i}$ is the $i^{\text{th}}$ feature and $|n_{i}\rangle$ is the
$i^{\text{th}}$ vector of the computational basis. The definition requires the
$x$ vector to be normalised,

$\sum_{i}|x_{i}|^{2}=1.$ (4)

If the number of features to encode is not a power of 2, the remaining
amplitudes can be padded with constant values. This model embeds the 16 (4)
variables of the _complete dataset_ (_muon dataset_) into the amplitudes of a
4-qubit (2-qubit) quantum state. The variational stage consists of a variable
number $L$ of strongly entangling layers. A strongly entangling layer consists
of trainable generic rotational gates $R(\alpha_{i},\beta_{i},\gamma_{i})$
applied to each qubit followed by a collection of CNOT gates applied to
neighbouring pairs of qubits, considering the last one as a neighbour of the
first one. The complexity of this kind of circuit can be tuned by changing the
number of strongly entangling layers $L$: for a generic $n$-qubit circuit, the
number of trainable parameters of the model $N_{\text{par}}$ is equal to

$N_{\text{par}}=3\times n\times L.$ (5)

On the final state, the expectation of the Pauli operator of the first qubit
$\langle\sigma^{0}_{z}\rangle\in[-1,+1]$ is measured and used to define the
probabilities $P_{b}$ and $P_{\bar{b}}$ of being a $b$-jet and a
$\bar{b}$-jet, respectively:

$\displaystyle P_{b}$
$\displaystyle=\frac{1}{2}(\langle\sigma^{0}_{z}\rangle+1)$ (6) $\displaystyle
P_{\bar{b}}$
$\displaystyle=\frac{1}{2}(1-\langle\sigma^{0}_{z}\rangle)=1-P_{b}.$ (7)

Figure 2: Circuit representation of the Amplitude Embedding model. In blue,
variables are embedded into the amplitudes of a quantum state. In red,
trainable generic rotational gates to be optimised during the training phase.
In green, CNOT gates entangling qubits with a circular topology.

###  3.2 Angle Embedding

The structure of the Angle Embedding model, represented in Fig. 3, differs
from the Amplitude Embedding model in the encoding used to embed the features
of the datasets into the quantum state: in this case, the embedding circuit
consists in a Angle Encoder that embeds 16 (4) features of the _complete
dataset_ (_muon dataset_) as rotation angles of 16 (4) $x$-axis rotational
gates $R_{x}(\theta_{i})$. Therefore, this circuit structure requires a one-
to-one correspondence between qubits and input features: that makes it
impractical to adopt with high-dimensionality datasets, due to computational
constraints of quantum simulators. The variational stage of the circuit is
identical to the Amplitude Embedding model, featuring a variable number $L$ of
strongly entangling layers that can be opportunely chosen to tune the number
of parameters $N_{\text{par}}$, defined in Eq. 5, and, therefore, the
complexity. The measurement of the expectation value of the Pauli $\sigma_{z}$
operator is mapped to the tagging probabilities $P_{b}$ and $P_{\bar{b}}$ as
expressed in Eq. 6 and Eq. 7, identically to the Amplitude Embedding model.

Figure 3: Circuit representation of the Angle Embedding model. In blue,
$x$-axis rotational gates used to embed the variables into the quantum
circuit. In red, trainable generic rotational gates to be optimised during the
training phase. In green, CNOT gates entangling qubits with a circular
topology.

##  4 $\mathbf{b}$-jets identification procedure

Quantum circuits are simulated by means of noiseless simulators (noise impact
is studied in Sec. 5.3) using Pennylane [32], a Python framework designed
specifically for QML applications. The quantum circuit is embedded into a
classical optimisation algorithm, using the Jax [33] Python library. Since the
quantum algorithms results are compared to classical DNN ones, the same
analysis is performed with a standard feed-forward network, implemented using
the Keras [34] framework with the TensorFlow [35] back-end. Additional details
on the structure and the optimisation of the DNN are reported in App. A.

###  4.1 Training and testing phases

The _muon_ and _complete_ datasets are both split into training and testing
sub-datasets: about $60\%$ of the samples are used in the training process
that includes also the validation and the remaining $40\%$ are used to test,
evaluate and compare the classifiers. In the _muon dataset_ analysis, 60000
jets are used for training and 40000 jets are used for testing. The _complete
dataset_ training is performed on 400000 jets and remaining 290000 are used
for testing and assessing performance. In the analysis of the _muon dataset_
the Angle Embedding structure is studied and compared to a DNN with the same
input variables. In the case of the _complete dataset_ , both Angle Embedding
and Amplitude Embedding classifiers are considered, and the 16 variables are
used as input for a DNN. The training process aims to find the values of the
model parameters $\theta$ that minimise the Mean Squared Error loss function

$\mathcal{L}(\theta)=\frac{1}{N}\sum_{i=1}^{N}(P^{i}_{b}(\theta)-T^{i})^{2},$
(8)

where $N$ is the number of training jets, $P^{i}_{b}$ is the predicted
probability, defined in Eq. 6, for the $i$-th jet, and $T^{i}$ is the target
probability for the $i$-th jet, i.e., $1$ for a $b$ jet and $0$ for a
$\bar{b}$ jet. Due to the large number of jets in the datasets, the quantum
models are trained implementing a mini-batch gradient descent [36] algorithm
using the ADAM optimiser [37] to minimise Eq. 8. The training dataset is split
in several mini-batches containing a fixed number of training samples. During
each training step, the gradient of Eq. 8 is evaluated, averaging over the
training samples of a mini-batch, and used to update model parameters. A
training epoch is completed when the whole training dataset is processed,
namely after a number of steps equal to the number of mini-batches. Unless
specified otherwise, the models are trained with learning rate111In Machine
Learning, the learning rate $\xi$ of an optimisation algorithm is usually
defined as the scaling factor applied to the gradient of the loss function
when updating the parameters. It can be tuned as a parameter of the learning
process. $\xi=0.01$ for 100 epochs, while the mini-batch size is fixed to the
maximum value allowed by memory constraints. The output of the classifier
gives the probability that a jet is generated by a $b$ or a $\bar{b}$ quark.
The label with the highest probability is assigned to the jet, i.e. if
$P_{b}>0.5$ ($P_{b}<0.5$) then it is classified as a $b$ jet ($\bar{b}$ jet).
In Fig. 4(a) the output distributions for the two classes ($b$ and $\bar{b}$
jets) after the training procedure are shown; a separation between the two
distributions around $0.5$ is visible leading to a good classification. It
should be noted that the $P_{b}$ distribution is shifted toward 1 for $b$
quarks, and toward 0 for $\bar{b}$ quarks, as expected. Fig. 4(b) shows the
output distribution for the Angle Embedding classifier on 16 qubits: in green
the jets that are correctly classified, in red the jets that are wrongly
classified and the sum of all jets in grey. As expected, correctly classified
jets tend to stay close to $0$ and $1$ while the wrongly classified jets are
peaked around $0.5$, where the prediction power is minimum.

(a)

(b)

Figure 4: Probability distributions for jet tagged to $b$ (blue) and $\bar{b}$
quarks (yellow), showing separation around 0.5 (a). Probability distribution
for the Angle Embedding circuit: jet correctly (wrongly) tagged are plotted in
green (red), showing around 0.5 worse classification. The probability
distribution for all jets is shown in grey (b).

Figure 5 shows the Receiver Operating Characteristic curve (ROC) and the Area
Under Curve (AUC) for the DNN and the quantum classifiers for the muon dataset
and the complete dataset.

(a)

(b)

Figure 5: ROC distributions and AUC score for DNN (green), Angle Embedding
(blue) and Amplitude Embedding circuits (yellow) for the muon dataset (a) and
the complete dataset (b). The dashed line represents a random classifier.

##  5 Results on $b$-jet tagging

The performance of the classifiers is evaluated by using the jet tagging
power, defined in Eq. 1. The tagging power is computed as a function of the
jet $p_{\mathrm{T}}$ and $\eta$ for both the quantum and the classical
classifiers. In order to optimise the tagging power, a region symmetric with
respect to 0.5 is defined, where no classification is performed. The width
$\Delta_{\text{cut}}$ of the excluded region is defined for each classifier by
maximising the tagging power evaluated using all the jets in the dataset. Such
an exclusion region reduces the tagging efficiency because less jets are
tagged, but enhances the identification probability. The probability
distributions and the excluded region are shown in Fig. 6: indeed the region
where the prediction power is minimum is excluded. The width
$\Delta_{\text{cut}}$ of the excluded region for each classifier and for
_muon_ and _complete dataset_ are summarised in Tab. 2. For comparison, the
unoptimised tagging power, obtained by identifying as $b$ ($\bar{b}$) the jets
with $P_{b}$ above (below) 0.5, is presented in App. B.

Figure 6: Probability distributions for jet tagged to (blue) $b$ and (yellow) $\bar{b}$ quarks  Table 2: Width $\Delta_{\text{cut}}$ for different classifiers and dataset. | Classifier
---|---
Dataset | DNN | Angle Embedding | Amplitude Embedding
Muon | 0.30 | 0.25 | 0.16
Complete | 0.21 | 0.19 | 0.12

The tagging power $\epsilon_{\text{tag}}$ as a function of jet
$p_{\mathrm{T}}$ and $\eta$ for the classical and quantum classifiers applied
to the _muon dataset_ is shown in Fig. 7. All the distributions have similar
behaviour demonstrating that no bias is created by any algorithm. The tagging
power dependence on the jet $p_{\mathrm{T}}$ is as expected, since at high
$p_{\mathrm{T}}$ the reconstruction and identification of the jet particle
content is more difficult, leading to a lower tagging power [38]. The Angle
Embedding algorithm results are slightly below the DNN ones, but they are
compatible within the statistical uncertainty, while both algorithms
outperform the Amplitude Embedding model and the muon tagging approach. The
muon tagging performance is expected, since it only uses the muon charge $q$
for the prediction. The other algorithms use also the muon $p_{\mathrm{T}}$,
the muon $\Delta R$ and the weighted jet charge $Q$. The simple 2-qubit
Amplitude Embedding model shows a slightly better performance with respect to
the muon tagging, but still worse than the DNN and the Angle Embedding model.

(a)

(b)

Figure 7: Tagging power $\epsilon_{\text{tag}}$ with respect to (a) jet
$p_{\mathrm{T}}$ and (b) jet $\eta$ for the _muon dataset_. The Angle
Embedding circuit and the DNN show similar performance.

The tagging power $\epsilon_{\text{tag}}$ for DNN and quantum classifiers on
the _complete dataset_ as function of jet $p_{\mathrm{T}}$ and $\eta$, is
shown Fig. 8. Also in this case, the usual dependence on the jet
$p_{\mathrm{T}}$ is visible. As expected, the tagging power of the QML and DNN
is higher in the _complete dataset_ with respect to the _muon dataset_ , since
a larger number of features is used. As before, for QML and DNN the
performance is above the muon tagging approach, given that these classifiers
use the information coming from the jet substructure. Differently from the
application to the _muon dataset_ , in the _complete dataset_ the QML
algorithms perform slightly worse than the DNN, with slightly better
performance for the Angle Embedding structure than the Amplitude Embedding. It
can be deduced that the DNN makes a better use of the features when a larger
number of them is used.

(a)

(b)

Figure 8: Tagging power $\epsilon_{\text{tag}}$ with respect to (a) jet
$p_{\mathrm{T}}$ and (b) jet $\eta$ for the _complete dataset_. The quantum
algorithms perform slightly worse than the DNN, with the Angle Embedding
circuit performing better than the Amplitude Embedding circuit.

###  5.1 Dependence of the results on number of training events and circuit
depth

The dependence of the quantum algorithms performance on the number of training
events and the circuit complexity has to be evaluated if near-term
applications on quantum hardware are considered. These parameters have an
impact on the execution times and therefore on the possibility to use it. The
performance dependence on the number of training samples is an interesting
parameter to compare QML and DNN methods, in order to assess the differences
between the two approaches. Given the high computational efforts of simulating
complex circuits with several qubits, only the _muon dataset_ is used. For
QML, the Angle Embedding structure is considered with different number of
strongly entangled layers and different number of training events. The results
are compared with the same DNN considered in the previous section; a
comparison with more complex networks is described in App. A. The metric used
to quantify the goodness of the quantum classifier is the accuracy on a test
subset of 40000 jets. The performance is calculated averaging over 10 training
rounds. In Fig. 9 (a) the accuracy of the Angle Embedding circuit is shown as
a function of the number of layers of the circuit. As expected, by increasing
the depth of the circuit, and therefore its complexity, the accuracy
increases. This behaviour stops at around 5 layers, where the accuracy is
saturating and no further improvement is evident. It is clear that, for a
given number of features and training data, the Angle Embedding model does not
profit of an arbitrarily large number of layers, therefore it is possible to
keep a low number of layers, and subsequently a lower complexity of the
circuit, to obtain the best performance. This would reduce also the computing
time and resources needed for the simulation.

(a)

(b)

Figure 9: (a) Accuracy of the Angle Embedding structure on the _muon dataset_
versus the number of layers. (b) Accuracy of the (red) Angle Embedding
structure and (blue) DNN on the _muon dataset_ versus the number of training
events.

The accuracy as a function of the number of training events for the Angle
Embedding circuit and the DNN is shown in Fig. 9 (b). Increasing the number of
training events the performance of the quantum algorithm is similar to the
DNN, but when the number of training events decreases the quantum algorithm
keeps very high performance, while the DNN is not able to perform a good
classification. This means that, with respect to the DNN, the QML method
reaches optimal performance with a lower number of events, which helps in
reducing the resources usage.

###  5.2 Time performance

Time performance is a fundamental figure of merit to understand the
feasibility of simulating such quantum algorithms. Here, the time performance
for the Angle Embedding and the Amplitude Embedding circuits is evaluated for
the _muon dataset_ , as a function of the number of strongly entangling layers
in the circuit structure. The quantum algorithms are trained using 4 NVIDIA
Volta V100 GPUs, and the training time to train 60000 jets for 100 epochs is
evaluated. Results are shown in Fig. 10.

Figure 10: Training time for 100 epochs for the Angle Embedding (red) and
Amplitude Embedding (blue) circuits on the _muon dataset_ with respect to the
number of circuit layers.

Results show that it takes less time to train the Amplitude Embedding circuit
than the Angle Embedding structure, and the training time increases as the
number of layers, although with a greater rate for the Angle Embedding
circuit: this may be expected since the complexity of the circuit increases
with the number of layers, therefore it takes more time to simulate the
quantum circuit.

###  5.3 Noise models

In order assess the performance of the algorithms in quantum hardware it is
important to understand the impact of noise on quantum circuits. Two kinds of
noise can affect quantum algorithms:

  * •

coherent noise: it originates from unitary errors in the application of
quantum gates. This lead to the construction of a different quantum state with
respect to the desired one. A typical source of this kind of noise is non-
ideal calibrations of the quantum hardware;

  * •

incoherent noise: this noise results from the interaction between the quantum
hardware and the environment. This noise gives quantum states that are not
pure anymore and are described by mixed states, _i.e._ probability
distributions over different states.

The simulations of noise contribution taking into account both sources of
noise in quantum circuit measurements have been performed using the
_pennylane-qiskit_ plugin [32, 39]. This plugin allows to simulate noise
models coming from different real IBM quantum computers [40] and to keep the
Pennylane syntax. The result is a simulation of a quantum algorithm on a real
device structure. Four IBM quantum computers are considered: ibmq-belem, ibmq-
santiago, ibmq-jakarta and ibmq-toronto, which have different numbers of
qubits (respectively 5, 5, 7 and 27 qubits), different quantum volumes222The
quantum volume is the maximum size of a quantum circuit that can be
effectively implemented on a noisy intermediate-scale quantum device. In this
paper, the definition from Ref. [41] is adopted. (respectively 16, 32, 16 and
32) and different qubits structure, as shown in Fig. 11 for the ibmq-santiago
and ibmq-belem which have the same number of qubits.

(a)

(b)

Figure 11: Qubit structure of ibmq-belem (a) and ibmq-santiago (b).

Studies are performed on the Angle Embedding circuit structure with three
strongly entangled layers. A small subset of the _muon dataset_ is used
because simulating circuits including noise contribution is more time and
computationally consuming; on the other hand, with a low number of events the
quantum algorithm performance is sufficiently high, as shown in Fig.9. In this
way, a subset of 1000 jets of the _muon dataset_ is selected for training
while validation is performed on a subset of 10000 jets. For each noise model
the training is performed for 50 epochs using ADAM with a learning rate
$\xi=0.01$ and batch size of 10 jets. The results are averaged over five
rounds of training, using five independent training subsets. The relevant
figure of merit to assess noise models performance is the accuracy on the
validation test. The results are shown in Fig. 12 and summarised in Tab. 3.
Models including noise need more epochs to reach convergence, but in the end
the results are consistent with those of noiseless simulations within error.
Such a result demonstrates that the proposed circuit model for the _muon
dataset_ is robust to noise.

Figure 12: Validation accuracy for noise models as a function of the number of epochs. Blue (light-blue) band represent $1\sigma$ ($2\sigma$) uncertainty bounds for the noiseless model. Table 3: Accuracy for noisy circuits obtained with the _muon dataset_. Noise model | accuracy
---|---
no noise | $0.640\pm 0.017$
ibmq-belem | $0.629\pm 0.047$
ibmq-santiago | $0.633\pm 0.038$
ibmq-jakarta | $0.637\pm 0.042$
ibmq-toronto | $0.631\pm 0.044$

##  6 Conclusions

The first application of QML algorithms to identify the charge of the $b$
hadron produced in the jet hadronization, at the LHCb experiment has been
presented. The results using the _muon dataset_ show that the Angle Embedding
structure reaches a performance consistent with the DNN while being better
than the muon tagging approach. The amplitude encoding structure presents no
evident improvement with respect to the muon tagging method. When the
_complete dataset_ is used, the best quantum algorithm performance is given by
Angle Embedding structure. The study of the performance dependence on the
circuit depth has shown that the number of layers is a parameter to be
optimised. More layers, which means more complex structures, do not
necessarily result in improved performance, since saturation is reached. The
impact of the noise on the results appears to be negligible, suggesting that
the proposed circuits could be implemented on a quantum hardware if available.
QML algorithms achieve performance consistent with classical methods like the
DNN with low-complexity circuits and a smaller number of training events. This
could have important implications for LHC experiments where often the training
phase is the most expensive in terms of resources. However, when a large
number of features is employed, the DNN performs better than QML algorithms.
Here huge improvements are expected when the hardware will be available. In
fact, the comparison of QML models to classical kernel methods [42] shows that
QML models achieve classification tasks separating data in Hilbert spaces
whose dimension scales exponentially with the number of qubits involved. As
this number increases, quantum simulations on classical computers become
unfeasible.

The full exploitation of QML in high energy physics experiments at colliders,
just began. As for any brand new tool, unexpected applications may manifest
themselves while using. For example, results obtained involving tensor network
methods [5] have shown that quantum algorithms can allow to study correlations
among the features. This is done by measuring the entanglement correlations
between qubits of an optimised multi-qubit system. This approach could allow
to extract information on jet constituents correlations starting from the
measurement of the entanglement between the qubits of a trained VQCs.

## Acknowledgements

The authors would like to thank the LHCb Data Processing & Analysis (DPA)
project colleagues for supporting this publication and reviewing the work. The
authors acknowledge support from INFN (Italy) and from the Department of
Physics and Astronomy of the University of Padova. The authors would also like
to thank CERN Openlab for useful and fruitful discussions. The authors thank
the CINECA Supercomputing Center for awarding access to Marconi100 cluster
based at CINECA (Bologna, Italy).

##  Appendix A Deep Neural Network model and further comparisons

The DNN is built following a classical feed-forward structure shown in Fig.
13: it starts with a Batch Normalisation Layer[43] and it applies several
Dense layers, each one followed by a Dropout[44] layer. Depending on the
number of input features the input vector for the DNN changes (4 variables for
the muon dataset and 16 for the complete one). For the hidden layers the
ReLu333ReLu is the rectifier activation function, defined as
$\mathrm{ReLu}(x)=\max(0,x)$ activation is used while a sigmoid444The sigmoid
activation function is defined as $\mathrm{sig}(x)=[1+\exp{(-x)}]^{-1}$
function is applied to the output node. The network is trained using the ADAM
optimiser[37]. The model is trained for 250 epochs with an early stopping of
25 epochs on the test loss function and a learning rate of 0.0001. Given the
classification task, we take cross entropy as loss function. A scheme of the
DNN model is shown in Fig. 13.

Figure 13: DNN structure.

Different DNN structures have also been considered, in order to have a fair
comparison between quantum and classical classifiers. Usual state-of-the-art
networks make use of LSTM [45] and convolutional layers, therefore the
following DNN models are studied:

  * •

“LSTM” model: starting from the DNN structure we apply a LSTM layer to the
particle features (therefore excluding the ”global” variable of the total jet
charge) before the first Dense layer. The Dense structure is identical to the
DNN structure.

  * •

“LSTM+Conv” model: from the LSTM model we firstly apply a convolutional layer
of the particle features (as before not applied to the total jet charge).

The two models are shown in Fig. 14.

(a)

(b)

Figure 14: DNN structure for (a) the “LSTM” model and for (b) the “LSTM+CONV”
model .

The DNN performance in terms of tagging power is compared to the LSTM and
LSTM+CONV models. Results for tagging power as function of the jet
$p_{\mathrm{T}}$ and $\eta$ are shown in Fig. 15, where results for different
models are comparable within the error and therefore allow us to consider just
the DNN model for the comparison with quantum algorithms.

(a)

(b)

Figure 15: Tagging power $\varepsilon_{\text{tag}}$ with respect to (a) jet
$p_{\mathrm{T}}$ and (b) jet $\eta$ for the complete dataset applied to
different DNN models.

##  Appendix B Unoptimised tagging power distributions

In this section the distributions for the unoptimised tagging power are shown,
both for the muon and the complete datasets. This is done in order to check if
there are performance biases when cutting on the probability distributions to
maximise the tagging power. In Figs. 16 and 17 the unoptimised tagging power
distributions as function of jet $p_{\mathrm{T}}$ and $\eta$ respectively for
the muon and the complete datasets are shown. No evident biases are present
and the same considerations done for the optimised distributions are valid.

(a)

(b)

Figure 16: Unoptimised tagging power $\varepsilon_{\text{tag}}$ with respect
to (a) jet $p_{\mathrm{T}}$ and (b) jet $\eta$ for the muon dataset. The angle
embedding circuit and the DNN show similar performance.

(a)

(b)

Figure 17: Unoptimised tagging power $\varepsilon_{\text{tag}}$ with respect
to (a) jet $p_{\mathrm{T}}$ and (b) jet $\eta$ for the complete dataset. The
angle embedding circuit and the DNN show similar performance.

## References

  * [1] D. Guest, K. Cranmer, and D. Whiteson, _Deep Learning and its Application to LHC Physics_ , Ann. Rev. Nucl. Part. Sci. 68 (2018) 161, arXiv:1806.11484
  * [2] A. J. Larkoski, I. Moult, and B. Nachman, _Jet Substructure at the Large Hadron Collider: A Review of Recent Advances in Theory and Machine Learning_ , Phys. Rept. 841 (2020) 1, arXiv:1709.04464
  * [3] P. Baldi et al., _Jet Substructure Classification in High-Energy Physics with Deep Neural Networks_ , Phys. Rev. D93 (2016) 094034, arXiv:1603.09349
  * [4] J. Shlomi, P. Battaglia, and J.-R. Vlimant, _Graph neural networks in particle physics_ , Machine Learning: Science and Technology 2 (2021) 021001
  * [5] T. Felser et al., _Quantum-inspired machine learning on high-energy physics data_ , npj Quantum Information 7 (2021) 111
  * [6] H. Qu and L. Gouskos, _ParticleNet: Jet Tagging via Particle Clouds_ , Phys. Rev. D101 (2020) 056019, arXiv:1902.08570
  * [7] J. Biamonte et al., _Quantum machine learning_ , Nature 549 (2017) 195–202
  * [8] K. Terashi et al., _Event Classification with Quantum Machine Learning in High-Energy Physics_ , Computing and Software for Big Science 5 (2021)
  * [9] S. L. Wu et al., _Application of quantum machine learning using the quantum kernel algorithm on high energy physics analysis at the LHC_ , Phys. Rev. Res. 3 (2021) 033221, arXiv:2104.05059
  * [10] A. Blance and M. Spannowsky, _Quantum Machine Learning for Particle Physics using a Variational Quantum Classifier_ , doi: 10.1007/JHEP02(2021)212 arXiv:2010.07335
  * [11] J. Y. Araz and M. Spannowsky, _Classical versus Quantum: comparing Tensor Network-based Quantum Circuits on LHC data_ , arXiv:2202.10471
  * [12] C. Tüysüz et al., _Particle track reconstruction with quantum algorithms_ , EPJ Web Conf. 245 (2020) 09013
  * [13] L. Funcke et al., _Studying quantum algorithms for particle track reconstruction in the LUXE experiment_ , in 20th International Workshop on Advanced Computing and Analysis Techniques in Physics Research: AI Decoded - Towards Sustainable, Diverse, Performant and Effective Scientific Computing, 2022, arXiv:2202.06874
  * [14] W. Guan et al., _Quantum machine learning in high energy physics_ , Machine Learning: Science and Technology 2 (2021) 011003
  * [15] LHCb collaboration, A. A. Alves, Jr. et al., _The LHCb Detector at the LHC_ , JINST 3 (2008) S08005
  * [16] LHCb collaboration, R. Aaij et al., _LHCb Detector Performance_ , Int. J. Mod. Phys. A 30 (2015) 1530022, arXiv:1412.6352
  * [17] ALEPH collaboration, D. Buskulic et al., _Performance of the ALEPH detector at LEP_ , Nucl. Instrum. Meth. A360 (1995) 481
  * [18] LHCb collaboration, R. Aaij et al., _Study of forward Z + jet production in pp collisions at $\sqrt{s}=7$ TeV_, JHEP 01 (2014) 033, arXiv:1310.8197
  * [19] M. Cacciari, G. P. Salam, and G. Soyez, _The anti- $k_{t}$ jet clustering algorithm_, JHEP 04 (2008) 063, arXiv:0802.1189
  * [20] LHCb collaboration, R. Aaij et al., _Identification of beauty and charm quark jets at LHCb_ , JINST 10 (2015) P06013, arXiv:1504.07670
  * [21] LHCb collaboration, R. Aaij et al., _First measurement of the charge asymmetry in beauty-quark pair production_ , Phys. Rev. Lett. 113 (2014) 082003, arXiv:1406.4789
  * [22] Particle Data Group, M. Tanabashi et al., _Review of Particle Physics_ , Phys. Rev. D98 (2018) 030001
  * [23] D0 Collaboration, _Measurements of $B_{d}$ mixing using opposite-side flavor tagging_, Phys. Rev. D74 (2006) 112002, arXiv:0609034v1
  * [24] CDF Collaboration, _Measurements of $B_{0}$ oscillations and calibration of flavor tagging in semileptonic decays_, CDF Note 8235 (2006)
  * [25] LHCb collaboration, R. Aaij et al., _Opposite-side flavour tagging of b mesons at the lhcb experiment_ , The European Physical Journal C 72 (2012)
  * [26] LHCb collaboration, M. Clemencic et al., _The LHCb simulation application, Gauss: Design, evolution and experience_ , J. Phys. Conf. Ser. 331 (2011) 032023
  * [27] T. Sjöstrand, S. Mrenna, and P. Skands, _A brief introduction to PYTHIA 8.1_ , Computer Physics Communications 178 (2008) 852
  * [28] D. J. Lange, _The EvtGen particle decay simulation package_ , Nucl. Instrum. Meth. A462 (2001) 152
  * [29] S. Agostinelli et al., _Geant4—a simulation toolkit_ , Nuclear Instruments and Methods in Physics Research Section A: Accelerators, Spectrometers, Detectors and Associated Equipment 506 (2003) 250
  * [30] M. Benedetti, E. Lloyd, S. Sack, and M. Fiorentini, _Parameterized quantum circuits as machine learning models_ , Quantum Science and Technology 4 (2019) 043001
  * [31] V. Havlíček et al., _Supervised learning with quantum-enhanced feature spaces_ , Nature 567 (2019) 209–212
  * [32] V. Bergholm et al., _Pennylane: Automatic differentiation of hybrid quantum-classical computations_ , 2020
  * [33] J. Bradbury et al., _JAX: composable transformations of Python+NumPy programs_ , http://github.com/google/jax, 2018
  * [34] F. Chollet et al., _Keras_ , https://keras.io, 2015
  * [35] M. Abadi et al., _TensorFlow: Large-scale machine learning on heterogeneous systems_ , https://www.tensorflow.org/, 2015
  * [36] L. Bottou, F. E. Curtis, and J. Nocedal, _Optimization methods for large-scale machine learning_ , SIAM Review 60 (2018) 223, arXiv:https://doi.org/10.1137/16M1080173
  * [37] D. P. Kingma and J. Ba, _Adam: A method for stochastic optimization_ , in 3rd International Conference on Learning Representations, ICLR 2015, San Diego, CA, USA, May 7-9, 2015, Conference Track Proceedings (Y. Bengio and Y. LeCun, eds.), 2015
  * [38] CDF, D0 collaboration, T. Wright, _$B$ -Tagging at CDF and D0, lessons for LHC_, in 17th Symposium on Hadron Collider Physics 2006 (HCP 2006), 2006, arXiv:0707.1712
  * [39] M. S. ANIS et al., _Qiskit: An open-source framework for quantum computing_ , 2021.  doi: 10.5281/zenodo.2573505
  * [40] IBM Quantum. https://quantum-computing.ibm.com/, 2021
  * [41] A. W. Cross et al., _Validating quantum computers using randomized model circuits_ , Phys. Rev. A 100 (2019) 032328
  * [42] M. Schuld, _Supervised quantum machine learning models are kernel methods_ , arXiv:2101.11020
  * [43] S. Ioffe and C. Szegedy, _Batch normalization: Accelerating deep network training by reducing internal covariate shift_ , in Proceedings of the 32nd International Conference on Machine Learning (F. Bach and D. Blei, eds.), 37, (Lille, France), 448–456, PMLR, 2015
  * [44] N. Srivastava et al., _Dropout: a simple way to prevent neural networks from overfitting_ , The journal of machine learning research 15 (2014) 1929
  * [45] S. Hochreiter and J. Schmidhuber, _Long short-term memory_ , Neural Computation 9 (1997) 1735