SciPost Phys. 14, 079 (2023)

Machine learning and LHC event generation
Anja Butter1,2 , Tilman Plehn1 , Steffen Schumann3 , Simon Badger4 , Sascha Caron5,6
Kyle Cranmer7,8 , Francesco Armando Di Bello9 , Etienne Dreyer10 , Stefano Forte11 ,Sanmay Ganguly12 , Dorival Gonçalves13 , Eilam Gross10 , Theo Heimel1 ,Gudrun Heinrich14 , Lukas Heinrich15 , Alexander Held16 , Stefan Höche17 ,Jessica N. Howard18 , Philip Ilten19 , Joshua Isaacson17 , Timo Janßen3 , Stephen Jones20 ,Marumi Kado9,21 , Michael Kagan22 , Gregor Kasieczka23 , Felix Kling24 , Sabine Kraml25 ,Claudius Krause26 , Frank Krauss20 , Kevin Kröninger27 , Rahool Kumar Barman13 ,Michel Luchmann1 , Vitaly Magerya14 , Daniel Maitre20 , Bogdan Malaescu2 ,Fabio Maltoni28,29 , Till Martini30 , Olivier Mattelaer28 , Benjamin Nachman31,32 ,Sebastian Pitz1 , Juan Rojo6,33 , Matthew Schwartz34 , David Shih25 , Frank Siegert35 ,Roy Stegeman11 , Bob Stienen5 , Jesse Thaler36 , Rob Verheyen37 ,Daniel Whiteson18 , Ramon Winterhalder28 , and Jure Zupan19

Abstract
First-principle simulations are at the heart of the high-energy physics research program. They link the vast data output of multi-purpose detectors with fundamental theory predictions and interpretation. This review illustrates a wide range of applicationsof modern machine learning to event generation and simulation-based inference, including conceptional developments driven by the specific requirements of particle physics.
New ideas and tools developed at the interface of particle physics and machine learning
will improve the speed and precision of forward simulations, handle the complexity ofcollision data, and enhance inference as an inverse simulation problem.
Copyright A. Butter et al.
This work is licensed under the Creative Commons
Attribution 4.0 International License.
Published by the SciPost Foundation.

Received 06-04-2022
Accepted 09-02-2023
Check forupdatesPublished 21-04-2023
doi:10.21468/SciPostPhys.14.4.079

1 Institut für Theoretische Physik, Universität Heidelberg, Germany
2 LPNHE, Sorbonne Université, Université Paris Cité, CNRS/IN2P3, Paris, France
3 Institut für Theoretische Physik, Georg-August-Universität Göttingen, Germany
4 Physics Department, Torino University and INFN Torino, Italy
5 IMAPP, Radboud Universiteit Nijmegen, The Netherlands
6 Nikhef Theory Group, Nikhef, Amsterdam, The Netherlands
7 Center for Cosmology and Particle Physics, New York University, New York, NY USA
8 Center for Data Science, New York University, New York, NY USA
9 Universita di Roma Sapienza, Roma, INFN, Italy
10 Weizmann Institute of Science, Rehovot, Israel
11 Dipartimento di Fisica, Università di Milano and INFN Sezione di Milano, Italy
12 ICEPP, University of Tokyo, Japan
13 Department of Physics, Oklahoma State University, Stillwater, OK, USA
14 Institut für Theoretische Physik, Karlsruher Institut für Technologie, Germany
15 Physik-Department, Technische Universität München, Germany
16 Department of Physics, New York University, New York, NY USA
17 Theoretical Physics Division, Fermi National Accelerator Laboratory, Batavia, IL, USA1

SciPost Phys. 14, 079 (2023)

18 Department of Physics & Astronomy, UC Irvine, Irvine, CA, USA
19 Department of Physics, University of Cincinnati, Cincinnati, OH, USA
20 IPPP, Physics Department, Durham University, Durham, UK
21 Université Paris-Saclay, CNRS/IN2P3, IJCLab, Orsay, France
22 Fundamental Physics Department, SLAC National Accelerator Laboratory, USA
23 Institut für Experimentalphysik, Universität Hamburg, Germany
24 Deutsches Elektronen-Synchrotron DESY, Germany
25 Univ. Grenoble Alpes, CNRS, Grenoble INP, LPSC-IN2P3, France
26 NHETC, Dept. of Physics and Astronomy, Rutgers University, Piscataway, NJ, USA
27 Department of Physics, TU Dortmund University, Germany
28 CP3, Université Catholique de Louvain, Louvain-la-Neuve, Belgium
29 Dipartimento di Fisica e Astronomia, Università di Bologna, Italy
30 Institut für Physik, Humboldt-Universität zu Berlin, Germany
31 Physics Division, Lawrence Berkeley National Laboratory, Berkeley, CA, USA
32 Berkeley Institute for Data Science, University of California, Berkeley, CA, USA
33 Department of Physics and Astronomy, Vrije Universiteit Amsterdam, The Netherlands
34 Department of Physics, Harvard University, Cambridge MA, USA
35 Institute of Nuclear and Particle Physics, Technische Universität Dresden, Germany
36 Center for Theoretical Physics, MIT, Cambridge, MA, USA
37 Department of Physics & Astronomy, University College London, UK

Contents
1 Introduction

3

2 Machine Learning in event generators
2.1 Phase space sampling
2.2 Scattering Amplitudes
2.3 Loop integrals
2.4 Parton shower
2.5 Parton distribution functions
2.6 Fragmentation functions

4
4
6
7
8
9
10

3 End-to-end ML-generators
3.1 Fast generative networks
3.2 Control and precision

11
12
14

4 Inverse simulations and inference
4.1 Particle reconstruction
4.2 Detector unfolding
4.3 Unfolding to parton level
4.4 MadMiner
4.5 Matrix element method

15
15
16
18
19
20

5 Synergies, transparency and reproducibility

21

6 Outlook

22

References

24
2

SciPost Phys. 14, 079 (2023)

1

Introduction

The defining goal of particle physics is to understand the fundamental nature of elementary
particles and their interactions. The outcome of a particle physics measurement is expressedin terms of a quantum field theory Lagrangian and its parameters. The great experimental
strength of collider-based particle physics is the availability of a huge amount of data andmeasurements in combination with a well-controlled environment. The theoretical and experimental poles are linked through precision simulations, starting from the Standard Modelor a hypothetical Lagrangian, generating particle-level events, and eventually simulating thedetector. The simulation chain realized by the standard LHC event generators [1–5] and illustrated in Fig. 1, should be based on first-principles physics rather than empiric modeling. For
these simulations precision and speed are essentially two sides of the same medal. A detailed
discussion of these traditional methods can be found in a parallel review, Ref. [6]. Adding
modern machine learning to the numerics toolbox has the potential to provide the simulations
needed for the LHC Run 3 and HL-LHC [7], as well as future energy frontier machines.
From a fundamental physics perspective there exist three distinctly different kinds of measurements at the LHC. First, basic and purely experimental measurements should be as independent of theory considerations and first-principle simulations as possible, to avoid expiration
dates. Their problem is that they provide no information about fundamental physics. These
basic measurements benefit from modern machine learning for instance in understanding thedata and calibrating the detectors. A second class of measurements is supplemented with afundamental theory interpretation framework. Examples are well-defined inclusive production
rates, like fiducial or total cross sections. They can be compared to predictions from perturbative quantum field theory. When we expect to find agreement with the Standard Model,
modern machine learning can help us in using these measurements to extract parton densities
or improve our Monte Carlo simulations. A third kind of measurement reflects our goal tofurther our understanding of fundamental physics by comparing data to predictions from perturbative or non-perturbative quantum field theory. We assume that interesting physics signals
hide in specific kinematic regions. Here, we can search for deviations between the Standard
Model predictions and experimental results, measure Standard Model parameters or higherdimensional Wilson coefficients, and aim for anomalies and eventually a proper discovery.
Such measurements of all possible features in the vast phase space of LHC collisions requireprecision simulations, specifically theory-based event generators. We will show how all ofthese aspects benefit significantly from the application of modern machine learning methods.
The challenges for event generators are, first of all, defined by the increase of the LHC
luminosity and the expected advances in experimental precision and reach. Going from theRun 2 dataset of 139 fb−1 to the projected HL-LHC dataset of 4 ab−1 suggests that experimental uncertainties at and below the percent level will become standard and need to be matchedby theory predictions, to allow for any kind of precision measurement. The same increase inrate will allow us to probe more and more exotic kinematic regions, with the hope of finding
forward

Theoryscattering

decay

QCD

showerfragmentation

detectors

EventsinverseFigure 1: Illustration of the LHC simulation chain. The forward direction is discussed
in Secs. 2 and 3, while the inverse simulation is the topic of Sec. 4.
3

SciPost Phys. 14, 079 (2023)

hints for new particles and interactions. The higher rates and an improved experimental understanding will also allow us to study more and more complex signatures with an increasing
multiplicity of backgrounds and potential signals. Each of these aspects poses a challenge tothe established event generators, and we will discuss ways modern machine learning can help
us meet them in Sec. 2. Next, we will introduce end-to-end (soup-to-nuts) ML-generators,
similar in structure to ML-detector simulations in Sec. 3. Finally, we discuss conceptual benefits from modern machine learning, for instance related to an invertible simulation chain andsimulation based inference, in Sec. 4. Because the main purpose of this report is to show new,
ML-driven developments in event generation, we refer to the main event generator Snowmass
white paper for a list of references and a detailed discussion of the physics background andthe classical approaches.

2

Machine Learning in event generators

Current multi-purpose event generators feature a modular structure, that reflects the factorization property of physics aspects at very different relevant energy scales [1–5]. While thehighest energy transfers, i.e. the hard process and QCD parton showers, can be treated byperturbative methods, phenomenological models are used to account for the hadronizationtransition, as well as non-trivial secondary interactions. The increase in perturbative precisionneeded to address the physics challenges posed by current and future collider experiments,
adds a sizeable number of more specialized numerical codes to the simulation toolbox. This
includes, for example, dedicated codes to construct and evaluate higher-order tree-level orloop amplitudes. Modern machine learning techniques can improve all aspects of event generation, ultimately making it more resource efficient and opening paths to yet more versatile
and accurate predictions. This includes important ingredients to precision predictions suchas parton densities and fragmentation functions, where neural network (NN) techniques areroutinely used already. First steps towards modeling the hadronization process with ML techniques have been presented in [8]. For the tuning of non-perturbative simulation parameters,
including an underlying event model, NN-based approaches have recently shown promise [9].

2.1Phase space sampling

The core of any scattering event simulation is the assumed hard process configuration or partonic scattering event. These are described by QFT transition amplitudes, where the physics
demands of the LHC experiments require us to consider high-multiplicity final states and oneor even two-loop QCD and/or EW corrections. The complexity of the resulting matrix elements
and the dimensionality of their phase space severely challenge the integration of cross sections
Table 1: Results for sampling the top decay width, the total cross section of top-pair
production and decay in e+ e− collisions, and g g → 3g and 4g production. Shownare the integral estimate, EN , and the unweighting efficiency, εuw , for a standard
importance sampler (Uniform), VEGAS, and NN-based optimization [10].
top decays
Sample εuw

EN [GeV]g g → 3g

top-pair production
εuw

EN [fb]εuw

Uniform 59 % 0.1679(2)
35 % 1.5254(8)VEGAS 50 % 0.16782(4) 40 % 1.5251(1)NN
84 % 0.167865(5) 78 % 1.52531(2)4EN [fb]g g → 4g
εuw

EN [fb]3.0 % 24806(55) 2.7 % 9869(20)27.7 % 24813(23) 31.8 % 9868(10)64.3 % 24847(21) 48.9 % 9859(10)SciPost Phys. 14, 079 (2023)

and the generation of partonic momentum configurations. Modern NN-techniques are ideallysuited to assist in these tasks. The standard technique used so far is based on importancesampling, employing mappings ⃗y : V → U ⊆ Rd for phase space integralsZZ∂ y(x)
f (x)
d
,
with
= g(x) .
(1)I=
d x f (x) =dd yg(x) x≡x( y)∂xVUIt can be chosen such that f /g ≈ const, to reduce the variance of the Monte Carlo integral
estimate. However, for complex matrix elements and high-dimensional phase spaces it is often
not possible to find a single function g that approximates the target function f sufficiently well.
Therefore, event generators use a multi-channel P
approach with independent
mappings ⃗yi forP
each channel i. Defining a total density g(x) = i βi g i (x), with i βi = 1 and 0 ≤ βi ≤ 1,
where βi are the channel weights, the phase space integral can be parametrized asZZXZ
f (x) X
f (x)
d
d
I=
.
(2)
d x f (x) ==d x βi g i (x)dd yi βi
g(x)
g(x)
x≡x(
y
)
V
V
U
i
i
i
i

Two ML-based approaches to phase space integration and event generation can be distinguished. The first directly hooks into existing phase space integrators and uses trainable maps
given for example by bijective normalizing flows to redistribute input random variables to themapping functions ⃗yi and better adapt to the integrand [10–16]. After an initial adaptation
phase these integrators can efficiently be used for generating weighted or unweighted events.
However, the very expressive NN-transformations can also deal with non-factorizable phasespace structures and correlations. Promising results in terms of efficiency improvements andspeed gains have been reported, see for example Tabs. 1 and 2. However, in particular for
high-multiplicity processes with non-trivial topologies the effective gains when comparing to
the established methods can fall below unity, cf. Tab 2. Therefore, next steps will be to better
combine NN-based approaches with multi-channel integrators [10, 14]. For example, one canallow the channel
weights to be phase space dependent, βi → αi (x), and solely start from thePcondition i αi (x) = 1 and 0 ≤ αi (x) ≤ 1,
Z
XZ
XZ
f (x)
d
d
dd yi αi (x)d x αi (x) f (x) =d x f (x) =I=.
(3)g i (x) x≡x( yi )
V
V
U
i
i
i

In fact, Eqs. (2) and (3) are mathematically equivalent, connected by αi (x) = βi g i (x)/g(x).
In NN-optimized event generation, Eq.(3) splits the optimization task into learning appropriate
phase space mappings for each channel and training another network to find optimal weights
αi (x) to connect all channels. This separation has two advantages: (i) possible missing correlations between the different channels can be described and recovered by the phase-space dependent channel weights, and (ii) the second network allows for a more flexible parametrization
as it does not need to be bijective.
A second approach to ML-assisted phase space sampling is based on directly learning thephase space distribution of events from input training samples, either weighted or unweighted.
Solutions employ autoregressive flows [17], generative adversarial networks (GANs) [18–21],
or variational autoencoders (VAEs) [20]. This motivates R&D to improve training throughdifferentiable programming; by merging matrix element codes with automatic differentiation [22], i.e. the automatic generation of derivatives of programs that is the backbone ofneural networks software frameworks. The gradients of matrix elements can be evaluatedand used as additional information for training generative models. Initial studies using differentiable matrix elements from MADJAX have explored extending normalizing flow training
with schemes uniquely enabled by the ability to automatically compute matrix element gradients [23], and show promise in terms of improving modeling and reducing the needed scaleof simulated datasets for training.
5

SciPost Phys. 14, 079 (2023)

10

Numerical
Analytical
Averaged models

3

normalized

Evaluation time (ms)

104

102
101
100
10

60 largest 100% ANN
50 largest 10% ANN
40 largest 1% ANN

boosted BNN
training

30
20
10

−1

0

10−2
4

5

6

0.04

0.02

0.00

(train)

0.02

0.04

Multiplicity

Figure 2: Left: comparison of the evaluation times for loop-induced amplitudes
for g g → γγ+ jets. NN-interpolation times include the averaging over ensembles for the uncertainty estimate. Figure from Ref. [29]. Right: precision∆train = |M|2NN /|M|2train − 1 for the process g g → γγ j using a Bayesian networkwith boosted training, ordered by the size of the amplitude. Figure from Ref. [31].
Closely related activities attempt to facilitate faster event unweighting and reweighting
methods using NN generative models [24, 25] or fast to evaluate NN surrogates for the transition amplitudes [26].

2.2Scattering Amplitudes

Perturbative precision calculations use interpolation methods to reduce the evaluation time for
expensive loop amplitudes, defining a task where appropriately designed neural networks canbe expected to outperform standard methods [26–30]. The challenge in NN-based surrogate
models for integrands and amplitudes is to ensure that all relevant features are indeed encodedin the network at sufficient precision and to establish a reliable uncertainty treatment of the
network training.
A relevant test case are loop-induced amplitudes such as those for
gg → ZZ

g g → γγ + jets .

and

(4)

The application of simple, gradient boosted machines to g g → Z Z highlights that fast interpolation times can lead to significant improvements in overall simulation times, if reliable models
Table 2: Unweighting efficiencies for V +jets production at the LHC. ‘SHERPA’ relies
on multi-channel importance sampling using VEGAS; ‘NN’ uses a normalizing flow;‘Gain’ shows the improvement of NN over SHERPA. Results from Ref. [14].
unweighting eff. εuw
process/sampling
W + + n jets

LO QCD
n =0

n =1

n =2

NLO QCD (RS)
n =3

n =4

n =0

n =1

SHERPA 2.8 · 10−1 3.8 · 10−2 7.5 · 10−3 1.5 · 10−3 8.3 · 10−4 9.5 · 10−2 4.5 · 10−3NNGain

6.1 · 10−1 1.2 · 10−1 1.0 · 10−3 1.8 · 10−3 8.9 · 10−4 1.6 · 10−1 4.1 · 10−3
2.2

3.3

Z/γ + n jets SHERPA 3.1 · 10
∗

−1

NN
Gain

3.6 · 101.4
−2

1.5 · 10

1.2
−2

1.1
−3

4.7 · 10

3.8 · 10−1 1.0 · 10−1 1.4 · 10−2 2.4 · 10−3
1.2

2.9

0.91

6

0.51

1.6
1.2 · 100.91
−15.3 · 10−3

1.8 · 10−3 5.7 · 10−3
1.5

1.1

SciPost Phys. 14, 079 (2023)

can be trained on fewer points than the original Monte Carlos. To control the features of theamplitude relevant for differential cross sections, separating soft and collinear regions enables
an ensemble of networks to reliably model full one-loop amplitudes for e+ e− →≤ 4 jets [28].
In Fig. 2 this scaling is shown for high-multiplicity scattering described by the NJET generator
for g g → γγ+ jets. Simulations for hadron colliders show overall improvements of around afactor Ninference /Ntraining [29] can be achieved. The right panel of Fig. 2 shows the achievable
precision on the γγ j loop amplitude from a single Bayesian network with boosted training to
improve the precision [32–35].
The reliability of the trained network is particularly at risk in divergent regions. However, these are precisely the phase space regions where the soft and collinear behavior of theamplitudes is universal and well known. Building the infrared factorization properties into
the NN-based model can lead to substantial improvements for the tree-level e+ e− → jets amplitude. Figure 3 shows that adopting a factorization-aware parametrization the achievable
precision is brought down to the per-mille level for 5-jet production [30]. While the shown
precision for this process does not translate into a clear improvement of higher-order LHC predictions, it illustrates how physics-informed network architectures can significantly improve
the network precision as the key criterion for an application in the LHC simulation chain. Perhaps of even greater interest would be the use of a single trained model to integrate over awide range of kinematic cuts, jet algorithms, PDF sets, scale choices, which could enable afurther order of magnitude in overall performance.

2.3

Loop integrals

Amplitudes beyond the leading order contain loop integrals, and machine learning can improve the calculation of (multi-)loop integrals by optimizing the integrands in Feynman parameter space [36]. When an analytic solution to these integrals is not feasible, they must be
evaluated numerically. Before attempting a numerical evaluation, the poles of the integrandneed to be controlled. In dimensional regularization, ultraviolet and infrared poles can befactorized efficiently with sector decomposition. After factorizing these poles, integrable singularities related, for example, to thresholds, remain. Such poles, located on the real axis in
Feynman parameter space, can be avoided by a deformation of the integration contour into thecomplex plane. An automated procedure to do this has already been implemented in standard
tools like SECDEC, FIESTA, and pySECDEC. The deformation of the integration contour is notunique and can be performed in many ways. In fact, the numerical precision of the integration
can vary by orders of magnitude depending on the chosen contour.

Figure 3: Accuracy for the full model from Ref. [30] for tree-level e+ e− → 5 jets
amplitudes. Figure adapted from Ref. [30].

7SciPost Phys. 14, 079 (2023)

10
10

°4

10

°5

elliptic2L (Sector 1)

pySecDec§-glob§-glob + Flow

Relative integration error

Relative integration error

triangle2L (Sector 1)
°3

10°6
10°7
10°8

101
10

10

°4

pySecDec§-glob§-glob + Flow

10°5
10°6
10°7
10°8
102

Ratio

Ratio

102

10

°3

mean: 5.20

mean: 1.6

10°1
10°2

10°1

100

101

102

103

101

mean: 2.1

10

10°1

104

s/m2

mean: 6.7

0

100

101

102

103

104

s/m2

Figure 4: Relative integration error for sector one of a 2-loop triangle integral (left)
and a 2-loop box integral known to contain elliptic functions (right) using the standard pySECDEC algorithm (green), the ML-assisted Λ-glob algorithm (blue) and including an additional normalizing flow (red). The lower panel shows the ratios tothe standard method. The figures are taken from Ref. [36].
For standard integrals, the contour deformation procedure implemented in pySECDEC
works fast and usually produces satisfactory contours in practice. However, for more complicated integrals and in specific phase-space regions, the chosen contour is sub-optimal and canbe optimized significantly, see Fig. 4. In this case, ML-assisted, or more specifically, NN-assisted
algorithms, offer great potential to amplify the precision. Like in the neural importance sampling methods [10, 13, 14, 17] for phase-space integrals, normalizing flows can be used tofind a better parametrization of the integration domain. As these contour integrals need tosatisfy certain boundary conditions, originating, for instance, from the Landau equations andCauchy’s theorem, the NN setup needs to be extended to obey these constraints. Furthermore,
the usage of complex-valued floats can entail the necessity to construct own implementations
for objects like gradients of complex determinants occurring during training and optimization.

2.4

Parton shower

The parton shower is an essential element of particle physics simulations. It describes theevolution of particles between the hard scale of the collision ∼ 100 GeV to the hadronization
scale ΛQCD . This evolution is typically modeled as a Markov process where partons evolve semiclassically, radiating gluons as they move with probabilities determined only by properties ofthe parton splitting and perhaps one or two spectator partons in the event. Although thesemi-classical approximation can be justified in the limit where the daughter particles are
emitted at small angles with respect to the mother, parton showers are used well outside of thisregime. The use of parton showers is thus justified not by physics but by necessity: computing
the full distribution from first principles is computationally intractable. This limitation is an
opportunity for machine learning; perhaps an improved parton shower could be learned rather
than built.
The simulation of parton showers offers an interesting structure compared to other generative tasks. When simulating entire collision events, as discussed in Sec. 3, commonly arepresentation encoding a small and often fixed number of 4-vectors is chosen. Simulating
showers all the way down to calorimeter sensors, or with calorimeter sensors themselves,

8

SciPost Phys. 14, 079 (2023)···

h(t 1)h(t)h(t+1)/m;?i2`b/m;?i2`b/m;?i2`b(t 1) (t 1)
kd1
kd2(t) (t)kd1 kd2(t+1) (t+1)
kd1
kd2

···
_LL, h(t)=b?Q`i?M//m;?i2`b(t)h(t 1)(t)

kd1 kd2

Figure 5: A promising approach to learning parton showers is to use a structure
inspired by the semi-classical approximation as a backbone for a general probability
estimator. In the JUNIPR approach, a recurrent neural network is used to emulated
the Markov-process aspect of a parton shower. Figure taken from Ref. [37].
yields a much larger number of particles in the final state. However, the output nodes of agenerative model can still be identified with different cells of the physical detector and therefore allow architectures that for example use convolutional layers.
Within the semi-classical approximation and even though the probability function at eachbranching in the shower is relatively simple, the overall distribution of particles produced isquite complex. It would be seriously challenging to learn this final distribution without somedomain knowledge of its structure [38]. One approach is to scaffold a learnable model overa semi-classical framework [37, 39], as sketched in Fig. 5. Additionally, network architectures
based on sets or graphs explicitly encoding permutation symmetry of the final state particles
have been investigated [40–45].
An alternative way of improving parton shower with ML-methods might be to stick to thefundamental splitting structure and measure the QCD splitting kernels in low-level observables. As before, the challenge of generating many particles covering several orders of magnitude in energy is taken care of by the usual Monte Carlo method. A modified and showerspecific form of the splitting kernels can be extracted from a combination of QCD predictions
and data using ML-based inference [46]. While this approach has practical advantages, it islimited by the applicability of the simple splittings picture.

2.5

Parton distribution functions

Parton distribution functions (PDFs), encoding the structure of colliding protons, are vitalfor the calculation of hard scattering cross sections at the LHC and appear in several stagesof the event simulation chain, in particular they guide the initial state parton showers andaffect the underlying event activity. The determination of PDFs is a classic pattern recognition
problem: it is known that an underlying law exists (the true analytic form of the PDFs, as
determined by QCD in the non-perturbative regime) but its explicit form is not known, andit must be inferred from discrete data (the cross-sections of PDF-dependent hard processes),that moreover are correlated to it indirectly and in a convoluted way. In comparison to more
standard pattern recognition problems, it has two peculiarities. First, the pattern – the setof PDFs — is a probability, rather than a deterministic outcome. Second, due to the noisy
nature of the input, which is affected by both experimental and theoretical uncertainties, witha complex correlation pattern, the final deliverable is a probability distribution of possible
results. Hence, one is delivering a probability distribution of probability distributions.

9SciPost Phys. 14, 079 (2023)

s at 1.65 GeV
Not Overfitted
Overfitted

0.08

xs(x)

0.06
0.04
0.02
0.00
0.2

0.4

x

0.6

0.8

Figure 6: Comparison between two results for the strange-quark PDF: one where
overfitting is clearly present, and one where this is not the case. Figure taken from
Ref. [47].
The way to approach this problem as a machine learning challenge was first suggested
long ago [48]: the basic idea is to deliver a Monte Carlo ensemble of machine learning models, such as neural networks, that provide the desired representation of a probability of probabilities. The successful implementation of this idea has led to the NNPDF family of proton
PDF determinations [47, 49–51] as well as to variants in the context of polarised PDF [52]and nuclear PDF [53, 54] global analyses. The current implementation frontier, which has led
to the recent NNPDF4.0 determination, involves a suite of contemporary machine learning
methods and tools, specifically cross-validation to avoid overtraining, hyperoptimization [55]combined with K-folding for the automatic selection of the methodology, feature scaling of the
input for the optimization of the neural networks used as basic underlying model [56], andGAN-enhanced compression for final efficient delivery [57, 58].
The current main challenge remains the maximal optimization of the extraction of available
information while avoiding overfitting, and the generalization to cases in which informationis scarce or altogether absent, such as extrapolation to kinematic regions where there are
no data. This is the physically most interesting case, as these are the regions where newphysics is being searched for, and also a challenge at the frontier of machine learning. While
several machine learning tools have been implemented with the aim of preventing overfitting,confirming whether the PDF resulting of a fit is indeed free of overfitting still relies – at leastin part – on the fitter’s accumulated knowledge of PDFs. To illustrate this point, Fig. 6 shows
a comparison of the strange-quark PDF xs(x, Q) at Q = 1.65 GeV, both for a good fit and aclearly overfitted alternative. The development of reliable quantitative measures of the degreeof overfitting is a challenge, both within the context of PDF determination and more in generalin machine learning, and it is a topic of ongoing research.

2.6

Fragmentation functions

Fragmentation functions (FFs) are the time-like equivalent of PDFs and encode the probabilities associated to the transition between partons produced in the hard-scattering and specifictypes of hadrons. Being based on the perturbative QCD factorization framework, FFs represent
an alternative strategy to model partonic hadronization as compared to the phenomenologicalmodels available in most MC event generators. FFs can be determined from a global analysis of10

SciPost Phys. 14, 079 (2023)

Figure 7: Fragmentation functions of up (left) and charm quarks (right) into charged
pions as a function of the time-like momentum fraction z, comparing the results of
two approaches (MAP and NNFF) based on machine learning techniques. Figure
taken from Ref. [59].
hard-scattering data from electron-positron collisions, semi-inclusive deep inelastic scattering,and proton-proton collisions (RHIC and LHC) with identified final-state hadrons.
A phenomenological analysis of FFs requires introducing a parametrization for their initial(h)scale (Q 0 ) dependence with the momentum fraction z, z Di (z, Q 0 ), where i is a partonic indexand (h) a hadronic label. To remove theory bias and model-dependence in the determinationof FFs, machine learning techniques can be adopted [59–62]. Feed-forward neural networks(h)are deployed as universal unbiased interpolants for z Di (z, Q 0 ), whose weight and thresholdparameters are obtained from a log-likelihood maximization by comparison with experimental data. This approach can be combined with the Monte Carlo replica method, originallydeployed for PDFs [63], to estimate and propagate the uncertainties from the input data to
the output FFs. The basic strategy is to generate Nrep replicas which sample the probability density associated to the data, and then train a separate neural network to each of these
replicas. The spread of the resulting networks (i.e. 68% CL intervals) provides then a robust
estimate of the uncertainties associated to the FFs.
Fig. 7 displays a comparison between FFs determined in two approaches (MAP and NNFF)
based on machine learning techniques. We show the FFs associated with the transition of up
and charm quarks into charged pions (π+ + π− ) as a function of the time-like momentum
fraction z. The bands represent the corresponding 68% CL ranges. It is worth emphasizingthat the resulting shapes, given the outcome of the NNs, are completely driven by the data,with no specific models (more or less inspired by QCD) assumed. The combination of the FFs(h)z Di (z, Q 0 ) obtained in this manner with higher-order perturbative QCD calculations provides
precise and accurate predictions for hard-scattering processes including identified hadrons in
the final state, which are important for many key phenomenological applications.

3

End-to-end ML-generators

In addition to applying a wide range of machine learning tools to improve the modules of classic event generators, we can train generative neural networks to directly generate events, atparton level and with or without detector effects [7]. End-to-end or, better, soup-to-nuts MLgenerators have to be developed together with the established generators and serve as studies
for phase space generators, enable inverted simulations, provide datasets for phenomenological analyses, and allow us to efficiently ship event samples. Their advantages include training

11

SciPost Phys. 14, 079 (2023)

on data combined with simulations, manipulation of event samples [64], or post-processing of
MC data for example to unweight events [17, 24, 26]. Finally, they define useful benchmarks
for conceptual work on uncertainty estimates for generative neural networks.
The work horses behind ML-generators are GANs, VAEs, optimal-transport-based probabilistic autoencoders, normalizing flows, and their invertible network (INN) variant. Given
the interpolation properties of neural networks and the benefits of their implicit bias in theapplications described in Sec. 2, we can quantify the amplification of statistics-limited training
data through generative networks [65, 66].

3.1

Fast generative networks

Theory-driven ML-generators at the parton level [18, 69] can be combined with experimentdriven fast detector simulations [70–81] into single generative networks [68,82–86], provided
we have sufficient control over the network and its uncertainties. Single, soup-to-nuts simulation networks are inspired by the fundamental goal of the detection process, namely toreconstruct parton-level information as accurately as possible.
Comparing different generative network architectures, we start with highly expressive
VAEs. They can be trained to generate events at the parton level, without or with fast detector simulation, by maximizing a lower bound of the data likelihood through variational
inference (ELBO). The model consists of a decoder p(x|z) which maps from a latent space
Z to the phase space X , and an encoder q(z|x) which is a variational approximation to theinverse of p(x|z). In practice, it is difficult to simultaneously optimize the separate components of the ELBO and the VAE performance can be improved by weighting the KL-divergencein the loss function term by a factor β. The B-VAE [20] is characterized by the limit β ≪ 1
and a strong preference for the reconstruction loss. After optimization, the Gaussian latent
distribution is replaced by a buffer which consists of the latent distribution derived from training events. This model simultaneously achieves a highly-optimized reconstruction loss, butwith a closely-matched and non-Gaussian latent distribution. While VAEs are very expressive
probabilistic models, the approximate nature of the ELBO and the need to balance the twocomponents of the loss function can become limiting factors.
Similarly, GANs can extract and reproduce the phase space density of LHC events. While
×10−2
True
No Weights
Flow

3
2
1

0.12
True
No Weights
Flow

0.10
dσ
dpT,tt [pb/GeV]dσdmtt [pb/GeV]

4

0.08
0.06
0.04
0.02

1.2

1.2

1.0

1.0

0.8

400

500

600

700
mtt [GeV]

800

900

0.8

1000

0

50

100

150
pT,tt [GeV]

200

250

Figure 8: Distributions of the invariant mass (left) and transverse momentum (right)of the t t̄ system in pp → t t̄ generated using MC@NLO. The true distribution (red)is compared with the normalizing flow distributions excluding (blue) or including(green) negative event weights. Figure from Ref. [67].

12

300

SciPost Phys. 14, 079 (2023)

Figure 9: Visualization of the transformation from parton-level theory-space to reconstructed data-space, Z → X . Left: learned transformation of OTUS’s decoder,p D (x | z). Right: true transformation from simulated samples, for comparison. Colors in the X projection indicate the source bin in Z for a given sample. Figure from
Ref. [68].
technically the difference between training on events without or with detector effects is negligible, parton-level events are more challenging when it comes to sharp kinematic features like
Breit-Wigner mass peaks. GANs generically do not achieve the necessary precision for suchfeatures, so they have to be enhanced, for example with a targeted MMD loss [18]. The main
challenge of GANs is the precision they can achieve in the underlying phase space density
while finding a Nash equilibrium.
Finally, normalizing flows avoid some of the limitations of the above architectures for LHC
event generation [67, 69, 87]. At the cost of some flexibility, they offer a direct evaluation of
the likelihood without having to resort to variational inference. They start from a latent space
Z and apply a series of bijective transforms, with tractable Jacobian, to the phase space X .
While the expressivity of the model may in some cases be limited, the advantage of a tractable
likelihood is significant. Flows can be trained on weighted events, including negative weights,
through a simple modification of the loss [67]. Figure 8 illustrates their performance applied to
pp → t t̄ events at the parton level, including shower evolution, generated with MC@NLO. In
addition, normalizing flows come with significant advantages in controlling their performanceand quantifying uncertainties, as discussed in the next section. Their invertible structure isuseful for many LHC-applications, including anomaly detection or related density estimation
tasks [88–91].
An attractive application of soup-to-nuts networks can be targeted using Optimal Transport-based probabilistic autoencoders [68]. Their structural advantage is that they learn the
mapping from parton-level information in theory space, Z , to detected and reconstructedobjects in data space, X , without requiring paired event samples, {z, x}. The probabilistic
autoencoder’s latent space is identified with a physically meaningful representation of partonlevel theory-space information, so the encoder and decoder networks define a simulator mapping, Z → X , and an unfolding mapping, X → Z . Properties of the OT-based method encourage the encoder and decoder to be conditional mappings, effectively sampling from theprobability distributions pE (z|x) and pD (x|z), respectively. Over many samples, these distributions will marginalize to the appropriate theory-space and data-space priors, p(z) and p(x),13

SciPost Phys. 14, 079 (2023)

normalized

Z + 1 jet exclusive
10−3
10−5

Reweighted
Train

Conditioned
Truth

BINN
Truth

wD

δ[%]

10.0
1.0
0.1
1.1
1.0
0.9
1.2
1.0

2

a ∈ [0, 6, 12]

1
0

50

100
pT,µ1 [GeV]150

Figure 10: Illustration of a complete control and uncertainty treatment for generative
networks applied to LHC event generation and simulation. Figure from Ref. [69].
respectively. Alternative methods to encode an unfolding mapping in neural networks are
discussed in Sec. 4.
Despite having no training pair information, OTUS’s learned mappings exhibit physicalbehavior, even picking up on invariant masses which were withheld during training. This suggests that further development in this direction should produce physically meaningful mappings, even if relations are missed or unknown, and therefore not included in the training
process. On the other hand, providing known relations as inductive biases on the data inputs,
network architectures, or loss functions will likely improve performance. Figure 9 depicts thejoint distribution and marginals of OTUS’ trained simulator as well as the true simulator for
one test-case. Despite OTUS only having information about marginal-matching during training, the decoder network learns a mapping which is qualitatively similar to the true simulator.

3.2

Control and precision

If we use neural networks to encode theory predictions for the LHC, we need to ensure that
all relevant phase space features are described with the required precision [92]. For neuralnetworks, this problem can be split into two distinct tasks: first, we need control over the
relevant phase space features, so the network does not interpolate over relevant, but narrow
phase space regions. Second, we have to estimate the precision with which the network haslearned the underlying phase space density.
Neural networks work much like a fit and not like an interpolation in the sense that they
do not reproduce the training data faithfully and instead learn a smooth approximation [65,66]. This is where we can gain some intuition for a NN-uncertainty treatment. For a fit,

14

SciPost Phys. 14, 079 (2023)

uncertainties on the training data are crucial information in the loss function. We then monitor
the fit quality and ensure that the fit is reliable over the entire phase space.
To guarantee that all relevant features are encoded in a generative network, we can follow
the GAN inspiration and train a simple discriminator network to identify and quantify deviations between training and generated data. As a post-processing step such a discriminator canbe used to reweight the events from the generative network [25, 69, 93]. In the GAN spirit wecan incorporate the discriminator into the generator training, either through adversarial training searching for a Nash equilibrium, or through alternative approaches for a normalizing flow
generator. Such a joint training will improve the generator, provide an uncertainty estimate,and prepare any remaining information in the discriminator for reweighting, as illustrated forZ+jets production at the parton level in Fig. 10.
Once we know that the neural network describes all features, we determine how wellit does. This can be done with Bayesian networks, where the learned network weights are
replaced by learned network weight distributions [32,33]. Bayesian network approaches havebeen shown to describe uncertainties in regression [35] and classification [34] tasks, and theconcept can be expanded to generative networks [87]. For generative networks we can assign
a training-related uncertainty in the underlying phase space density to the (unit) weight ofeach event. In Fig. 10 we see, for instance, the increasing uncertainty in the kinematic tail,driven by a lack of training data.
We often know systematic or theoretical limitations of describing certain kinematic
regimes. In that case we augment the training data, representing this uncertainty throughan additional parameter in event weights. We train the generative network conditionally on
this parameter, either in a deterministic or a Bayesian setup, and generate events either for agiven parameter or sampling over it. Again, this approach is illustrated in Fig. 10, where a directly affects the p T -distribution of the leading jet and enters many other observables throughkinematic correlations. We see that its effect is larger than the uncertainty from the Bayesian
network for the individual a-values. This first attempt of a comprehensive uncertainty treatment for generative networks will allow us to build confidence in the applications of generative
networks to LHC simulation and inference.

4

Inverse simulations and inferenceMonte Carlo simulations based on first principles have allowed us to properly understandessentially all aspects of LHC data. The price to pay for an extremely fast and reliable forward
Monte Carlo simulation chain is that the corresponding inverse simulation is not feasible inpractice. ML-based simulations can be built symmetrically, for instance INNs encode a bijective
mapping between two physics spaces linking different levels of the simulation chain illustrated
in Fig. 1 [94, 95]. Similarly, we can relate different levels of the simulation chain througha reweighting procedure working on the full respective phase spaces and accounting for all
correlations [96]. Moreover, as ML-based simulations are often differentiable, we can use
their gradients to probe and learn about distributions on phase space [97]. Finally, we canconstruct generative inverse simulations with conditional versions of the respective forward
generative networks [95, 98, 99]. This last approach is based on progress with soup-to-nuts
ML-generators and their essentially identical network architectures.

4.1

Particle reconstruction

The first stage of the inverse problem uses the set of energy deposits in the detector to reconstruct the set of particles present at the first interaction with the detector, that is, following
hadronization. In its fullest sense, reconstruction also involves the prediction of the particles’
15

SciPost Phys. 14, 079 (2023)

properties, in particular, their class and momenta. The difficulty of this task stems from thebusy LHC environment caused by pileup interactions and the inherently collimated signaturesassociated with jets. Traditional particle flow (PF) algorithms rely on parameterized schemes
for merging and splitting to disentangle overlapping calorimeter cell clusters as well as trackbased subtraction to infer the contribution from neutral particles.
A series of publications [100–102] have established the potential for ML-based reconstruction to go beyond traditional PF algorithms. In Ref. [100], particle reconstruction was recast
as a computer vision problem using state-of-the-art ML architectures including U-net, graphneural network (GNN) and DeepSets. A simplified dataset was used comprising overlapping
pairs of charged and neutral pions in a 6-layer calorimeter block. In comparison to a traditional PF algorithm, the ML models regress the component of neutral energy better by afactor of two to four in terms of resolution. The study also finds significant improvements via
a super-resolution approach (see also Ref. [103]), where the network is trained to predict a
corresponding calorimeter signature with higher granularity.
This proof of concept has been extended to particle reconstruction in more realistic environments resembling multiple pileup interactions in a full-coverage simulated detector [101,102].
In both cases, GNN architectures are employed for their ability to handle the complexity ofdetector data: variable numbers of input and target entities, lack of ordering, irregularity ofdetector components, and sparsity of “pixels”. Moreover, GNNs are able to leverage the spatial
relationships between calorimeter cells alongside their input features to optimize the prediction tasks.
Based on these developments, it can already be anticipated that ML methods will take a keyrole in particle reconstruction at future runs of the LHC, especially to handle HL-LHC conditions. GNN-based models in particular show potential to outperform current PF algorithms forparticle identification and regression while opening new possibilities such as super-resolution
and resolving neutral particles inside of jets. Finally, the learned deep latent representation of
detector information, which underlies the prediction tasks, should serve as a more expressive
input format for both event classification and downstream tasks in the inverse problem.

4.2

Detector unfolding

While the physical processes behind an LHC collision are described by fundamental physics
and are therefore universal, the observed data depend in an intimate way on the technical details of the detector. Detector effects like phase space coverage, detection thresholds, particle
reconstruction, efficiencies, or calibration induce not only resolution smearing in the measurements, but can lead to systematic deviations between the properties of particles reaching thedetector and the objects reconstructed from actually measured data. For individual experiments, these detector effects differ greatly and can only be estimated by the collaboration. It
is therefore essential for future interpretations of a measurement to unfold detector effects
so that we can compare measurements by different experiments to each other and to theory
predictions.
Traditional approaches to unfolding are based on matrices connecting binned particle-level
distributions at truth level with histograms of corresponding detector-level observables. While
the folding or convolution of detector effects with kinematic distributions at particle level ispossible with Monte Carlo simulations, the inverse direction often suffers from instabilities andscales poorly for high-dimensional phase spaces. The limitation to low-dimensional representations requires an unwanted pre-selection of interesting observables. Finally, the matrix-based
approach requires fixed bin sizes, which limits the re-optimization options for future analyses.
ML-approaches establish high-dimensional and binning-independent unfolding. We candistinguish two fundamentally different concepts [104]: a classification-based approach toreweight a Monte Carlo simulation with the learned likelihood ratio of data and simula16

10
8
6

Truth: Pythia 8.243
Truth (alt.): Herwig 8.1.5
Detector: Delphes 3.4.2 (CMS)
Z+jet: pTZ > 200 GeV, R = 0.4

4
2

0
1.50
1.25
1.00
0.75
0.00 0.05 0.10 0.15 0.20 0.25 0.30
0.00

Ratio to Truth

Number of Events / 104

Truth
OmniFold
cINN

0.05

0.10

0.15

0.20

0.25

7
6
5
4
3
2
1
0

Truth
OmniFoldcINN
Truth: Pythia 8.243
Truth (alt.): Herwig 8.1.5
Detector: Delphes 3.4.2 (CMS)
Z+jet: pTZ > 200 GeV, R = 0.4

0.00

0.25

0.50

0.75

1.00

1.25

1.50

1.75

2.00

0.30

Ratio to Truth

Number of Events / 104

SciPost Phys. 14, 079 (2023)

2

1.50
1.25
1.00
0.75
0.0

0.5

1.0

21 = 2/ 1

1.5

2.0

Figure 11: An illustration of classifier-based (OmniFold) and density-based (cINN)
unfolding for the two-dimensional space of N -subjettiness variables τ1 and τ2 . The
right plot shows that since the unfolding is done simultaneously and unbinned, we
can produce a measurement of the widely-used N -subjettiness ratio for free. Figure
from Ref. [104].
tion [96]; and a complementary approach that learns directly the probability density atparticle-level [95]. Figure 11 illustrates that both approaches can perform unbinned unfolding
in multiple dimensions.
Classification-based approaches start by learning the likelihood ratio between data andsimulation at reconstruction level [96, 105, 106]. Using matched event pairs at the truth andreconstruction levels, the resulting weights are pulled to the particle level. Next, a classifier
learns the likelihood ratio of the weighted and unweighted distributions at particle level toreplace an event-based weight with a generalized weighting function. After several iterationsof weight updates, the algorithm converges to an unfolded distribution which is compatible
with the observed measurement.
Density-based approaches build on generative networks that predict probability configurations of truth-level events given a detector-level measurement. They are trained on pairsof reconstruction- and particle-level events from Monte Carlo simulations, to learn a direct
mapping between both levels. Unfolding built on generative adversarial networks has beenshown to work on kinematic distributions [94, 99]. Event-wise unfolding requires a meaningful probabilistic treatment, which can be achieved with conditional normalizing flows [95,97],
the kind of generative networks which also allows for the uncertainty treatment discussed inSec. 3.2. This unfolding method yields calibrated probability distributions for each measured
event. It admits multiple approaches; one approach frames unfolding in terms of learning a
conditional density of particle-level quantities conditioned on reconstructed inputs [95], while
another approach frames unfolding as an empirical Bayes / maximum marginal likelihood /data-informed prior learning problem [97].
Because classification-based and density-based unfolding techniques have distinct
strengths and weaknesses, the natural next step will be to combine the two methods to benefit
from both. While there is an extensive R&D program required to integrate both methodologies and to achieve precision, these tools are starting to be applied to data analysis in collider
physics [107]. Looking ahead, it is clear that future versions of these tools will play an important role in the data analysis of future colliders. Unfolded differential cross sections are
one of the main data products from collider experiments. By performing the unfolding with asmuch information as possible, we ensure that the measurements achieve the maximal preci17

SciPost Phys. 14, 079 (2023)

×10−2

2 jet incl.
Parton Truth
Parton cINN
Detector Truth

1 dσ
σ dpT,W[GeV−1 ]

1.0
0.8
0.6
0.4

cINN
Truth

0.2

10
1
dσ
−1
σ dMW,reco [GeV ]

1.2

0.0
1.2
1.0
0.8
0

50

100

150 200 250
pT,W [GeV]

300

350

10−2

10−3

10−4

400

cINN vs Truth
2 jet
3 jet
4 jet

−1

60

65

70

75
80
85
MW,reco [GeV]

90

95

100

Figure 12: Unfolded parton-level distributions for the process pp → Zℓℓ W j j +jets
using a cINN. The unfolding covers detector effects as well as additional jet radiation.
Figure from Ref. [95].
sion, making the best use of the data. Furthermore, high-dimensional and unbinned unfolding
ensures that these data products are ‘future proof’ in the sense that binning and even observables can be chosen post-hoc [104]. This enables downstream data analysis long after the data
were published, including when new theoretical insights are available.

4.3

Unfolding to parton level

Once we control ML-unfolding of detector effects, we can target other parts of the simulation
chain shown in Fig. 1 and invert them for a given LHC analysis. To probe the kinematics ofa hard scattering process we can use neural networks to unfold QCD jet radiation and heavy
particle decays to study the production kinematics of top quarks, electroweak gauge bosons, orthe Higgs without binning and with full correlations. Such measurements are standard in top
physics and provide the ideal input to global SMEFT analyses. Once we know the parton-level
configurations for a given observed event, we can use NN-techniques to evaluate observables
like CP-sensitive angular correlations in their original reference frames.
The inversion of QCD radiation or decays relies on the same classification or generative
networks as detector unfolding. For instance, we can train a normalizing flow to map random
numbers to the parton-level phase space, under the condition of a given detector-level event.
The underlying model is encoded in the forward simulation chain used to train the network.
Part of it is the assumed hard process, including the number of jets which are part of the hard
scattering and do not get unfolded. When analyzing an event and sampling into parton-level
phase space, we extract a probability distribution of parton-level configurations [95], which
we can use to define observables suitable for standard analyses.
One challenge for such analyses are combinatorics. For the hard scattering qq̄ → Zℓℓ W j j
and up to two additional QCD jets we ask how well cINN-unfolding extracts the W -kinematics.
In the left panel of Fig. 12 we illustrate how the network reproduces the momentum of the
decaying W -boson. The relation between the up to four jets and the two partonic quarks fromthe W -decay is learned by the network. In the right panel of Fig. 12 we show the reconstructed W -mass stacked for different numbers of jets. The network resolves the underlying
combinatorics such that the W -widths for the different jet multiplicities are identical, all by byaccessing correlations combined with the truth information from the forward simulation. This
corresponds to results from a systematic study which shows that deep networks outperform
classical approaches to solving the combinatorics in the reconstruction of top-quark final statessignificantly [108].
18

SciPost Phys. 14, 079 (2023)

High-level observables encoded into neural networks will find their way into standard
experimental analyses. They are motivated by existing top-sector measurements, and using
NN-techniques will simplify their use considerably. Moreover, the comprehensive uncertainty
treatment discussed Sec. 3.2 and the merged classification-based [96] and density-based techniques from Sec. 4.2 can be applied to any part of an inverted or unfolded simulation chain.

4.4

MadMiner

The relation between data x and physics parameters θ is, fundamentally, described by the
likelihood function or normalized fully differential cross section, which we can predict in afactorized form,p(x|θ ) =1 dσ(x|θ ).
σ(θ ) dx

(5)

While we can predict this likelihood at the detector level using the standard, forward simulation tools, we can only compute it in a closed form at the parton level. This challenge in therelation of simulations and inference is where neural networks might lead to transformative
progress.
Inspired by the standard simulation chain we can assume that the likelihood of Eq.(5)approximately factorizes into the form [109, 110]ZZZp(x|θ ) =

dzd

dzs

dz p p(x|zd ) p(zd |zs ) p(zs |z p ) p(z p |θ ) .
|
{z
}

(6)

p(x, z|θ )Here we integrate over latent variables z, where zd characterize the detector effects, zs the parton shower and hadronization, and z p the partonic phase space including helicities, charges,
and flavors, etc. Given the typically large number of latent variables, it is unrealistic to integrate over them or evaluate the joint-likelihood p(x, z|θ ). However, it is possible to calculate
the joint likelihood ratio relative to a reference point in terms of the ratio of squared matrix
elements from parton-level generators [109–113],p(z p |θ )|M |2 (z p |θ ) σ(θref )p(x, z|θ )r(x, z|θ ) ==∼.
p(x, z|θref )p(z p |θref ) |M |2 (z p |θref ) σ(θ )(7)The starting point to new ML-methods is to construct functionals in terms of the joint likelihood
ratio r(x, z|θ ), which are minimized by the true likelihood or likelihood ratio function [114,115]. The result of this training are neural networks that approximate the true likelihood
ratio r(x|θ ). Given such a neural network, established statistical techniques can be used toconstruct confidence limits in parameter space.
Note that here simulation-based inference provides the primary statistical model, i.e. theprobability model p(x|θ , ν) that describes the dependence on the data x, the parameters ofinterest θ , and the nuisance parameters ν, even when the data is high-dimensional and traditional modeling approaches are inadequate. The publication of the trained network in a
re-usable form, as discussed below, can thus be of great benefit for an optimal use of experimental results [117]. This approach is separate from tools like DNNLikelikood [118], whichaims at approximating likelihoods derived from traditional approaches to model building, using libraries like RooFit.
Instead of the full likelihood function, one can also use the score t(x|θ ) = ∇θ log p(x|θ )to define statistically optimal observables at the detector level. This approach is motivated by
an expansion of the log likelihood ratio around θref ,log r(x|θ ) = log r(x|θref ) + t(x|θref ) (θ − θref ) + · · · .
19(8)SciPost Phys. 14, 079 (2023)

Rate
pT,
SALLY
ALICES

1.4

0.3

1.2

0.0

1.0

Searches in pp tth @ HL-LHC (68% CL)t

100 cG

0.6

0.3

0.8

0.6

0.6
0.4

0.2

0.0
cu

0.2

Di-leptonic tth
Semi-leptonic tth

60

0.4

30

Hadronic tth
Combination

0

+30

+60

Figure 13: Expected sensitivity of a MadMiner based analysis for t t̄ H production,
probing SMEFT coefficients (left) and the CP-structure of the top Yukawa coupling(right). Figures taken from Ref. [113] and Ref. [116], respectively.
For parameter points close enough to θref the score components are the sufficient statistics,
so for measuring θ knowing t(x|θref ) is as powerful as the full likelihood. Since the scoreis defined through the likelihood function, it is also intractable. However, similarly to theapproach discussed above, we can train a neural network on a suitable loss function such that
it will converge to the score. The trained network will now represent the optimal observable. In
a next step, the likelihood can be determined for instance with simple histograms of the scorecomponents [113, 119]. This approach requires only minor changes to established analysis
pipelines. Alternatively, the scores can be used to evaluate the Fisher Information and set
limits based on the Cramer-Rao bound [113]. One challenge with training using the score isthat the relevant gradient information of the matrix elements must be accessible for training
the neural network, but this information is typically only accessible for a subset of parameterswith analytic dependence that facilitates easy gradient estimation. One approach to enable
score based training for any parameter is through differentiable programming; when matrix
elements are merged with automatic differentiation frameworks, the required gradients canbe computed automatically with relatively small additional computational overhead. Casestudies using differentiable matrix elements from MADJAX for score based training successfully
trained networks for inference on parameters that were inaccessible without differentiablematrix elements [23].
The previously outlined inference strategy has been fully automated in the MadMiner
tool [113, 120]. The increase in physics sensitivity relative to a total rate or single kinematic
distribution is illustrated in Fig. 13. In the left panel we consider t t̄ H production to constrain
the two SMEFT Wilson coefficients cu and cG . In the right panel we consider the same process
to constrain CP-violation in the top-Higgs coupling, as parameterized by the magnitude κ t and
CP-phase α of the top Yukawa coupling [116, 121].

4.5

Matrix element method

Inverting the entire simulation chain in Fig. 1 allows us to extract the transition amplitude foran observed event and relate it to the theory prediction. This so-called matrix element method(MEM) can be used to estimate fundamental physics parameters from individual events andhas, for instance, been applied for measure the top mass. Being defined on the single-event
level, it is in particular suitable for low-statistics signals, where an optimal exploitation of all
20

SciPost Phys. 14, 079 (2023)

kinematic features is critical.
The MEM relies on our ability to extract the likelihood for detector-level events as a function of a model parameter θ , as an ML-application through density-based unfolding or an
inverted simulation. Extending the discussion in Sec. 4.4, the transition amplitude as a function of detector-level phase space is unknown, but it can be calculated at the parton level. The
two phase spaces can be related by transfer functions T (⃗x , z⃗), probabilities to observe partonlevel configurations z⃗ as detector-level signatures x⃗ , as part of the forward simulation. In the
N -event likelihood they appear asZNY
dl σ(θ )1
dm σ(θ )1=
L(θ ) =
p(⃗x |θ ) =dm z
T (⃗x (i) , z⃗) . (9)
σ
(θ
)
dx
...
dx
σ
(θ
)
dz
...
dz
fid
1
l x⃗ (i)
fid
1
m
i=1
i=1
i=1
N
Y

(i)

N
Y

The dimensionality of the parton-level and detector-level phase spaces is different. For instance, longitudinal neutrino momenta are unobservable, while additional jets have to be included with higher-order QCD corrections. Existing approaches model the transfer functions
heuristically, and for non-trivial cases the numerical convolution is impossible. The form ofEq.(9) indicates ways of enhancing the accuracy of the matrix element method: first, higherorder corrections can be included at parton level, for instance using the MEM@NLO program.
Second, general and highly non-Gaussian transfer functions can account for parton shower,
hadronization, detector resolution, acceptance, and efficiency, as well as a possible mismatch
between theoretically described and actually measured quantities, event by event.
The transfer function is defined as a probability density T (⃗x , z⃗) = p(⃗x |⃗z , θ ). This allows us
to learn it directly from simulated data using a conditional normalizing flow or INN as a density estimator. Because the matrix element spans several orders of magnitude and the transfer
function usually is a narrow peak in phase space, the integral in Eq.(9) is numerically challenging for a regular Monte-Carlo integration. However, we know from Bayes’ theorem thatthe integrand becomes trivial when the parton level samples are drawn from the distributionp(⃗z |⃗x , θ ). In analogy to Secs. 4.2 and 4.3, the partonic configurations z⃗ for a given detector
event x⃗ (i) can be sampled by another conditional INN. Then the likelihood can be expressed
as
N
Y

­


·
1
∂ z⃗(⃗r ; x⃗ (i) , θ ) dm σ(θ )(i)L(θ ) =T (⃗x , z⃗).
σfid (θ )∂ ⃗r
dz1 ... dzmz⃗(⃗r ;⃗
x (i) , θ ) ⃗r ∼p(⃗r )
i=1

(10)

This way, density estimation of the transfer function in combination with density-based unfolding will allow us to make optimal use of the statistical power of the MEM, exploiting thefull and correlated event kinematics event by event for critical LHC observables like the top
mass, the Higgs self-coupling, or CP-violating phases.

5

Synergies, transparency and reproducibility

A key paradigm in the development of simulation tools for high-energy collider experiments ispublicly accessible open source software. The versioning of code releases and the reproducibility of predictions is vital for a reliable analysis and interpretation of collider data. As we haveseen in the previous sections, ML-methods are entering all aspects of the simulations chain athigh pace. They range from initial proof-of-concept applications to well established use caseswith largely consolidated techniques, for example in the determination of parton densities.
Machine learning models efficiently encode arbitrary decision functions of a given set of
inputs, and thus offer a chance to easily exchange complex relations. This might correspondto the value of a scattering matrix element given a set of momenta, or probability models
21

SciPost Phys. 14, 079 (2023)

from simulation-based inference, like MEM or MadMiner. The sharing of neural networks
used for various generative or discriminative tasks will be of central importance and should
be further extended. This will allow researchers to critically examine and build upon previous
results more easily, enable synergies between different use cases, and facilitate reproducibilityof results.
Successfully sharing a machine learning model entails two challenges: (i) sharing themodel itself, including architectures, software versions, and weights; and (ii) sharing data
it can be used on. Exchanging models is technically relatively straightforward and several
corresponding tools exist, for example Open Neural Network Exchange (ONNX). It allows theexchange of neural networks and BDTs between training frameworks.
Suitable input data poses the more difficult problem. On the side of results by large collaborations, this adds additional weight to the ongoing move towards publishing open data
along with measurements. Containerization, as enabled by software tools like Docker can be
useful in bundling the correct versions of different software packages used for data processing
and machine learning in a coherent fashion.
An opportunity exists in the realm of phenomenological studies based on the DELPHES
detector simulation [122]. Here a common specification on how quantities are translated into
the inputs to machine learning algorithms might — together with publishing the ML models —boost sharing and meaningful exchange. Another interesting angle are generative models. As
these do not need data to evaluate, sharing the architecture and weights is already sufficient.
Generative networks themselves can even be used as an efficient alternative way of sharing
simulated data.
Publication of ML models for their reuse is not yet standard in the particle physics community. Examples where trained networks have been published in ONNX format for futurereuse are the DNNLikelihood [123], a package for cross-disciplinary training of discriminator
networks [124], and the ATLAS search for R-parity-violating supersymmetry [125, 126], thelatter also being available in the ATLAS SimpleAnalysis framework. However, detailed documentation for instance of the input variables is missing. Further development is strongly
encouraged for, e.g., the purpose of analysis preservation [127, 128], and in general for the
implementation of the Findable, Accessible, Interoperable, and Reusable (FAIR) principles forscientific data management [129] of ML models. An example for a dataset with special emphasis on these aspects can be found in Ref. [130]. Making the newly developed simulation
and analysis tools, along with the required data, accessible to other scientists and future users
forms an essential element of open and thriving science.

6

OutlookAs a field combining vast datasets with excellent, first-principle simulations, particle physics
is benefiting tremendously from developments in data science and machine learning. While
new AI-inspired methods will not magically solve all challenges in LHC simulations and analysis, they are providing a crucial and transformative contribution to our numerical toolbox.
Moreover, given the quality of the LHC datasets, simulations, and simulation-based analysis
methods, we expect particle physics to eventually contribute to broader machine learning research.
Event generation, or the simulation of signals for the LHC detectors from QFT Lagrangians,
is the main link between experimental and theoretical particle physics. It has stringent requirements when it comes to first principles vs modeling, control, precision, speed, and flexibility. In
this review we have shown that even within the physics-motivated modular structure of standard event generators, there is no aspect that cannot be improved through modern machine

22

SciPost Phys. 14, 079 (2023)

learning. This includes phase space sampling, scattering amplitudes, loop integrals, parton
showers, parton densities, and fragmentation. Some of these ML-applications have a long
history and are accepted as standard approaches, other ML-based improvements of physics
modules are currently under rapid development and are finding their way into standard generators. All of them will be key to address the needs for example of the HL-LHC.
In addition to ML-enhanced event generators, an interesting application of generative neural networks are ML-generators at parton level and fast ML-detector simulations. They provide
an excellent testing ground for phase space generators, precision networks, and inverted simulations. This includes conceptual developments in the field of generative networks, driven by
LHC-specific requirements of controlling precision-generative networks as numerical tools andproviding a full range of uncertainties. They allow us to define, produce, and encode datasetsfor phenomenological studies and serve as a compression for data entering experimental analyses.
The main conceptual advantage of ML-event generation is that simulations with generative
networks are symmetric: given a fundamental physics model we can predict the probability
distributions of LHC events over phase space, or we can predict the probability distributions ofmodel parameters given observed LHC events. Different ML-approaches to simulation-based
inference include classification-based methods, conditional generative networks as a direct
inversion, or indirect ways of learning likelihood ratios. In combination, they will allow us tosystematically use unfolding or inverted simulations at the HL-LHC, from particle identificationand detector unfolding all the way to an event-wise matrix element method analysis.
Finally, there are many simulation-related questions in fundamental physics, where AImethods allow us to make significant progress. Examples going beyond immediate applications to event generation include symbolic regression [131], sample and data compression [58, 132], detection of symmetries [133–136], and many other fascinating new ideas
and concepts.

Acknowledgments
Anja Butter, Gudrun Heinrich, and Tilman Plehn are supported by the Deutsche Forschungsgemeinschaft under grant 396021762 - TRR 257 Particle Physics Phenomenology after the
Higgs Discovery. Lukas Heinrich is supported by the Excellence Cluster ORIGINS, which isfunded by the Deutsche Forschungsgemeinschaft under Germany’s Excellence Strategy - EXC2094-390783311. Alexander Held is supported by the U.S. National Science Foundation (NSF)
Cooperative Agreement OAC-1836650. Stefan Höche and Joshua Isaacson are supported by
the Fermi National Accelerator Laboratory (Fermilab), a U.S. Department of Energy, Office of
Science, HEP User Facility. Fermilab is managed by Fermi Research Alliance, LLC (FRA), acting
under Contract No. DE–AC02–07CH11359. Jessica N. Howard is supported by the U.S. National Science Foundation under the grant DGE-1839285. Michael Kagan is supported by the
US Department of Energy (DOE) under grant DE-AC02-76SF00515. Gregor Kasieczka and Felix Kling are supported by the Deutsche Forschungsgemeinschaft under Germany’s Excellence
Strategy – EXC 2121 Quantum Universe – 390833306. Sabine Kraml acknowledges support
by the IN2P3 master project Théeorie – BSMGA and the joint ANR-FWF project PRCI SLDNP
grant no. ANR-21-CE31-0023. Claudius Krause and David Shih are supported by DOE grant
DOE-SC0010008. Rahool Kumar Barman and Dorival Gonçalves thank the U.S. Departmentof Energy for the financial support under grant number DE-SC 0016013. Benjamin Nachman
is supported by the U.S. Department of Energy, Office of Science under contract DE-AC0205CH11231. Tilman Plehn is supported by the Deutsche Forschungsgemeinschaft under Germany’s Excellence Strategy EXC 2181/1 - 390900948 (the Heidelberg STRUCTURES Excel-

23

SciPost Phys. 14, 079 (2023)

lence Cluster). Steffen Schumann acknowledges support from the German Federal Ministryof Education and Research (BMBF, grant 05H21MGCAB) and the German Research Foundation (DFG, project number 456104544). Rob Verheyen is supported by the European ResearchCouncil (ERC) under the European Union’s Horizon 2020 research and innovation programme(grant agreement No. 788223, PanScales). Ramon Winterhalder is supported by FRS-FNRS(Belgian National Scientific Research Fund) IISN projects 4.4503.16.

References[1] T. Sjöstrand et al., An introduction to PYTHIA 8.2, Comput. Phys. Commun. 191, 159(2015), doi:10.1016/j.cpc.2015.01.024.
[2] E. Bothmann et al., Event generation with Sherpa 2.2, SciPost Phys. 7, 034 (2019),
doi:10.21468/SciPostPhys.7.3.034.
[3] J. Alwall et al., The automated computation of tree-level and next-to-leading order differential cross sections, and their matching to parton shower simulations, J. High Energy
Phys. 07, 079 (2014), doi:10.1007/JHEP07(2014)079.
[4] J. Bellm et al., Herwig 7.0/Herwig++ 3.0 release note, Eur. Phys. J. C 76, 196 (2016),
doi:10.1140/epjc/s10052-016-4018-8.
[5] W. Kilian, T. Ohl and J. Reuter, WHIZARD: Simulating multi-particle processes at LHC and
ILC, Eur. Phys. J. C 71, 1742 (2011), doi:10.1140/epjc/s10052-011-1742-y.
[6] J. M. Campbell et al., Event generators for high-energy physics experiments, (arXiv
preprint) doi:10.48550/arXiv.2203.11110.
[7] A. Butter and T. Plehn, Generative networks for LHC events, in Artificial Intelligencefor High Energy Physics, World Scientific, Singapore, ISBN 9789811234026 (2022),
doi:10.1142/9789811234033_0007.
[8] P. Ilten, T. Menzo, A. Youssef and J. Zupan, Modeling hadronization using machine learning, SciPost Phys. 14, 027 (2023), doi:10.21468/SciPostPhys.14.3.027.
[9] M. Lazzarin, S. Alioli and S. Carrazza, MCNNTUNES: Tuning shower Monte Carlogenerators with machine learning, Comput. Phys. Commun. 263, 107908 (2021),
doi:10.1016/j.cpc.2021.107908.
[10] E. Bothmann, T. Janßen, M. Knobbe, T. Schmale and S. Schumann, Exploring phase space with neural importance sampling, SciPost Phys. 8, 069 (2020),
doi:10.21468/SciPostPhys.8.4.069.
[11] J. Bendavid, Efficient Monte Carlo integration using boosted decision trees and generative
deep neural networks, (arXiv preprint) doi:10.48550/arXiv.1707.00028.
[12] M. D. Klimek and M. Perelstein, Neural network-based approach to phase space integration, SciPost Phys. 9, 053 (2020), doi:10.21468/SciPostPhys.9.4.053.
[13] C. Gao, J. Isaacson and C. Krause, i-flow: High-dimensional integration and sampling with normalizing flows, Mach. Learn.: Sci. Technol. 1, 045023 (2020),
doi:10.1088/2632-2153/abab62.

24

SciPost Phys. 14, 079 (2023)[14] C. Gao, S. Höche, J. Isaacson, C. Krause and H. Schulz, Event generation with normalizing flows, Phys. Rev. D 101, 076002 (2020), doi:10.1103/PhysRevD.101.076002.
[15] I.-K. Chen, M. D. Klimek and M. Perelstein, Improved neural network Monte Carlo simulation, SciPost Phys. 10, 023 (2021), doi:10.21468/SciPostPhys.10.1.023.
[16] S. Pina-Otey, F. Sánchez, T. Lux and V. Gaitan, Exhaustive neural importance sampling applied to Monte Carlo event generation, Phys. Rev. D 102, 013003 (2020),
doi:10.1103/PhysRevD.102.013003.
[17] B. Stienen and R. Verheyen, Phase space sampling and inference from weighted events withautoregressive flows, SciPost Phys. 10, 038 (2021), doi:10.21468/SciPostPhys.10.2.038.
[18] A. Butter, T. Plehn and R. Winterhalder, How to GAN LHC events, SciPost Phys. 7, 075(2019), doi:10.21468/SciPostPhys.7.6.075.
[19] R. Di Sipio, M. F. Giannelli, S. K. Haghighat and S. Palazzo, DijetGAN: A generativeadversarial network approach for the simulation of QCD dijet events at the LHC, J. HighEnergy Phys. 08, 110 (2019), doi:10.1007/JHEP08(2019)110.
[20] S. Otten, S. Caron, W. de Swart, M. van Beekveld, L. Hendriks, C. van Leeuwen, D.
Podareanu, R. Ruiz de Austri and R. Verheyen, Event generation and statistical samplingfor physics with deep generative models and a density information buffer, Nat. Commun.
12, 2985 (2021), doi:10.1038/s41467-021-22616-z.
[21] S. Choi and J. H. Lim, A data-driven event generator for hadron colliders using
Wasserstein generative adversarial network, J. Korean Phys. Soc. 78, 482 (2021),
doi:10.1007/s40042-021-00095-1.
[22] A. G. Baydin, B. A. Pearlmutter, A. A. Radul and J. M. Siskind, Automatic differentiation in machine learning: A survey, J. Mach. Learn. Res. 18, 1 (2018),
doi:10.48550/arXiv.1502.05767.
[23] L. Heinrich and M. Kagan, Differentiable Matrix Elements with MadJax, J. Phys.: Conf.
Ser. 2438, 012137 (2023), doi:10.1088/1742-6596/2438/1/012137.
[24] M. Backes, A. Butter, T. Plehn and R. Winterhalder, How to GAN event unweighting,
SciPost Phys. 10, 089 (2021), doi:10.21468/SciPostPhys.10.4.089.
[25] R. Winterhalder, M. Bellagente and B. Nachman, Latent space refinement for deep generative models, (arXiv preprint) doi:10.48550/arXiv.2106.00792.
[26] K. Danziger, T. Janßen, S. Schumann and F. Siegert, Accelerating Monte Carlo event generation – Rejection sampling using neural network event-weight estimates, SciPost Phys.
12, 164 (2022), doi:10.21468/SciPostPhys.12.5.164.
[27] F. Bishara and M. Montull, Machine learning amplitudes for faster event generation, Phys.
Rev. D 107, L071901 (2023), doi:10.1103/PhysRevD.107.L071901.
[28] S. Badger and J. Bullock, Using neural networks for efficient evaluation ofhigh multiplicity scattering amplitudes, J. High Energy Phys. 06, 114 (2020),
doi:10.1007/JHEP06(2020)114.
[29] J. A.-Bullock, S. Badger and R. Moodie, Optimising simulations for diphoton production at hadron colliders using amplitude neural networks, J. High Energy Phys. 08, 066(2021), doi:10.1007/JHEP08(2021)066.
25

SciPost Phys. 14, 079 (2023)[30] D. Maître and H. Truong, A factorisation-aware Matrix element emulator, J. High Energy
Phys. 11, 066 (2021), doi:10.1007/JHEP11(2021)066.
[31] S. Badger, A. Butter, M. Luchmann, S. Pitz and T. Plehn, Loop amplitudes from a boosted
precision network, (arXiv preprint) doi:10.48550/arXiv.2206.14831.
[32] Y. Gal, Uncertainty in Deep Learning, Ph.D. thesis, Cambridge, UK (2016).
[33] A. Kendall and Y. Gal, What uncertainties do we need in Bayesian deep learning for computer vision?, (arXiv preprint) doi:10.48550/arXiv.1703.04977.
[34] S. Bollweg, M. Haussmann, G. Kasieczka, M. Luchmann, T. Plehn and J. Thompson, Deep-learning jets with uncertainties and more, SciPost Phys. 8, 006 (2020),
doi:10.21468/SciPostPhys.8.1.006.
[35] G. Kasieczka, M. Luchmann, F. Otterpohl and T. Plehn, Per-object systematics using deeplearned calibration, SciPost Phys. 9, 089 (2020), doi:10.21468/SciPostPhys.9.6.089.
[36] R. Winterhalder, V. Magerya, E. Villa, S. Jones, M. Kerner, A. Butter, G. Heinrich and T.
Plehn, Targeting multi-loop integrals with neural networks, SciPost Phys. 12, 129 (2022),
doi:10.21468/SciPostPhys.12.4.129.
[37] A. Andreassen, I. Feige, C. Frye and M. D. Schwartz, Binary JUNIPR: An interpretable probabilistic model for discrimination, Phys. Rev. Lett. 123, 182001 (2019),
doi:10.1103/PhysRevLett.123.182001.
[38] Y. S. Lai, D. Neill, M. Płoskoń and F. Ringer, Explainable machine learning of the underlying physics of high-energy particle collisions, Phys. Lett. B 829, 137055 (2022),
doi:10.1016/j.physletb.2022.137055.
[39] A. Andreassen, I. Feige, C. Frye and M. D. Schwartz, JUNIPR: A framework for
unsupervised machine learning in particle physics, Eur. Phys. J. C 79, 102 (2019),
doi:10.1140/epjc/s10052-019-6607-9.
[40] E. Bothmann and L. del Debbio, Reweighting a parton shower using a neural network: The
final-state case, J. High Energy Phys. 01, 033 (2019), doi:10.1007/JHEP01(2019)033.
[41] L. de Oliveira, M. Paganini and B. Nachman, Learning particle physics by example:
Location-aware generative adversarial networks for physics synthesis, Comput. Softw. Big
Sci. 1, 4 (2017), doi:10.1007/s41781-017-0004-6.
[42] J. W. Monk, Deep learning as a parton shower, J. High Energy Phys. 12, 021 (2018),
doi:10.1007/JHEP12(2018)021.
[43] K. Dohi, Variational autoencoders
doi:10.48550/arXiv.2009.04842.

for

jet

simulation,(arXivpreprint)[44] B. Orzari, T. Tomei, M. Pierini, M. Touranakou, J. Duarte, R. Kansal, J.-R. Vlimant andD. Gunopulos, Sparse data generation for particle-based simulation of hadronic jets in the
LHC, (arXiv preprint) doi:10.48550/arXiv.2109.15197.
[45] S. Tsan, R. Kansal, A. Aportela, D. Diaz, J. Duarte, S. Krishna, F. Mokhtar, J.-R. Vlimant
and M. Pierini, Particle graph autoencoders and differentiable, learned energy mover’s
distance, (arXiv preprint) doi:10.48550/arXiv.2111.12849.

26

SciPost Phys. 14, 079 (2023)[46] S. Bieringer, A. Butter, T. Heimel, S. Höche, U. Köthe, T. Plehn and S. T. Radev,
Measuring QCD splittings with invertible networks, SciPost Phys. 10, 126 (2021),
doi:10.21468/SciPostPhys.10.6.126.
[47] R. D. Ball et al., The path to proton structure at 1% accuracy, Eur. Phys. J. C 82, 428(2022), doi:10.1140/epjc/s10052-022-10328-7.
[48] S. Forte, L. S. Garrido, J. I. Latorre and A. Piccione, Neural network parametrization of deep-inelastic structure functions, J. High Energy Phys. 05, 062 (2002),
doi:10.1088/1126-6708/2002/05/062.
[49] R. D. Ball, L. del Debbio, S. Forte, A. Guffanti, J. I. Latorre, J. Rojo and M. Ubiali, A first
unbiased global NLO determination of parton distributions and their uncertainties, Nucl.
Phys. B 838, 136 (2010), doi:10.1016/j.nuclphysb.2010.05.008.
[50] R. D. Ball, Parton distributions for the LHC run II, J. High Energy Phys. 04, 040 (2015),
doi:10.1007/JHEP04(2015)040.
[51] R. D. Ball et al., An open-source machine learning framework for global analyses of parton
distributions, Eur. Phys. J. C 81, 958 (2021), doi:10.1140/epjc/s10052-021-09747-9.
[52] E. R. Nocera, R. D. Ball, S. Forte, G. Ridolfi and J. Rojo, A first unbiased global determination of polarized PDFs and their uncertainties, Nucl. Phys. B 887, 276 (2014),
doi:10.1016/j.nuclphysb.2014.08.008.
[53] R. Abdul Khalek, R. Gauld, T. Giani, E. R. Nocera, T. R. Rabemananjara and J. Rojo,nNNPDF3.0: Evidence for a modified partonic structure in heavy nuclei, Eur. Phys. J. C
82, 507 (2022), doi:10.1140/epjc/s10052-022-10417-7.
[54] R. Abdul Khalek, J. J. Ethier, J. Rojo and G. van Weelden, nNNPDF2.0: Quark
flavor separation in nuclei from LHC data, J. High Energy Phys. 09, 183 (2020),
doi:10.1007/JHEP09(2020)183.
[55] S. Carrazza and J. M. Cruz-Martinez, Towards a new generation of parton densities withdeep learning models, Eur. Phys. J. C 79, 676 (2019), doi:10.1140/epjc/s10052-0197197-2.
[56] S. Carrazza, J. M. Cruz-Martinez and R. Stegeman, A data-based parametrization ofparton distribution functions, Eur. Phys. J. C 82, 163 (2022), doi:10.1140/epjc/s10052022-10136-z.
[57] S. Carrazza, J. I. Latorre, J. Rojo and G. Watt, A compression algorithm for the combination of PDF sets, Eur. Phys. J. C 75, 474 (2015), doi:10.1140/epjc/s10052-015-3703-3.
[58] S. Carrazza, J. M. Cruz-Martinez and T. R. Rabemananjara, Compressing PDF
sets using generative adversarial networks, Eur. Phys. J. C 81, 530 (2021),
doi:10.1140/epjc/s10052-021-09338-8.
[59] R. A. Khalek, V. Bertone and E. R. Nocera, Determination of unpolarized pion fragmentation functions using semi-inclusive deep-inelastic-scattering data, Phys. Rev. D 104,
034007 (2021), doi:10.1103/PhysRevD.104.034007.
[60] V. Bertone, S. Carrazza, N. P. Hartland, E. R. Nocera and J. Rojo, A determination of
the fragmentation functions of pions, kaons, and protons with faithful uncertainties, Eur.
Phys. J. C 77, 516 (2017), doi:10.1140/epjc/s10052-017-5088-y.

27

SciPost Phys. 14, 079 (2023)[61] V. Bertone, N. P. Hartland, E. R. Nocera, J. Rojo and L. Rottoli, Charged hadronfragmentation functions from collider data, Eur. Phys. J. C 78, 651 (2018),
doi:10.1140/epjc/s10052-018-6130-4.
[62] M. Soleymaninia, H. Hashamipour and H. Khanpour, Neural network QCD analysis ofcharged hadron fragmentation functions in the presence of SIDIS data, Phys. Rev. D 105,
114018 (2022), doi:10.1103/PhysRevD.105.114018.
[63] T. NNPDF collaboration, L. Del Debbio, S. Forte, J. I. Latorre, A. Piccione and J. Rojo,p
Unbiased determination of the proton structure function F2 with faithful uncertainty estimation, J. High Energy Phys. 03, 080 (2005), doi:10.1088/1126-6708/2005/03/080.
[64] A. Butter, T. Plehn and R. Winterhalder, How to GAN event subtraction, SciPost Phys.
Core 3, 009 (2020), doi:10.21468/SciPostPhysCore.3.2.009.
[65] A. Butter, S. Diefenbacher, G. Kasieczka, B. Nachman and T. Plehn, GANplifying event
samples, SciPost Phys. 10, 139 (2021), doi:10.21468/SciPostPhys.10.6.139.
[66] S. Bieringer et al., Calomplification — the power of generative calorimeter models, J.
Instrum. 17, P09028 (2022), doi:10.1088/1748-0221/17/09/P09028.
[67] B. Stienen and R. Verheyen, Phase space sampling and inference from weighted events withautoregressive flows, SciPost Phys. 10, 038 (2021), doi:10.21468/SciPostPhys.10.2.038.
[68] J. N. Howard, S. Mandt, D. Whiteson and Y. Yang, Learning to simulate high energy
particle collisions from unlabeled data, Sci. Rep. 12, 7567 (2022), doi:10.1038/s41598022-10966-7.
[69] A. Butter, T. Heimel, S. Hummerich, T. Krebs, T. Plehn, A. Rousselot and S.
Vent, Generative networks for precision enthusiasts, SciPost Phys. 14, 078 (2023),
doi:10.21468/SciPostPhys.14.4.078.
[70] M. Paganini, L. de Oliveira and B. Nachman, Accelerating science with generative adversarial networks: An application to 3D particle showers in multilayer calorimeters, Phys.
Rev. Lett. 120, 042003 (2018), doi:10.1103/PhysRevLett.120.042003.
[71] M. Paganini, L. de Oliveira and B. Nachman, CaloGAN: Simulating 3D high energy particle showers in multilayer electromagnetic calorimeters with generative adversarial networks, Phys. Rev. D 97, 014021 (2018), doi:10.1103/PhysRevD.97.014021.
[72] P. Musella and F. Pandolfi, Fast and accurate simulation of particle detectors using generative adversarial networks, Comput. Softw. Big Sci. 2, 8 (2018), doi:10.1007/s41781018-0015-y.
[73] M. Erdmann, L. Geiger, J. Glombitza and D. Schmidt, Generating and refining particle detector simulations using the Wasserstein distance in adversarial networks, Comput.
Softw. Big Sci. 2, 4 (2018), doi:10.1007/s41781-018-0008-x.
[74] M. Erdmann, J. Glombitza and T. Quast, Precise simulation of electromagnetic calorimeter
showers using a Wasserstein generative adversarial network, Comput. Softw. Big Sci. 3, 4(2019), doi:10.1007/s41781-018-0019-7.
[75] ATLAS Collaboration, Energy resolution with a GAN for fast shower simulation in ATLAS,
Https://atlas.web.cern.ch/Atlas/GROUPS/PHYSICS/PLOTS/SIM-2019-004/

28

SciPost Phys. 14, 079 (2023)[76] D. Belayneh et al., Calorimetry with deep learning: Particle simulation and reconstruction
for collider physics, Eur. Phys. J. C 80, 688 (2020), doi:10.1140/epjc/s10052-020-82519.
[77] E. Buhmann, S. Diefenbacher, E. Eren, F. Gaede, G. Kasieczka, A. Korol and K. Krüger,Getting high: High fidelity simulation of high granularity calorimeters with high speed,Comput. Softw. Big Sci. 5, 13 (2021), doi:10.1007/s41781-021-00056-0.
[78] E. Buhmann, S. Diefenbacher, E. Eren, F. Gaede, G. Kasieczka, A. Korol and K. Krüger,Decoding photons: Physics in the latent space of a BIB-AE generative network, Eur. Phys.
J. Web Conf. 251, 03003 (2021), doi:10.1051/epjconf/202125103003.
[79] C. Chen, O. Cerri, T. Q. Nguyen, J. R. Vlimant and M. Pierini, Analysis-specific fast
simulation at the LHC with deep learning, Comput. Softw. Big Sci. 5, 15 (2021),
doi:10.1007/s41781-021-00060-4.
[80] C. Krause and D. Shih, CaloFlow: Fast and accurate generation of calorimeter showers
with normalizing flows, (arXiv preprint) doi:10.48550/arXiv.2106.05285.
[81] C. Krause and D. Shih, CaloFlow II: Even faster and still accurate generation of calorimeter
showers with normalizing flows, (arXiv preprint) doi:10.48550/arXiv.2110.11377.
[82] S. Otten, S. Caron, W. de Swart, M. van Beekveld, L. Hendriks, C. van Leeuwen, D.
Podareanu, R. R. de Austri and R. Verheyen, Event generation and statistical samplingfor physics with deep generative models and a density information buffer, Nat. Commun.
12, 2985 (2021), doi:10.1038/s41467-021-22616-z.
[83] B. Hashemi, N. Amin, K. Datta, D. Olivito and M. Pierini, LHC analysis-specific datasetswith Generative Adversarial Networks, (arXiv preprint) doi:10.48550/arXiv.1901.05282.
[84] R. di Sipio, M. F. Giannelli, S. K. Haghighat and S. Palazzo, DijetGAN: A generativeadversarial network approach for the simulation of QCD dijet events at the LHC, J. HighEnergy Phys. 08, 110 (2019), doi:10.1007/JHEP08(2019)110.
[85] J. Arjona Martínez, T. Q. Nguyen, M. Pierini, M. Spiropulu and J.-R. Vlimant, Particle
generative adversarial networks for full-event simulation at the LHC and their application
to pileup description, J. Phys.: Conf. Ser. 1525, 012081 (2020), doi:10.1088/17426596/1525/1/012081.
[86] Y. Alanazi et al., Simulation of electron-proton scattering events by a feature-augmented
and transformed generative adversarial network (FAT-GAN) in Proceedings of the thirtiethinternational joint conference on artificial intelligence, ISBN 9780999241196 (2021),
doi:10.24963/ijcai.2021/293.
[87] M. Bellagente, M. Haussmann, M. Luchmann and T. Plehn, Understanding
event-generation networks via uncertainties, SciPost Phys. 13, 003 (2022),
doi:10.21468/SciPostPhys.13.1.003.
[88] S. Caron, L. Hendriks and R. Verheyen, Rare and different: Anomaly scores from a combination of likelihood and out-of-distribution models to detect new physics at the LHC,
SciPost Phys. 12, 077 (2022), doi:10.21468/SciPostPhys.12.2.077.
[89] A. Hallin, J. Isaacson, G. Kasieczka, C. Krause, B. Nachman, T. Quadfasel, M. Schlaffer,
D. Shih and M. Sommerhalder, Classifying anomalies through outer density estimation,
Phys. Rev. D 106, 055006 (2022), doi:10.1103/PhysRevD.106.055006.
29

SciPost Phys. 14, 079 (2023)[90] B. Nachman and D. Shih, Anomaly detection with density estimation, Phys. Rev. D 101,
075042 (2020), doi:10.1103/PhysRevD.101.075042.
[91] T. Buss, B. M. Dillon, T. Finke, M. rämer, A. Morandini, A. Mück, I. Oleksiyuk and T.
Plehn, What’s anomalous in LHC jets?, (arXiv preprint) doi:10.48550/arXiv.2202.00686.
[92] B. Nachman, A guide for deploying deep learning in LHC searches: How toachieve optimality and account for uncertainty, SciPost Phys. 8, 090 (2020),
doi:10.21468/SciPostPhys.8.6.090.
[93] S. Diefenbacher, E. Eren, G. Kasieczka, A. Korol, B. Nachman and D. Shih, DCTRGAN:
Improving the precision of generative models with reweighting, J. Instrum. 15, P11004(2020), doi:10.1088/1748-0221/15/11/P11004.
[94] K. Datta, D. Kar and D. Roy, Unfolding with generative adversarial networks, (arXiv
preprint) doi:10.48550/arXiv.1806.00433.
[95] M. Bellagente, A. Butter, G. Kasieczka, T. Plehn, A. Rousselot, R. Winterhalder, L. Ardizzone and U. Köthe, Invertible networks or partons to detector and back again, SciPost
Phys. 9, 074 (2020), doi:10.21468/SciPostPhys.9.5.074.
[96] A. Andreassen, P. T. Komiske, E. M. Metodiev, B. Nachman and J. Thaler, OmniFold: A
method to simultaneously unfold all observables, Phys. Rev. Lett. 124, 182001 (2020),
doi:10.1103/PhysRevLett.124.182001.
[97] M. Vandegar, M. Kagan, A. Wehenkel and G. Louppe, Neural empirical Bayes: Source
distribution estimation and its applications to simulation-based inference, in Proceedingsof the 24th international conference on artificial intelligence and statistics, AAAI Publications, Palo Alto, USA, ISBN 9781577357384 (2021).
[98] L. Ardizzone, C. Lüth, J. Kruse, C. Rother and U. Köthe, Guided image generation with conditional invertible neural networks, (arXiv preprint)
doi:10.48550/arXiv.1907.02392.
[99] M. Bellagente, A. Butter, G. Kasieczka, T. Plehn and R. Winterhalder, How to GAN away
detector effects, SciPost Phys. 8, 070 (2020), doi:10.21468/SciPostPhys.8.4.070.
[100] F. A. di Bello, S. Ganguly, E. Gross, M. Kado, M. Pitt, L. Santi and J. Shlomi, Towards acomputer vision particle flow, Eur. Phys. J. C 81, 107 (2021), doi:10.1140/epjc/s10052021-08897-0.
[101] J. Pata, J. Duarte, J.-R. Vlimant, M. Pierini and M. Spiropulu, MLPF: Efficient machinelearned particle-flow reconstruction using graph neural networks, Eur. Phys. J. C 81, 381(2021), doi:10.1140/epjc/s10052-021-09158-w.
[102] S. R. Qasim, K. Long, J. Kieseler, M. Pierini and R. Nawaz, Multi-particle reconstructionin the high granularity calorimeter using object condensation and graph neural networks,
Eur. Phys. J. Web Conf. 251, 03072 (2021), doi:10.1051/epjconf/202125103072.
[103] P. Baldi, L. Blecher, A. Butter, J. Collado, J. N. Howard, F. Keilbach, T. Plehn, G. Kasieczka
and D. Whiteson, How to GAN higher jet resolution, SciPost Phys. 13, 064 (2022),
doi:10.21468/SciPostPhys.13.3.064.
[104] M. Arratia et al., Publishing unbinned differential cross section results, J. Instrum. 17,
P01024 (2022), doi:10.1088/1748-0221/17/01/P01024.

30

SciPost Phys. 14, 079 (2023)[105] A. Andreassen, P. T. Komiske, E. M. Metodiev, B. Nachman, A. Suresh and J. Thaler,Scaffolding simulations with deep learning for high-dimensional deconvolution, (arXiv
preprint) doi:10.48550/arXiv.2105.04448.
[106] P. Komiske, W. P. McCormack and B. Nachman, Preserving new physics while
simultaneously unfolding all observables, Phys. Rev. D 104, 076027 (2021),
doi:10.1103/PhysRevD.104.076027.
[107] V. Andreev et al., Measurement of lepton-jet correlation in deep-inelastic scattering with theH1 detector using machine learning for unfolding, Phys. Rev. Lett. 128, 132002 (2022),
doi:10.1103/PhysRevLett.128.132002.
[108] J. Erdmann, T. Kallage, K. Kröninger and O. Nackenhorst, From the bottom to the
top–reconstruction of t t̄ events with deep learning, J. Instrum. 14, P11015 (2019),
doi:10.1088/1748-0221/14/11/P11015.
[109] J. Brehmer, K. Cranmer, G. Louppe and J. Pavez, Constraining effective
field theories with machine learning, Phys. Rev. Lett. 121, 111801 (2018),
doi:10.1103/PhysRevLett.121.111801.
[110] J. Brehmer, K. Cranmer, G. Louppe and J. Pavez, A guide to constraining effective field theories with machine learning, Phys. Rev. D 98, 052004 (2018),
doi:10.1103/PhysRevD.98.052004.
[111] J. Brehmer, G. Louppe, J. Pavez and K. Cranmer, Mining gold from implicit models to improve likelihood-free inference, Proc. Natl. Acad. Sci. 117, 5242 (2020),
doi:10.1073/pnas.1915980117.
[112] M. Stoye, J. Brehmer, G. Louppe, J. Pavez and K. Cranmer, Likelihood-free inference with
an improved cross-entropy estimator, (arXiv preprint) doi:10.48550/arXiv.1808.00973.
[113] J. Brehmer, F. Kling, I. Espejo and K. Cranmer, MadMiner: Machine learning-based inference for particle physics, Comput. Softw. Big Sci. 4, 3 (2020), doi:10.1007/s41781-0200035-2.
[114] K. Cranmer, J. Pavez and G. Louppe, Approximating likelihood ratios with calibrated
discriminative classifiers, (arXiv preprint) doi:10.48550/arXiv.1506.02169.
[115] P. Baldi, K. Cranmer, T. Faucett, P. Sadowski and D. Whiteson, Parameterized neural networks for high-energy physics, Eur. Phys. J. C 76, 235 (2016), doi:10.1140/epjc/s10052016-4099-4.
[116] R. K. Barman, D. Gonçalves and F. Kling, Machine learning the Higgs boson-top quark
CP phase, Phys. Rev. D 105, 035023 (2022), doi:10.1103/PhysRevD.105.035023.
[117] K. Cranmer et al., Publishing statistical models: Getting the most out of particle physics
experiments, SciPost Phys. 12, 037 (2022), doi:10.21468/SciPostPhys.12.1.037.
[118] A. Coccaro, M. Pierini, L. Silvestrini and R. Torre, The DNNLikelihood: Enhancing likelihood distribution with deep learning, Eur. Phys. J. C 80, 664 (2020),
doi:10.1140/epjc/s10052-020-8230-1.
[119] J. Brehmer, K. Cranmer and F. Kling, Improving inference with matrix elements and machine learning, Int. J. Mod. Phys. A 35, 2041008 (2020),
doi:10.1142/S0217751X20410080.

31

SciPost Phys. 14, 079 (2023)[120] J. Brehmer, F. Kling, I. Espejo, S. Perez and K. Cranmer, MadMiner: ML based inference
for particle physics, Zenodo (2019), doi:10.5281/zenodo.1489147.
[121] H. Bahl and S. Brass, Constraining CP -violation in the Higgs-top-quark interaction using machine-learning-based inference, J. High Energy Phys. 03, 017 (2022),
doi:10.1007/JHEP03(2022)017.
[122] J. de Favereau, C. Delaere, P. Demin, A. Giammanco, V. Lemaître, A. Mertens and M.
Selvaggi, DELPHES 3: A modular framework for fast simulation of a generic collider experiment, J. High Energy Phys. 02, 057 (2014), doi:10.1007/JHEP02(2014)057.
[123] A. Coccaro, M. Pierini, L. Silvestrini and R. Torre, The DNNLikelihood: Enhancing likelihood distribution with deep learning, Zenodo (2019), doi:10.5281/zenodo.3567822.
[124] L. Benato et al., Shared data and algorithms for deep learning in fundamental physics,
Comput. Softw. Big Sci. 6, 9 (2022), doi:10.1007/s41781-022-00082-6.
[125] G. Aad et al., Search for R-parity-violating supersymmetry in a final state containing lepptons and many jets with the ATLAS experiment using s = 13 TeV proton-proton collision
data, Eur. Phys. J. C 81, 1023 (2021), doi:10.1140/epjc/s10052-021-09761-x.
[126] G. Aad et al., “SUSY-2019-04-ONNX.tgz” of “Search for R-parity-violating supersymmetrypin a final state containing leptons and many jets with the ATLAS experiment using s =13
TeV proton-proton collision data” (2021), doi:10.17182/hepdata.104860.v1/r3.
[127] W. Abdallah et al., Reinterpretation of LHC results for new physics: Status and recommendations after Run 2, SciPost Phys. 9, 022 (2020), doi:10.21468/SciPostPhys.9.2.022.
[128] S. Bailey et al., Data and analysis preservation, recasting, and reinterpretation, (arXiv
preprint) doi:10.48550/arXiv.2203.10057.
[129] M. D. Wilkinson et al., The FAIR guiding principles for scientific data management andstewardship, Sci. Data 3, 160018 (2016), doi:10.1038/sdata.2016.18.
[130] Y. Chen et al., A FAIR and AI-ready Higgs boson decay dataset, Sci. Data 9, 31 (2022),
doi:10.1038/s41597-021-01109-0.
[131] A. Butter, T. Plehn, N. Soybelman and J. Brehmer, Back to the formula – LHC edition,(arXiv preprint) doi:10.48550/arXiv.2109.10414.
[132] A. Butter, S. Diefenbacher, G. Kasieczka, B. Nachman, T. Plehn, D. Shih and R. Winterhalder, Ephemeral learning - Augmenting triggers with online-trained normalizing flows,
SciPost Phys. 13, 087 (2022), doi:10.21468/SciPostPhys.13.4.087.
[133] P. Betzler and S. Krippendorf, Connecting dualities and machine learning, Fortschr. Phys.
68, 2000022 (2020), doi:10.1002/prop.202000022.
[134] S. Krippendorf and M. Syvaeri, Detecting symmetries with neural networks, Mach. Learn.:
Sci. Technol. 2, 015010 (2021), doi:10.1088/2632-2153/abbd2d.
[135] G. Barenboim, J. Hirn and V. Sanz, Symmetry meets AI, SciPost Phys. 11, 014 (2021),
doi:10.21468/SciPostPhys.11.1.014.
[136] K. Desai, B. Nachman and J. Thaler, Symmetry discovery with deep learning, Phys. Rev.
D 105, 096031 (2022), doi:10.1103/PhysRevD.105.096031.

32