Accelerating HEP simulations with Neural Importance Sampling
Authors: Nicolas Deutschmann
Niklas Götz
Venue: ICLR 2022 Submitted (2022)
Keywords: Importance Sampling
Normalizing Flows
High-Energy-Physics

Abstract: Virtually all high-energy-physics (HEP) simulations for the LHC rely on Monte Carlo using importance sampling by means of the VEGAS algorithm. However, complex high-precision calculations have become a challenge for the standard toolbox.
As a result, there has been keen interest in HEP for modern machine learning to power adaptive sampling. Despite previous work proving that normalizing-flow-powered neural importance sampling (NIS) sometimes outperforms VEGAS, existing research has still left major questions open, which we intend to solve by introducing ZüNIS, a fully automated NIS library.
We first show how to extend the original formulation of NIS to reuse samples over multiple gradient steps, yielding a significant improvement for slow functions. We then benchmark ZüNIS over a range of problems and show high performance with limited fine-tuning. This is crucial for ZüNIS to be a mature tool for the wider HEP public. We outline how the the library allows for non-experts to employ it with minimal effort, an essential condition to widely assess the value of NIS for LHC simulations.

ACCELERATING HEP SIMULATIONS WITH NEURAL IMPORTANCE SAMPLING
Anonymous authors Paper under double-blind review
ABSTRACT
Virtually all high-eneryy-phystee (HEP) simulations for the LHC rely Monte Carlo using importance sampling by means the VEGAS algorithm. How- ever, complex high-precision calculations have become challenge for the stan- dard toolbox. As a result, there has been keen interest in HEP for modern ma- chine learning to power adaptive sampling. Despite previous work proving that normalizing-foo flow-powered neural importance sampling (NIS) sometimes outper- forms VEGAS, existing research has still left major questions open, which we intend to solve by introducing ZüNIS, a fully automated NIS library. We first show how to extend the original formulation of NIS to reuse samples over multi- ple gradient steps, yielding a significant improvement for slow functions. We then benchmark ZüNIS over a range of problems and show high performance with limited fine finettuiin. The library can be used by non-experts with minimal effort, which is crucial to become a mature tool for the wider HEP public.
INTRODUCTION
Beghh-Energy-Pyyeee (HEP) simulations are at the heart of the Large Hadron Collider (LHC) pro- gram for studying the fundamental laws of nature. Most HEP predictins are expressed as expectation values, evaluated numerically as Monte Carlo (MC) integrals. This permits both the integration of the very complex functions and the reproduction of the data selection process by experiments.
Most HEP simulations tools (Alwall et al., 2014) perform MC integrals using importance sampling, which allows to adaptively sample points to speed up convergence while keeping independent and identically distributed samples, crucial to reproduce experimental analyses which can only ingest uniformly-weighted data, typically produced by rejection sampling (see appendix A).
The most popular tool to optimize importance sampling in HEP by far the VEGAS algorithm (Lep- age, 1980; 2021), which fights the curse of dimensionality by assuming correlations between the variables. While this is rarely the case in general, a good understanding of the integrand func- tion can help significantly. Indeed optimized parametrizations using multichannelling (Kleiss al., 1986; Ohl, 1999; Kleiss & Pittau, 1994) have become bread-and -nd-eettter tools for HEP event genera- tion simulators, with good success for leading-order (I.O) calculations However, simulations get more complex, either by having more complex final states or by including higher orders in pertur- bation theory, performance degrades fast.
While newer approaches to adaptive importance sampling have been developped, like Population Monte Carlo (Bugallo et al., 2017: 2015; Cappé et al., 2004; Iba, 2001) and its extensions (Cappé et al., 2008; Koblents & Míguez, 2015; Elvira et al., 2017; Douc et al., 2007; Cornuet et al., 2012), and have been successfully applied other fields, VEGAS remains the standard tool used by essen- tially all modern simulation tools (Alwall et al., 2014; Bothmann et al., 2019; Bellm et al., 2016). We suspect that the reasons for this are twofold. First of all the lack of exchange between the two scientific communities is probably at cause and should probably not be neglected. Second of all however, the impetus in the HEP community is to move toward more "black box" approaches where little is known of the structure of the integrand, which seems to be in tension with the situation with PMC, which is sensitive to to the initial proposal density (Beaujean & Caldwell, 2013; Cappé et al., 2008). As the main practical goal of this paper is reduce the reliance of HEP simulations the careful tuning of integrators, we will focus on comparing our work with the de facto HEP standard: the VEGAS algorithm.

---

There is much room for investing computational time into improving sampling (ATLAS Collabora- tion, 2020): modern HEP theoretical calculations are taking epic proportions and can require hours for a single function evaluation (Jones, 2018). Furthermore, unweighting samples can be extremely inefficient, with upwards of 90% sampled points discarded (HEP Software Foundation, 2020). More powerful importance sampling algorithms would therefore be a welcome improvement (Buckley, 2020; WG et al., 2021).
First attempts to use machine learning (ML) to address this challenge explored using classical neural networks to sample (Bendavid, 2017; Klimek & Perelstein, 2020; Chen al., 2021) but typically suffer from excessive computational costs. Another avenue of research has been leverage gener- ative models successful in other fields such generative adversarial networks (Butter et al., 2019; Di Sipio et al., 2019; Butter et al., 2020; Ahdida et al., 2019; Hashemi et al., 2019; Carrazza & Dreyer, 2019). While such approaches do improve sampling speed by large factor, they have ma- jor limitations. In particular, they have no theoretical guarantees of providing correct answer on average (Matchev al., 2021) and poor control of uncertainties.
To avoid these disadvantages, our work exploits Neural Importance Sampling (NIS) (Müller et al., 2019 Zheng & Zwicker, 2019), which relies on normalizing flows and has strong theoretical guar- antees.
number of exploratory papers have been published on using NIS for LHC simulations (Gao al., 2020b; Bothmann et al., 2020; Gao et al., 2020a), as well as closely related variations (Bellagente et al., 2021; Stienen & Verheyen, 2021), but most studies have focused on preliminary investigation of performance without much concern for the practical usability of the method. Indeed, training requires function evaluations, which we are trying to minimize and data-efficiency training is there- fore important but under-appreciated concern. Furthermore. few authors have provided usable open source code, making the adoption of the technique in the HEP community difficult.
Our contribution to improve this situation can be summarized in two items:
The introduction of a new loss function and an associated new training algorithm for NIS. This permits the re-use of sampled points over multiple gradient descent steps, therefore making NIS much more data efficient.
The introduction of ZUNIS, a PyTorch-basee based library providing robust and usable NIS tools, usable by non-experts. It implements previously-developped ideas as well our new training procedure and is extensively documented.
BACKGROUND
2.1 IMPORTANCE SAMPLING AS AN OPTIMIZATION PROBLEM
Importance sampling relies on the interpretation of integrals as expectation values. Indeed, let us consider an integral over finite volume:
dx dxf(x), where V( / dx is finite.
(1)
Let be strictly positive probability distribution over Q, we can re-express our integral as the expectation random variable
/ dx f(x) = p(x)
N
f(X,)
I I=
(2)
N iml
p(X.)*
a (f.p))
whose mean is indeed 7 and whose standard deviation is deviation f(X)/p(X) for ~
where o(.p) is the standard
VN
o'(f.p) = f(())
12.
(3)
The problem statement of importance sampling is to find the probability distribution function p that minimizes the variance of our estimator for given N. In Neural Importance Sampling, rely
2

---

Normalizing Flows to approximate the optimal distribution, which we can optimize using stochastic gradient descent.
2.2 NORMALIZING FLOWS AND COUPLING CELLS
Normalizing flows (Tabak & Vanden-Eijnden, 2010; Tabak & Turner, 2013; Rippel & Adams, 2013; Rezende & Mohamed, 2015) provide way to generate complex probability distribution functions from simpler ones using parametric changes variables that can be learned to approximate a tar- get distribution. As such, normalizing flows are diffeomorphisms: invertible, (nearly-everywherc)) differentiable mappings with a differentiable inverse.
Indeed, if p(u), then T(u) (2x) where
(x= = p(u) [r('')',,
(4)
where Jr is the Jacobian determinant of T:
0T: Jr(u) det Ouj (u).
(5)
In the world of machine learning, T is called a normalizing flow and typically part of a parametric family of diffeoomopphimm (I) such that gradients VaJr are tractable.
Coupling cell mappings perfectly satisfy this requirement (Dinh et al., 2015; 2017; Müller et al., 2018): they neural-network- parrallntt bijections whose Jacobian factor can be obtained an- alytically without backpropagation or expensive determinant calculation. As such, they provide a good candidate for importance sampling as long they can be trained to learn an unnormalized target function, which is exactly what neural importance sampling proposes.
2.3 NEURAL IMPORTANCE SAMPLING
Neural importance sampling was introduced in the context of computer graphics (Müller al., 2018) and proposes to use normalizing flows as a family of probability distributions over which to solve the minimization problem of importance sampling.
/ dr P(x)
(())
(6)
p (x
Of course, to actually do so, one needs to find a way to explicitly evaluate L(0) and the original neural importance sampling approach proposes to approximate it using importance sampling. One needs to be careful that the gradient of the estimator the loss need not be the estimator of the gradient of the loss. The gradient of the loss can be expressed as
fx(x) Vo logp(x, VeL(O) = = - S dx (xx,0)
for which an estimator is proposed as
vec())- Σ
f(Xi) p (X.,0)
V log (X;,0), X
(8)
The authors also observed that other loss functions possible which share the same global mini- mum as the variance based loss: for example, the Kullback-Leibler divergence Dr between two functions is also minimized when they are equal. Such alternative loss functions not guaranteed to work for importance sampling, but they prove quite successful in practice. After training to mini- mize the loss estimator of eq. (8), the normalizing flows provides tractable probability distribution function from which to sample points and estimate the integral.
3 CONCEPTS AND ALGORITHMS
In this section we describe the original contributions of this paper. The major conceptual innova- tion we provide in ZONIS is more flexible and data-efficient way of training normalizing flows in

---

the context of importance sampling. This relies on more ngorous formulation of the connection between the theoretical expression of ideal loss functions in terms of integrals and their practi- cal realizations as random estimators than previous literature. We describe this improvement in section 3.1. We also give a high-level overview of the organization of the ZUNIS library, which implements this new training procedure.
3.1 EFFICIENT TRAINING FOR IMPORTANCE SAMPLING
In this section, we describe how we train probability distributions within ZUNIS using gradient- based optimizers. While the solution proposed the original formulation of NIS defined eq. (8) works, its main drawback is that it samples points from the same distribution that it tries to optimize. As result, new points X, must be sampled from p after every gradient step, which is very inefficient for slow integrands.
Our solution this problem is remember that the loss function is an integral, which can be evaluated by importance sampling using any PDF, not only We will therefore define an auxiliary probability distribution function ((z), independent from 0, from which we sample to estimate our loss function:
f(0) I dx p(z,0) a0a ()(())
f(x)2
This is the basis for the general method we use for training probability distributions within ZÜNIS, described in algorithm 2. Because the sampling distribution separated from the model train, the same point sample can be reused for multiple training steps, which is not possible when using eq. (8). This is particularly important for high-precision particle physics predictions that involve high-order perturbative calculations or complex detector simulations because function evaluations can be extremely costly. We show in section 4, in particular in fig. that reusing data indeed has very significant impact on data efficiency.
Algorithm 1: Backward sampling training in ZONIS
Data: Parametric PDF p(x,0 0)
Result: Trained PDF p(x.
for M steps do
Sample a batch
from g;
Compute the sampling PDF values ((<;);
Compute function values (xx);
Start tracking gradients with respect to 0;
for N steps
Compute the PDF values from the parametric PDF p(a;, 0);
f(x,)2
Estimate the loss f.
N 1=1 p(z,0)a(z))
Compute VoL using backpropagation;
Set θ t 0 VoL; Reset gradients;
10
end end
14 return p(x, 8):
After training, is discarded and the integral is estimated from the optimized
The only constraint model. Much like in other applications of neural networks, we have found that stochastic gradient descent can yield good results despite noisy loss estimates. We propose three schemes for q:
is that it yields good enough estimate that gradient steps improve the
A uniform distribution
4

---

A frozen copy of the model, which can be updated once in while' (survey stttsssts_ sattegy
forwardd")
An adaptive scheme starting with a uniform distribution and switching to sampling from frozen model when it is expected to yield a more accurate loss estimate variance").
An important point to notice is that the original NIS formulation in eq. (8) can be restated limiting case of our approach. Indeed, we take to frozen version of the model p(æ; 0o). which we update everytime sample points (setting N 1 in algorithm 2), the gradient update in line 9
f(x)² x--((r,00 p 0) ()p(2, 80) 10=8
f(x) Vo logp(r,0) 0).
V
E
(10)
xm p(r
3.2 THE ZÜNIS LIBRARY
On the practical side, ZONIS a PyTorch-based library which implements many ideas formulated in previous work but organizes them in the form of modular software library with casy-to-uss interface and well-defined building blocks. We believe this structure will help non-specialist use without understanding all the nuts and bolts, while experts easily add new features to responds to their needs. The ZUNIS library relies on three layers abstractions which steer the different aspects of using normalizing flows to learn probability distributions from un-normalized functions and compute their integrals:
F 1 ows, which implement bijective mapping which transforms points and computes the corresponding Jacobian are described appendix F.1
Trainers, which provide the infrastructure to perform training steps and sample from flow models are described in appendix F.2
Integrators which use trainers steer the training model and compute integrals are described in appendix F.3
4 RESULTS
In this section, we evaluate ZÜNIS on a variety of test functions to assess its performance and com- pare it to the commonly used VEGAS algorithm (Peter Lepage, 1978; Ohl, 1999). We first produce a few low dimensional examples for illustrative purposes, then move on to integrating paramet- ric functions in various dimensions and finally evaluate performance on particle scattering matrix elements.
4.1 LOW-DIMENSIONAL EXAMPLES
Let start by illustrating the effectiveness of ZONIS in low dimensional setting where readily visualize results. We define three functions on the 2D unit hypercube which each illustrate some failure mode of VEGAS (see appendix C).
We ran the ZUNIS Integrator with default arguments over ten repetitions for each function and report the performance of the trained integral estimator compared to flat estimator and VEGAS in table 1. Overall, ZÜNIS Integrators learn to sample from their target function extremely well: we outperform VEGAS by a factor 100 for the camel and the slashed circle functions and factor 30 for the sinusoidal function and VEGAS itself provides no advantage over uniform sampling for the latter two.
We further illustrate the performance of our trained models by comparing the target functions and density histogram for points sampled from the normalizing flows in fig. 1, which shows great quali- tative agreement.
"Tii is inspired by deep-Q. @ learning, where two copies of the value model are used: frozen one used to sample actions, and trainable one used to estimate values in the loss function. Here the frozen copy is used to sample points, and the trainable model is used compute PDF values used the loss function
5

---

Variance Reduction
Camel
Slashed Circle
Sinusoidal
uniform
1.8 0.4 10
8.9 0.9
10*
2.0 + 0.5 10a
VS. VEGAS
7 10a
8.8 + 0.9 x 10
1.6 + 0.5 10
Table I: Variance reduction (high is good) for the camel, slashed circle and sinusoidal functions compared uniform sampling and to VEGAS over 10 repetitions.
Target Function
Target Function
Target Function
0.2
0.2
0.2
0.4
0.4
Sampling Dersity Histagram
Sampling Dereity Histogram
Semp ing Density Histogrem
0.6
0.5
0.4
0.4
0.4
0.2
0.2
0.4
0.6
0.8
0.2
0.4
0.6
0.8
0.2
0.4
0.6
0.8
Sampled Points
Samp ed Points
Sampled Points
(a) Camel function
(b) Sinusoidal function
(c) Slashed circle function
Figure 1: Comparison between target functions and point sampling densities for la the camel func- tion, 1b the sinusoidal function, Ie the slashed circle function. Supplementary fig. 7 shows how points are mapped from latent to target space.
4.2 SYSTEMATIC BENCHIARRSS
Let us now take a more systematic approach to benchmarking ZÜNIS. We compare ZUNIS Inte- grators against uniform integration and VEGAS using the following metrics: integrand variance (a measure of convergence speed, see section 2. unweighting efficiency (a measure of the efficiency of exact sampling with rejection, see appendix A) and wall-clock training and sampling".
ZÜNIS improves convergence rate compared to VEGAS
this experiment, we focus on camel function defined in eq. (12) and scan 35 configurations spanning from 2 to 32 dimensions over function variances between
and 102 as shown in table
Except in the low variance limit, ZÜNIS reduce the required number of points sampled to attain a given precision on integral estimates without any parameter tuning, attaining speed-ups of up to 1000 both compared to uniform sampling and VEGAS-based importance sampling, shown in fig. 2a-2b and table 4. Unweighting efficiencies are also boosted significantly, although more mildly than variances, as shown fig. 2c-2d, which could attribute to PDF underestimation in regions with low point density; the nature of the veto algorithm makes it very sensitive to few bad behaving points in the whole dataset.
ZÜNIS is slower than VEGAS ZÜNIS does not, however, outclass VEGAS on all metrics by far: as shown in fig. 3, training is few hundred times slower than VEGAS and sampling is 10-50 times slower, all while ZÜNIS runs on GPUs. This be expected given the much increased com-
We provide details hardware in appendix

---

benchmark
Camel function benchmark
10F
Dim
Dim
Sperc-up vs. uniferm
L07
Speok Vs. VEGAS
L07
16
10
1.0
10"
10
2.5
10.0
12.5
2.5
7.5
10.0
12.5
Function Standard Deviation (eU//BUU)
Function Standard Devistion (rt///UU))
(b)
Camel Function banchmark
Camel function banchmark
Unweighting ERt. Robio va Uniterm
102
10a
Dim
Unweighting us
18
10
209
2.5
5.0
1.5
10.0
12.5
2.5
5.0
1.5
10.0
12.5
Function Standard Deviation (eU))EL))
Function Standard Devistion (oU)EEU))
(c)
(d)
Figure 2: Benchmarking ZONIS against uniform sampling and VEGAS with default settings. In (2a-2b), we show the sampling speed-up (ratio of integrand variance) as a function of the relative standard deviation of the integrand, while we show the unweighting speed-up (ratio of unweighting efficiencies) in (2c-2d).
putational complexity of normalizing flows compared the VEGAS algorithm. As such, ZONIS not a general replacement for VEGAS, but provides clear advantage for integrating time-intensive functions, where sampling is negligible overhead, such as precision gech--nergy--hhssie simula- tions.
Sampling Speed Comparison: ZuNIS VEGAS
Training Speed Comparison: ZuNis VEGAS
0.10
0.003
Samping spaed ratlo
0.004
0.06
Speed
16
0.003
0.04
0.002
0.001
0.00
10a
10" Batch Size
10"
0.000
8
Dim
(a)
(b)
Figure 3: Comparison of the training and sampling speed of ZUNIS and VEGAS. As can be ex- pected, ZUNIS is much slower than VEGAS, both for training and sampling, although larger batch sizes can better leverage the power of hardware accelerators.
The new loss function introduced in ZÜNIS improves data efficiency We have shown that ZÜNIS is very performant importance sampling and event generation tool and provides significant improvements over existing tools, while requiring little fine tuning from users. Another key result that the approach training introduced in section 3.1 has large positive impact

---

performance. Indeed, as can be seen in fig. 4, re-using samples for training over multiple epochs provides a 2- to 10-fold increase in convergence speed, making training much more data-efficient.
For this experiment, we use forward sampling, where the frozen model is used sample batch of points which are then used for training over multiple epochs before resampling from an update of the frozen model. As a result, we reproduce the original formulation of NIS in eq. (8) when we use single epoch as shown in eq. (10).
Camel Benchmark: ZONIS vs. Uniform
Camel Benchmark: 5 Epochs vs. 1 Epoch
Dim
0.8
Speed -pppaxxpxxxxxxxx
0.6
8
16
Spoat
0.4
32
6
2
0.0
0
10
20
30
40
50
4
8
16
Epochs
Dim
(a)
(b)
Figure 4: Figure 4a: Effect of repeatedly training on the same sample of points multiple epochs. For all settings, there is large improvement when going from one to moderate epoch counts, with a peak around 5-10. Larger number of epochs lead overfitting, which impacts performance nega- tively. Figure 4b: Comparison between optimal data reuse (5 epochs) and the original NIS algorithm epoch).
4.3 MADGRAPII CROSS SECTION INTEGRALS
Cross-sections are integrals of quantum transition matrix elements for a scattering process such as a LHC collision and express the probability that specific particles are produced. Matrix ele- ments themselves are un-normalized probability distributions for the configuration of the outgoing particles: therefore both valuable integrate them to know the overall frequency given scattering process, and sample from them to understand how particles will be spatially distributed as they fly off the collision point.
We study the performance of ZONIS in comparison to VEGAS by studying three simple processes at leading order in perturbation theory. e 14 + e via Z, dd - dd via Z and uc - ueg (with 3-jet the first process
cuts based on AR), table 2 and fig. 5. We two other, quark-initiated processes are used to illustrate specific points. Indeed, both feature narrow regions of their integration space with large enhancements, due respectively to Z-boson resonances and infra-red divergences.
very easy reference while the
www
www.
000000o
Z
0gg
d
(a)
(b)
(c)
Figure 5: Sample Feynman Diagrams for e f + e via Z, dd - dd via and uc - ueg
We evaluate the matrix elements for these three processes by using the FORTRAN standalone in- terface MADGRAPH5 LARAAAH55SMCCO AMC@NA (Alwall et al., 2014). The two hadronic processes are con- volved with parton-distribution functions from LHAPDF6 (Buckley et al., 2015). We parametrize
8

---

14 C via dd dd via Z uc -
ueg
dimensions
normalized standard deviation
1.45 x 10
6.57 x 10-2
0.96
Table 2: Comparison of the three test processes.
phase space (the particle configuration space) using the RAMBO diet algorithm (Plätzer, 2013) implemented for PyTorch in TORCHPS (Götz, 2021).
We report benchmark results in fig. in which trained 500,000 points for each process using near-default confi aafiguuaation, scanning only over variance and Kullback-Leiin losses.
Speed-up rilativ VECAS
Unweighting Efficinccy raBo relative VEGAS
107
10A
10a
10
10
Survey Strategy
Survey Strategy Maptiee Des Adaptive
Mapeive Dec
10
C o
CT of vi
di
sig
(b)
Figure 6: Average performance of ZONIS over 20 runs relative to VEGAS, measured by the relative speed-up in fig. and by the relative unweighting efficiency in fig. 6b.
As previously observed, little convergence acceleration achieved for low variance integrands like uu e μ., but unweighting still benefits from NIS. The two hadronic processes illustrate typical features for cross sections: training performance is variable and different processes are optimized by different loss function choices².
The performance dd + dd shows nice improvement with ZÜNIS while that of uc - ucg is more limited. This be understood by comparing to importance sampling (see appendix D.3): it is in fact VEGAS, which performs significantly better on uc - ueg compared to dd - dd because the parametrization of RAMBO is based splitting invariant masses, making them aligned with the enhancements in the ucg phase space and allowing great VEGAS performance. This drives a key conclusion for the potential role of ZUNIS in the HEP simulation landscape: not to replace VEGAS, but to fill in the gaps where it fails due to inadequate parametrizations, we illustrate here by using non- non-multichanneled dd - dd as proxy for more complex processes.
5 CONCLUSION
We have showed that ZÜNIS outperform VEGAS both in terms of integral convergence rate and unweighting efficiency on specific cases, at the cost of a significant increase in training and sampling time, which is acceptable tradeoff for high-precision HEP computations with high costs. In this context, the introduction of efficient training power of neural importance sampling where function evaluation costs further testing is required to ascertain how far NIS can fill the gaps left by VEGAS for integrating complex functions, there is already good evidence that ZUNIS can provide needed improvements in specific cases. We hope that the publication of usable toolbox for NIS such as ZÜNIS will stir wider audience within the HEP community to apply the method so that the exact boundaries its applicability can be more clearly ascertained.
key element to making the most of the
major concern. While
GGenrallly, smoother functions are better optimized with the Kullback-Leibler loss while functions with peaks benefit from using the variance loss. As show in appendix D.4, choosing an adaptive strategy is generally advisable whatever the loss

---

REPRODUCIBILITY STATEMENT
An anonymized version of the library is available Github.
The recommended procedure to reproduce the experiments to clone the repository and install the requirements. txt.
Python requirements using pip install
The data to reproduce the experiments can be generated using scripts provided in the repository experi innersmmestss... /bwnchharkss......ee.. in which Jupyter notebooks are also available reproduce the figures of the paper. The following scripts are available:
Badddaattt./////mmmm//aaaaaaaaat to generate camel integration data
btadhhaak...///mmmmmmmmmooaaaaaa to generate camel sampling speed data
benoobaakk_mdggggpp////ebeeaaann generate e - 6 a via integra- tion data
benohhartL_magggggp///eeenennnrd to generate dd dd via Zintegration data
sh to generate we + weg integration data
These scripts assume that CUDA GPUs available and run benchmarks in parallel. If fewer GPUs are available, it is recommended to modify the scripts to run the benchmarking scripts se- quentially (by removing the ampersand) and adapt the --cuda=N option.
REFERENCES
C. Ahdida, R. Albanese, Alexandrov, A. Anokhina, S. Aoki, G. Arduini, E. Atkin, N. Azorskiy, J.J. Back, A. Bagulya, and et al. Fast simulation of muons produced the ship experiment using generative adversarial networks. Journal of Instrumentation, (4(11):P1100 1012)))P11 1028, Nov 2019. ISSN 1748-0221. doi: 1088/1748- 01177748-222//1 1wt.//////////1///02200000011110
URL http: //tt:://dxxoot org/10.
Alwall, R. Frederix, S. Frixione, V. Hirschi, F. Maltoni, Mattelaer, H.-S. Shao, T. Stelzer, P. Torrielli, and M. Zaro. The automated computation of tree-level and next-to-leading order differentiaal cross sections, and their matching to parton shower simulations. Journal of High Energy Physics, 2014(7), Jul 2014. ISSN 1029-8479. doi: 10.1007/hhep077004)07999 URL
httpp:/dde doi. orgoorrg/// .. JHEP (2014) 079.
The ATLAS Collaboration. ATLAS HL-LHC Computing Conceptual Design Report. 2020.
Frederik Beaujean and Allen Caldwell. Initializing adaptive importance sampling with markov chains, arXiv preprint arXiv:1304.7800 2013.
Marco Bellagente, Manuel Haußmann, Michel Luchmann, and Tilman Plehn. Understanding Event- Generation Networks via Uncertainties. wxXiv:2104004543 [hep-ph], April 2021.
Johannes Bellm et al. Herwig 7.0/Herwig++ 3.0 release note. Eur. Phys. J. C, 76(4):196, 2016. doi: 10.11400epps3100222014441118
Joshua Bendavid. Efficient Monte Carlo Integration Using Boosted Decision Trees and Generative Deep Neural Networks. 2017.
Enrico Bothmann, Timo Janßen, Max Knobbe, Tobias Schmale, and Steffen Schumann. Exploring phase space with Neural Importance Sampling. SciPost Physics, 8(4):069, April 2020. ISSN 2542-4653. doi: 10.2 4ww.//SiiossPys.8....66.
Enrico Bothmann et al. Event Generation with Sherpa 2.2. SciPost Phys., 7(3):034, 2019. doi: 10.21 26144///1oss111yy77..33333
10

---

Johann Brehmer and Kyle Cranmer. Flows for simultaneous manifold learning and density estima- tion. 2020.
Andy Buckley. Computational challenges for MC event generation. Journal of Physics: Conference Series, 1525:012023, April 2020. ISSN 1742-6596. doi: 10.1088/1742-65900 1255////22//////222.
Andy Buckley, James Ferrando, Stephen Lloyd, Karl Nordström, Ben Page, Martin Rüfe- nacht, Marek Schönherr, and Graeme Watt. LHAPDF6: parton density access in the LHC precision era. The European Physical Journal C. 75(3), Mar 2015. ISSN 1434-6052. doi: 10. s s10052- 015 3318-8.
URL http:/ /wtp/dxx.ooooooooooo/ooooooo/ooo www.ooorg/o//00 111100.11//ppppo
Monica F. Bugallo, Victor Elvira, Luca Martino, David Luengo, Joaquin Miguez, and Petar M. Djurie, Adaptive importance sampling: The past, the present, the future. IEEE Signal Pro- cessing Magazine, 34(4):60-79, 2017. doi: 10.1 (0.109///SS.221122299922
Mónica Bugallo, Luca Martino, and Jukka Corander. Adaptive importance sampling in signal processing. Digital Signal Processing, 47:36-49, 49, 2015. ISSN 1051-2004. doi: https://doi. orgp//0..0166f.ddp,2055........0 URL https: / /www. sc encedi rect. com scc./cc... com/scienee////// 864. Special Issue in Honour of William J. (Bill) Fitzger-
wtt..cllepp../S155100044.. ald.
Anja Butter, Tilman Plehn, and Ramon Winterhalder. How to Ihe events. SciPost Physics, (6), Dec 2019. ISSN 2542-4653. doi: 10.21 10.21 440/sccipostphys....0000 URL http: doi. org/10. 211446//ciiostPhyya 7.6 7.70.7
Anja Butter, Tilman Plehn, and Ramon Winterhalder. How gan event subtraction. SciPost Physies Core, 3(2), Nov 2020. ISSN 2666-9366 doi: 10211468/scppsiphyycco..33 3.000..... URL http: / //dxx.oof org/ C 214 68/ wwws/sssttsssoor......o http://...1..//eesssaspeeeeeeeee
009.
Olivier Cappé, Amaud Guillin, Jean-Miche Marin, and Christian P Robert. Population monte carlo. Journal of Computational and Graphical Statisties, 13(4):907-929 2004.
Olivier Cappé, Randal Douc, Arnaud Guillin, Jean-Michel Marin, and Christian P Robert. Adaptive importance sampling in general mixture classes. Statistics and Computing, 18(4):447-45) 2008.
Stefano Carrazza and Frédéric A. Dreyer. Lund jet images from generative and cycle-consistent adversarial networks. The European Physical Journal C, 79(11), Nov 2019. ISSN 1434- 6052. doi: 10. 11144//pp//100050011999797 sp:c/s1c052-011--
URL http:/ /dxx: doi. 11/111/1...111
I.-Kai Chen, Matthew D. Klimek, and Maxim Perelstein. Improved Neural Network Monte Carlo Simulation. SciPost Physics, 10(1):023, January 2021. ISSN 2542-4653. doi: 10.21468/
SciPostPhys.10.1.023.
Jean-Marie Cornuet, Jean-Michel Marin, Antonietta Mira, and Christian Robert. Adaptive multiple importance sampling. Scandinavian Journal of Statistics, 39(4):798-81) 2012. doi: https://doi. org/10.1111f..161 1467 9464447779999649 2011.00776..... URL https: /tt//onniilaarr.../////// wwwww. com/ doi/ abl//b///tt 10. 1111111111111111 1467 - 9469. 2011 00756.
Riccardo Di Sipio, Michele Faucci Giannelli, Sana Ketabehi Haghighat, and Serena Palazzo. Dijetgan: a generative-adveesariall network approach for the simulation of qed dijet events at the lhe. Journal of High Energy Physics, 2019(8), Aug 2019. ISSN 1029-8479. doi: 10. 0011007/hhe008(211)1111 URL http: //xp:.///x////// doi. org/10....... 1007 JHEP (2019) 110.
Laurent Dinh, David Krueger, and Yoshua Bengio. NICE: Non-linear Independent Components Estimation, 2015.
Laurent Dinh, Jascha Sohl-Dickstein, and Samy Bengio. Density estimation using Real NVP, 2017.
Randal Douc, Arnaud Guillin, J-M Marin, and Christian Robert. Minimum variance importance sampling via population monte carlo. ESAIM: Probability and Statistics, 11:427-447, 2007.
11

---

Integration and Victor Elvira, Luca Martino, David Luengo, and Mónica F. Bugallo. Improving population monte carlo: Alternative weighting and resampling schemes. Signal Processing. 131:77-91, 2017. ISSN 0165-1684.
doi: http:///do..or//1001111........
URL. https:
/ h/www. www.eeiieeciie.rrcccoo http:///.....ccao.ocaoa/aa//aa/a
Christina Gao, Stefan Höche, Joshua Isaacson, Claudius Krause, and Holger Schulz. Event genera- tion with normalizing flows. Physical Review D. 101(7):07002,, April 2020a. ISSN 2470-0010, 2470-0029. doi: 10.1 10.110//hhy.Ree..00.....
Christina Gao, Joshua Isaacson, and Claudius Krause. I-flow:
Sampling with Normalizing Flows. Machine Learning: Science and Technology, 144:0045023, November 2020b. ISSN 2632-2153. doi: 11-008/266222255/abbbbb
Niklas Götz. NGoet//TorchhS:.v11..01 March 2021. URL
zenodo... 1eroodo.463300 4639109. Available at hug.tubb..tbbtuuo/oo con///ootta/Toocc/////ez v1.0.1, version 1.0.1.
Bobak Hashemi, Nick Amin, Kaustuv Datta, Dominick Olivito, and Maurizio Pierini. Lhe analysis- specific datasets with generative adversarial networks, 2019.
HEP Software Foundation. HL-LHC Computing Review: Common Tools and Community Soft- ware. arXiv:2008. 13636 [hep-ex, physies:pyysics) August 2020. doi: 10.5281/zemodo..ooo 40090 14.
Yukito Iba. Population based Monte Carlo algorithms. Trans. Jap. Soc. Artif. Intell., 16:279, 2001.
doi: 10. 14.1577/jsaai.6622799
F. James. Monte Carlo Theory and Practice. Rept. Prog. Phys., 43:1145, 1980. doi: 10.1088/ 0034-4885/43/9/002.
S.P. Jones. Higgs Boson Pair Production: Monte Carlo Generator Interface and Parton Shower, Acta Physica Polonica B Proceedings Supplement, 11(2):295, 2018. ISSN 1899-2358, 2082-7865. doi:
J0.500//APhysPooBSSpp. 11.295.
R. Kleiss and R. Pittau. Weight optimization in multichannel Monte Carlo. Computer Physics Com- munications, 53(2-3):141-146 December 1994. ISSN 00104655. doi: 10.101600010-4655094)) 90043-4.
Kleiss, W. J Stirling, and S. Ellis. A new Monte Carlo treatment of multiparticle phase space at high energies. Computer Physics Communications, (0(2):359-373 June 1986. ISSN 0010-4655. doi: 10.10160010-445556699011 19-0.
Matthew Klimek and Maxim Perelstein. Neural basedd approach to phase space integration. SciPost Physics, 9(4):053, October 2020. ISSN 2542-4653. doi: (0.21466/SciPootthyss.... 0,500.5
Eugenia Koblents and Joaquín Míguez. A population monte carlo scheme with transformed weights and its application to stochastic kinetic models. Statistics and Computing, 25(2):407-425, 2015,
G. Peter Lepage. VEGAS: AN ADAPTIVE DOAATIVEEMUTTDDDMEENINNON INTEGRATION PROGRAM. March 1980.
G. Peter Lepage. Adaptive Multidimensional Integration: VEGAS Enhanced. Journal of Computa- tional Physics, 439: 10386, August 2021. ISSN 00219991. doi: 20.101.....p.2222.1.1 10386.
Konstantin T. Matchev, Alexander Roman, and Prasanth Shyamsundar. Uncertainties associated with GAN-generated datasets in high energy physics. arXiv:2002.06307 [hep-ex, physics:hep- ph, physics:physics), June 2021.
Thomas Müller, Brian MiWilliams, Fabrice Rousselle, Markus Gross, and Jan Novák. Neural Im- portance Sampling, 2018
Thomas Müller, Brian Mcwilliams, Fabrice Rousselle, Markus Gross, and Jan Novák. Neural Impor- tance Sampling. ACM Transactions on Graphics, 38-5):1-19, November 2019. ISSN 0730-0301, 1557-7368. doi: 10. 145/3341156
12

---

Thorsten Ohl. Vegas revisited: Adaptive Monte Carlo integration beyond factorization. Comput.
Phys. Commun., 120:13-19, 1999. doi: 10.1016/500110445559990000000000
Peter Lepage. new algorithm for adaptive multidimensional integration. Journal of
Computational Physics, 17(2):12-2033, 1978. ISSN 0021-9991. doi: https://oi.org//
1016/0021-9991 (7990002199991189000444 URL https://wwet http://wwcciiiccccccc http://wwwwwiiinnieeeee......... com/ coom/science
Simon Plätzer. RAMBO on diet, 2013.
Danilo Jimenez Rezende and Shakir Mohamed. Variational inference with normalizing flows. 32nd International Conference on Machine Learning, ICML 2015, 2:1530-1538, 2015. ISSN
9781510810587.
Oren Rippel and Ryan Prescott Adams, High-dimensional probability estimation with deep density models. arXiv preprint urXiv: 2013.
Bob Stienen and Rob Verheyen. Phase Space Sampling and Inference from Weighted Events with Autoregressive Flows. SciPost Physics, 10(2):038, February 2021. ISSN 2542-4653, doi: 10.
21468/SciPoxtPhyy.02...33
Esteban G Tabak and Cristina V Turner. A family of nonparametric density estimation algorithms. Communications on Pure and Applied Mathematics, 66(2):145-164, 2013.
Esteban G Tabak and Eric Vanden-Eijndes Density estimation by dual ascent of the log-like logliikl Communications in Mathematical Sciences, 8(1):217-233, 2010.
The HSF Physics Event Generator WG, Andrea Valassi, Efe Yazgan, Josh McFayden, Simone Amoroso, Joshua Bendavid, Andy Buckley, Mattco Cacciari, Taylor Childers, Vitaliano Ciulli, Rikkert Frederix, Stefano Frixione, Francesco Giuli, Alexander Grohsjean, Christian Gütschow, Stefan Höche, Walter Hopkins, Philip Iiten, Dmitri Konstantinov, Frank Krauss, Qiang Li, Leif Lönnblad, Fabio Maltoni, Michelangelo Mangano. Zach Marshall, Olivier Mattelaer, Javier Fer- nandez Menendez, Stephen Mrenna, Servesh Muralidharan, Tobias Neumann, Simon Plätzer, Stefan Prestel, Stefan Roiser, Marek Schönherr, Holger Schulz, Markus Schulz, Elizabeth Sexton- Kennedy, Frank Siegert, Andrzej Siódmok, and Graeme A. Stewart. Challenges in Monte Carlo event generator software for Luginosttm LHC. Computing and Software for Big Science, (1):12, December 2021. ISSN 2510-2036, 2510-2044. doi: 10. 10.10077/417781-02--000
Quan Zheng and Matthias Zwicker. Learning to Importance Sample in Primary Sample Space. Computer Graphics Forum, 38(2): 179, 2019. ISSN 1467-8659, doi: 10. 10.1111/cgf.1662..
A EVENT GENERATION WITH THE VETO ALGORITHM
Generating i.id points following an arbitrary probability distributions in high dimensions is not priori a trivial task. straightforward way to obtain such data use rejection sampling, which can be based on any distribution q from which we can sample. Given an i.i.d sample 21....,IN %. we can define weights wxx;)) = p(x;)//(t /q(x)) and keep/reject points with probability w(x:)) Wmmx.
The main metric for evaluating the performance of this algorithm is the unweighting efficiency: how much data kept from an original sample of size N on average, which is expressed as
Cynw
E 2r0q W maxx
(11)
B BACKWARD SAMPLING ALGORITHM
13

---

Algorithm 2: Backward sampling training in ZONIS
Data: Parametric PDF p(x, 0)
Result: Trained PDF p(x,
for M steps do
Sample batch 21:--..: Tuaa from d:
Compute the sampling PDF values ((r;);
Compute function values f (xz);;
Start tracking gradients with respect to 0;
for N steps do
Compute the PDF values from the parametric PDF p(a;, 0);
f(1;)²
Estimate the loss i.
p(x,0)q(z)
Compute VoL using backpropagation;
Set 0 + 0 VVL; Reset gradients;
10
end end
return p(x, 8):
C FUNDAMENTAL LIMITATIONS OF THE VEGAS ALGORITHM
We define three functions over the two dimensional hypercube:
x 112) 2
(1)
Seamet (x:) exp
exp
(12)
(-(1a))
[10 Oz )
2)
2
fo(x) min 1, exp
+ exp
(13)
fsi(x)) cos (A: xx²
(14)
which we refer respectively as the camel, slashed circle and sinusoidal target functions. We set their parameters follows
(0.55)
0.75 0.75
k - (1)) 0 0.1, r 0.3 (15)
f1 =
0.255
#22 =
We chose these functions because they illustrate different failure modes of the VEGAS algorithm:
Because of the factorized PDF of VEGAS, the camel function leads to 'phantom peaks" in the off diagonal. This problem grows exponentially with the number of dimensions but can be addressed by change of variable to align the peaks with an axis of integration.
The sinuoidal function makes it nearly impossible for VEGAS to provide any improvement: the marginal PDF over each variable of integration nearly constant. Again this type of issue can be addressed by change of variable provided one knows how to perform it.
The slashed circle function is an example a hard-for- VEGAS function that cannot be improved significantly by change of variables. One can instead use multi-channeling but this requires lot of prior knowledge and has computational costs since each channel is its integral.
D SUPPLEMENTARY RESULTS
D. QUALITATIVE EXAMPLES
14

---

Latent Space
0.6
Torget Space
3.0
Target Space
Target Space
0.5
0.6
0.2
0.0
0.2
0.8
0.6
0.8
1.0
Figure 7: Mapping between the uniform point density in the latent space and the target distribution for the camel function, the sinusoidal function, the slashed circle function. Points are colored based on their position in latent space.
D.2 SYSTEMATIC BENCHMARK DETAILS
Dimension
2
8
16
32
did
rel. st.d
2.00 X 10-2
8.93 x 10
2.04 X 10
3.61
10
5.06
10-1
0.001
14.1
6.32 x 10-2
1.64 x 10-1
3.14 x 10-1
4.64 x 10-1
6.06 x 10-1
0.01
4.41
Camel param.
2.21 x 10-1
3.78 10-1
5.23 10-1
6.67
10-1
8.25
10
0.1
1.22
4.51 X 10-1
5.93 x 10-1
7.43 x 10-1
9.11
10-1
1.10
0.3
0.56
6.44 x 10-1
7.99 x 10-1
9.74 x 10-1
1.18
1.41
0.5
0.32
8.62 x 10
1.05
1.26
1.51
1.81
0.7
0.19
1.21
1.45
1.73
2.07
2.47
1.0
0.10
Table 3: Setup of the 35 different camel functions considered to benchmark ZONIS. We scan over function relative standard deviations, which correspond to different parameters for each dimen- sion(eq. (12)). We provide the corresponding width of ID gaussian (σia) with the same variance for reference.
15

---

Dimension
Tid
rel. st.d
2
4
8
16
32
7.6
x 10
6.0 X 10th
4.3
x 10a
9.0
x 10
3.5 x 10
0.001
14.11
2.4
VEGAS speed-up
2.0
10
5.8
10
5.1
10
2.8103
10
4.5
10
0.01
4.41
3.4 x 10
1.0 X 10
3.5
x 10a
1.9201
X 10
1.2 X 10a
0.1
1.22
6.3 10
2.0
10'
7.2 100
2.8+0.1
10
1.3 100
0.3
0.56
1.2 x 10
3.5
1
9.1 x 10-1 3.7
x 10-1
2.0 X 10-1
0.5
0.33
2.4 0.2 x 10
5.3+0.5
10
1.3
x 10 5.5
10
3.0
X 10
0.7
0.20
0.2
4.0 0.5
10
7.8+06
10
1.6
10 7.8
10
4.8 X 10
1.0
0.11
Table 4: Variance reduction factor compared VEGAS for each of the 35 different camel setups defined in table 3.
D.3 COMPARING ZONIS WITH UNIFORM SAMPLING MATRIX ELEMENTS
Speed-up relative Uniforn Sampling
Unweighting Efficiency ratio relabve Uniform Sampling
10P
Survey Strategy
Adaptive Dez Mapptive
10a
Spod--p
18
10"1
Survey Strategy
Mapeive Dec
10
C -e0
CT -0
di
sig
Figure 8: Average performance of ZÜNIS over 20 runs relative to flat sampling, measured by the relative speed-up 8a and by the relative unweighting efficiency in 8b.
D.4 EFFECT OF SURVEY STRATEGIES
In the following, want investigate how the integration of the three example processes in 4.3 with ZUNIS behaves relative to VEGAS dependence of the choice of the loss function, the survey strategy and the number of epochs during training. For all other options, the default values are chosen again, except for the number of points during the survey phase, which is set to 500,000.
Figure 9a shows that for simple process like J e I via Z where no correlations exist, Z0- NIS camnot reach the speed-up achieved by VEGAS. Variance loss seems to lead to higher variance improvements than DE loss. Contrary, 9b shows that ZÜNIS can greatly improve the unweighting efficiency for this process. The effect is again consistently stronger when using variance loss. Using a flat survey strategy suffers for both loss functions from overfitting, whereas adaptive sampling in average performs slightly better and does not show overfitting.
dd + dd via Z presents more realistic use case, as the parton distribution functions introduce correlations between the integration variables which present challenge to the VEGAS algorithm.
For this process, both the speed-up and the unweighting efficiency ratio clearly favor the variance loss again, which outperforms in both metrics the DxL loss when multiple epochs are used, as can be seen in 10. The richer structure the integrand reduces the effect overfitting. Therefore, the performance increases stays approximately constant except for the combination of D loss and the forward survey strategy. For the variance loss, it becomes apparent that the flat survey strategy increases much slower in performance than alternative strategies function of the number of epochs. However, for combination of losses and survey strategies, VEGAS could outperformed substantially in terms of integration speed.
16

---

M"" -> Jr Speed-up relative to VEGAS
Number Epocht
20
-0.70
S
6.11
6M
67X
-0.85
Forword (flt init)
Forward
Flat
Adaptive
Foowwrrd (flat init)
Forword
Flat
Adaptive
Variance
Dau Loss
I pow and Survey Strategy
Unweighting Efficiency Ratio VEGAS
Number of Epochs
10 20
5.00
8.5
51
Forward (fla init)
Adaptive
Forward (fla: init)
Adaptivs
Variance
Dua.
Loss and Survey Strategy
(h)
Figure 9: Median of the performance of ZONIS over 20 runs relative to VEGAS for the process eT e J via Z depending the loss function, the survey strategy and number of epochs, measured by the relative speed-up in 9a and by the relative unweighting efficiency in 9b.
dd Al Speed up relative VEGAS
Number Epochs
20
144422
LiKe
12
100
10
5
15
-00
Forward (flat
Forward
Flat
Adaptive
Forward
Forward
Flat
Adaptive
Variance Loss
Dgg Loss
and Survey Strategy
(a)
dd - di Unweighting Efficiency Ratio relative to VEGAS
Number Lpochs
20
011
10
5.0
5
36a
144
0%5
2.5
Forward (flt init)
Forward
Flat
Adaptive
Forward (flas init)
Forward
Flat
Adaptive
Variance Loss
Pigz Do
Loss and Survey Strategy
(b)
Figure 10: Median of the performance of ZUNIS over 20 runs relative to VEGAS for the process dd + dd via Z depending on the loss function, the survey strategy and number epochs, measured by the relative speed-up in 10a and by the relative unweighting efficiency in 10b.
An opposite picture is drawn by the process uc - ucg in figure 11, for which, apart from the flat survey strategy, DEL loss is in general favored both for integration speed well unweighting efficiency ratio. The adaptive survey strategy is here giving the best results, although for a high number of epochs causes overfitting for the unweighting efficiency.
17

---

MC + ucy Speed-up relative to VEGAS
Numbbe Cpochs
20
6
10
4
5
6:5
CAS
632
52
Forward init)
Forvard
Flat
Adaptive
Forward init)
Forvard
Flat
Adaptive
Variance av
Fic Loss
Loss and Survey Strategy
(a)
mog Unweighting Efficiency Ratio relative to VEGAS
Number Loochs
20
1,004
3.477
L1
30
13
5
3000
Ca61
1.2
wits
Forward (flat init)
Adaptive
Forward (flat init]
Adaptive
Variance
Pys
Loss and Survey Strategy
(b)
Figure 11: Median the performance of ZÜNIS over 20 runs relative to VEGAS for the process we ueg depending the loss function, the survey strategy and number of epochs, measured by the relative speed-up in laa and by the relative unweighting efficiency in 1 lb.
The take-home message of this section is, one the one hand, that the flat survey strategy is in general not recommended. Apart from this, the most important mean to improve the quality of importance sampling are testing whether, independent of the survey strategy, the loss function should be chosen differently.
E EXACT MINIMZATION OF THE NEURAL IMPORTANCE SAMPLING ESTIMATOR VARIANCE
Let us show that the optimal probability distribution for importance sampling is the function itself. Explicitly, as discussed in section 2.1, want to find the probability distribution p defined over some domain Ω which minimizes the variance of the estimator f ((X)/p(X) for X p(X). We showed that this amounts to solving the following optimization problem:
min C(p) = / dx fxx)²
such that / dx p(x) = 1,
(16)
p(z)*
which we can encode using Lagrange multipliers
= arg min C(p, = dx f(x)2
+
p((x) V(
(17)
p(x)
We can now solve this problem by finding extrema with functional derivatives
5L(p,>) - (
f(x)2
(18)
8p(x)
p(x)2'
which indeed zero if p(x) & [f(z)) Furthermore, this extremum is certainly a minimum because the loss function is positive and unbounded. Indeed, if we separate Q into two disjoint measurable subdomains 22 and Q2 and define Pa(x)) such that points are drawn uniformly over 22 with prob- ability and uniformly over Gh with probability 1 a, then the resulting loss function would be
((pa) = V(2+) / dr f(x)²
v( / da f(x)²
(19)
1 - Ja,
which can be made arbitrarily large by sending to O.
18

---

F HIGH-LEVEL CONCEPTS OF THE ZUNIS API
F.1 NORMALIZING FLOWS WITH FLOW CLASSES
Flows map batches of points and their densities. 'he ZONIS library implements normalizing flows by specifying general interface defined as a Python abstract class: GeneralFlow. All flow models ZONIS are child classes GeneralF Flow which itself inherits from Pytorch Modall interface.
As defined in section 2.2, a normalizing flow in ZÜNIS is a bijective mapping between two dimensional spaces, which in practice always the unit hypereube [0, 1] or RR with tractable Jacobian so that can be used to map probability distributions. To this end, the Genera ra 1F1 Flow interface defines normalizing flows callable Python object which acts batches of point drawn from a known PDF p. A batch of points 2; with their PDF values is encoded as Pytorch Tensor object X organized follows
x= (X1,.... bbach E R batch Rd+1
(20)
where each X; corresponds to points stacked with its negative log density
Fi,1
X,
(21)
Fi,d log p(x,))
Encoding point densities by their negative logarithm makes their transformation under normalizing flows additive. Indeed if we have mapping Q with Jacobian determinant jQ. then x ~ p(z) is mapped to = (x) ~ p(y) such that
-log((y) logp(x)) logjq(z))
(22)
Coupling Cells are flows defined by an element-wise transform and a mask. All flow mod- els used in ZÜNIS in practice are implemented as sequence of coupling cell transformations acting on subset of the variables. The abstract class GeneralCoupl heppinggell and its child In- vertibles eCoupl vertibleeoupl ngCell specifies the general organization of coupling cells as needing to be instantiated with
dimension
mask defined list of boolean specifying which coordinates transformed or not a transform that implements the mapping of the non-masked variables
In release v1.0 of ZUNIS two such classes are provided: PWLinearCoupl ing and PlQuadraticCo cenddraatccopppingg which respectively implement the piecewise linear and piecewise quadratic coupling cells proposed in (Müller al., 2018). New coupling cells can easily be im- plemented, as explained in appendix F.4. Both existing classes rely on dense neural networks for the prediction of the shape of their dimeensiiional piecewise polyynomia mappings, whose parameters are set at instantiation.
Here how one can use piecewise-linear coupling cell for sampling points
import torch
from zunis wtt.ss.odess... www.s. wttp://wwww.in.n..ccc.. piecewise__padratic. import iinpaarratccoupllin
leslinng_elllss cells. piecewise htteecewwsse_coppllnn.._. pew.we.eeeeewwssseeeeee........e coupling.
d-2
N_batch=10
mask = False,
x = torch.zeros (N_ mu_batch,ede d+1))
Sample the d first entries uniformly, keep 0. for the negative log jacobian
x al....:-1).unifomm..
uniform
(xint # 10.3377, 0.4362, 0.] f = FRQuadraticCoupl. (d=d, mask=mask)
f(x)
print (y[o)) [0.3377, 0.4411, -0.0314]
19

---

We provide further details of the use and possible parameters of flows in the documentation of ZONIS: https::t///// http://uniisnooonnnoaaaaao a-an g thub. 1 / /uzii/n.......////....////
F.2
TRAINING WITH THE TRAINER CLASS
The design of the ZONIS library intentionally restricts Flow models to being bijective mappings instead of being ways to evaluate and sample from PDFs
not restrict their applicability (see Brehmer Cranmer (2020) for an example). The specific application in which one uses a how precisely
normalizing flow, and in such a model trained. As a result, ZUNIS bundles together the low-level training tools for F ow models together with sampling tools inside the Trainer classes.
samples from is intimately linked how
The general blueprint for such classes is defined in the the main implementation for users is provided as StatefulTraines At instantiation, all trainers expect Flow model and flow_pr_ These elements together define a probability distribution over the target space from which one can sample.
abstract class while
which samples point from fixed PDF in latent space.
There are two main ways one interacts with Trainers:
One can sample points from the PDF defined by the model and the prior using the sam- ple_forward method.
One train a pre-sample batch of points, their sampling PDF and the corresponding function values using train_on_batch (self, px, fx)
In practice, we expect that the main way users will use Trainers is for sampling pre-trained models. In the context of particle physics simulations for example, unweighting is common task, which aims sampling exactly from known function f. common approach the Hit-or-miss algorithm (James, 1980), whose efficiency is improved by sampling from PDF approaching f- This is how one would use trainer trained on f::
[...] import or train trainer exeetraine_ttraiee
import torch
Sampling points
xj
x - x[:, :-1]
px (-x(:,-1]).ee
fx (xx)
Applying the veto algorithm fmax = fx.max
veto (fx/fmax torch. zeros (fx) (0 1.)) 0.
unwweighted x[veto]
x unweightedd follows the PDF obtained by normalizing f.
F.3 INTEGRATING WITH THE INTEGRATOR CLASS
Integrators are intended as the main way for standard users to interact with ZUNIS. They provide high-level interface the functionalities of the library and only optionally require users to know to what lower levels of abstractions really entail and to what their options correspond. From practical point of view, the main interface of ZUNIS for integration implemented as the Integrator: which is factory function that instantiates the appropriate integrator class based on choice of options.
All integrators follow the same pattern, defined in the SurveyyeefieeInneggaaaaai BaseIntegrator abstract classes. Integration start by performing survey phase, in which optimizes the way it samples points and then refine phase, in which it computes the integral by using its learned sampler. Each phase proceeds through a number steps, which can be set at instantiation when integrating:
and
20

---

Setting values at instantiation time
integrator Integrator I=f, deiter_surveyye n __iter_rrre retirrr)ff # Override at integration time
entegral_data
interrator. (n_survey=10, n_refine=10)
For both the survey and the refine phases, using multiple steps is useful monitor the stability of the training and of the integration process: if one step is not within a few standard deviations of the next, either the sampling statistics are too low, or something is wrong. For the refine stage, this is the main real advantage of using multiple steps. On the other hand, at each new survey step, new batch of points re-sampled, which can be useful to mitigate overfitting.
By default, only the integral estimates obtained during the refine stage are combined to compute the final integral estimate, and their combination is performed by taking their average. Indeed, because the model is trained during the survey step, the points sampled during the refine stage correlated an uncontrolled way with the points used during training. Ignoring the survey stage makes all estimates used in the combination independent random variables, which permits to build formally correct estimator of the variance of the final result.
F.4 IMPLEMENTING NEW COUPLING CELLS
To implement a new invertible coupling cell inheriting from Liverriblecougpinnggee.. one must provide an InvertibleTee form object and define callable attribute T computing the parameters of the transform. For example, consider a very simple linear coupling cell over
yA= = x y = exp (T(xa) x B
== Qxx):
(23)
where T(*A) is scalar strictly positive value. This can be defined in the following way in ZONIS
import torch
from zunis, wwn.is.mmeels..m flows.
/ow. ott.://w..coopooo.oo.oooogeeaaaa www.geeeenanopplliiie import
Invertiblecouplingee..
from zunis cumis.models..c flows. cowpling cells. transforms import
Invertibleerannsfom
from ounis.models. aun...noeell..lyyrrr.r.raaaaaa import AiiitaaryyhhpeRRcc.gggllaaaa
alpha torch.exp(1) class
def forward
logj shx.ehape [-1] return x*alpha, logj det back
alpha torch.exp logj -Tax.shape return x*alpha, logj
21

---

class LinearCouplingCell (LwvertiblecouppLnggeee.
def init (self, d, mask, nn width, nn :
transferm - LinearTransform/
super Lineerppplll self) init_ (d=d, mask-nask, transform= transform)
i_in sun (mask)
self. = ArhtrarrshapeeettggglaaaMMM (d (d_indd__ in,
out out_shape- (1,), d_hiddea=an_wide n n_hhidden-n__eepth
G HARDWARE SETUP DETAILS
The computations presented in this work were performed on computing cluster using Intel(R) Xeon(R) Gold 5218 CPU @ 2.30GHz with 376 GB RAM. Processes which could be performed on the GPU were done on GeForce RTX 2080 having 12 GB memory and running on CUDA 11.0.
22


Reviews:
[{"review_id": "aUkzj2ydz9y", "reviewer": "Reviewer_dDs7", "summary": "", "questions": "", "limitations": "", "rating": 3, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "The work contains a new importance sampling scheme. The state-of-the-art discussion  is quite poor, ignoring all the recent literature on importance sampling and adaptive importance sampling schemes (except the VEGAS algorithm). See my comments below, please.\n\n- The differential is missed in Equations (7) and (8).\n\n- The state-of-the-art discussion that the authors do not take into account is, for instance,\n\nM. F. Bugallo et al, \"Adaptive Importance Sampling: The Past, the Present, and the Future\", IEEE Signal Processing Magazine, Volume 34, Issue 4, Pages: 60-79, 2017.\n\nM. Bugallo et al, \"Adaptive Importance Sampling in Signal Processing\", Digital Signal Processing, Volume 47, Pages: 36-49, 2015.\n\nJ. M. Cornuet, et al, “Adaptive multiple importance sampling,” Scandinavian Journal of Statistics, vol. 39, no. 4, pp. 798–812, December 2012.\n\nO. Cappe et al, “Population Monte Carlo,” Journal of Computational and Graphical Statistics, vol. 13, no. 4, pp. 907–929, 2004.\n\nV. Elvira et al, \"Improving Population Monte Carlo: Alternative Weighting and Resampling Schemes\", Signal Processing Volume 131, Pages: 77-91, 2017\n\nG.R. Douc et al, “Minimum variance importance sampling via population Monte Carlo,” ESAIM: Probability and Statistics, vol. 11, pp. 427–447, 2007.\n\nL. Martino et al, \"Layered Adaptive Importance Sampling\", Statistics and Computing, Volume 27, Issue 3, Pages: 599-623, 2017.\n\nO. Cappe et al, “Adaptive importance sampling in general mixture classes,” Statistics and Computing, vol. 18, pp. 447–459, 2008.\n\nL. Martino et al, \"An Adaptive Population Importance Sampler: Learning from the Uncertanity\", IEEE Transactions on Signal Processing, Volume 63, Issue 16, Pages 4422-4437, 2015.\n\n\n", "summary_of_the_paper": "The authors propose a new importance sampling scheme based on a transformation and relying on Normalizing Flows.", "summary_of_the_review": "The paper could contain interesting material but requires additional work on order to be published.\n", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "4: All of the claims and statements are well-supported and correct.", "technical_novelty_and_significance": "2: The contributions are only marginally significant or novel.", "empirical_novelty_and_significance": "2: The contributions are only marginally significant or novel."}, "extra_text": {}}, {"review_id": "Uk7lMBvG7IY", "reviewer": "Reviewer_wAEt", "summary": "", "questions": "", "limitations": "", "rating": 5, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "# Pros\n\n1. The paper introduces ZüNIS, a fully automated software library for easy NIS. \n\n2. ZüNIS offers flexible and data-efficient training of normalizing flows in the context of importance sampling. \n\n3. The paper establishes a rigorous link between the true loss function and its Monte Carlo approximation (which can be used for training flows). \n\n# Cons\n\n1. The paper is in most parts a description of the ZüNIS package but does not contain much original research. \n\n2. The approach is focusing on a single target $f$, whereas other sampling algorithms such as MCMC typically aim to generate samples that can be used for multiple targets. \n\n3. Comparison with other sampling approaches (MCMC or quasi Monte Carlo) is missing. \n\n# Comments and Questions\n\n- In case of expensive integrands $f(x)$, the authors propose to sample from an auxiliary pdf $q(x)$ to approximate the intractable variance of the importance sampling estimator (eq. 10). This separates sampling (to estimate the intractable integral) from the trained model. However, how to choose the auxiliary pdf $q$? Hasn't the problem of solving an intractable integral just been cast into another intractable integral? Why should it be easier to choose a suitable $q$ than to choose a suitable proposal for the main task (i.e. solving $\\int_\\Omega f(x) dx$)?\n\n- Also the related remark about deep $Q$-learning on page 4 remains quite obscure to me. Are you using a $Q$-learning strategy? \n\n- Against which implementation of VEGAS do you compare ZüNIS? There seems to be an enhanced version, VEGAS+, which performs superior to the original formulation (see Lepage, 2021). \n\n- I could not reproduce your results by following the instructions in \"Reproducibility statement\" due to import errors. After fixing some import errors, new import errors occurred such that I finally gave up. \n\n- It seems that figures 3, 4 and 6 are not discussed in the main text. \n\n# Typos\n\n- page 1: \"predictins\" \n\n", "summary_of_the_paper": "The paper describes ZüNIS, a neural importance sampling (NIS) library. The main application area is high-energy physics (HEP) simulation which is used to, for example, interpret experiments at the LHC. Traditionally, integrals arising in HEP simulations are solved numerically with the VEGAS algorithm. VEGAS is an adaptive importance sampler that neglects correlations between variables and optimizes a coordinate-wise discrete transformation of the integration domain. VEGAS potentially suffers from large variances in the Monte Carlo approximation and violations of the underlying assumptions. NIS offers a remedy in that it can adapt itself more flexibly to the integrand and thereby achieves smaller Monte Carlo errors. The authors propose a simple extension of NIS that allows them to reuse samples when the evaluation of the integrand is very complex and expensive (as is the case in HEP simulations). They show that ZüNIS outperforms the VEGAS algorithm on various toy examples and simple HEP processes. \n", "summary_of_the_review": "The paper mainly describes the ZüNIS package but does not contain much original research. ", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "4: All of the claims and statements are well-supported and correct.", "technical_novelty_and_significance": "2: The contributions are only marginally significant or novel.", "empirical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work."}, "extra_text": {}}, {"review_id": "RQClelT_lB8", "reviewer": "Reviewer_9xQC", "summary": "", "questions": "", "limitations": "", "rating": 3, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "It is unclear what is the major contribution of this paper. The library ZÜNIS seems an important tool in the field of HEP, but is discussed very briefly. The loss function (10) seems new, but it is unclear for readers why such an auxiliary probability distribution q is required and q can be simply selected as a uniform distribution.", "summary_of_the_paper": "The paper investigate sampling problems in HEP simulations and neural importance sampling (NIS) based methods. In addition, an NIS algorithm developed based on the ZÜNIS, a PyTorch library, is introduced. The performance of the proposed algorithm is demonstrated by toy examples and some numerical experiments of HEP simulations.", "summary_of_the_review": "The innovation is insufficient.", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "3: Some of the paper’s claims have minor issues. A few statements are not well-supported, or require small changes to be made correct.", "technical_novelty_and_significance": "1: The contributions are neither significant nor novel.", "empirical_novelty_and_significance": "1: The contributions are neither significant nor novel."}, "extra_text": {}}, {"review_id": "LkrPFin2GAW", "reviewer": "Reviewer_TzbE", "summary": "", "questions": "", "limitations": "", "rating": 3, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "# Strengths\n\n* The background on IS and NIS is very didactic and has tutorial value.\n* The toolbox is modular and straightforward use.\n\n# Weaknesses\n* There is no discussion regarding activation functions and their role in the behavior of the resulting distribution (IS proposal). \n* The technical novelty is very limited. The only innovation is a gradient trick in Section 3.2 to make training more stable.\n\n\n# Observations\n* In Section 2.3 onwards, the authors use the notation $p(x, \\theta)$ — and $q(x, \\theta)$. Since $\\theta$ parameterizes $p$ and $x$ is the only random variable here, the usual notation is something like $p(x;\\theta)$ or $p_{\\theta}(x)$.", "summary_of_the_paper": "This manuscript proposes a framework to compute integrals of arbitrary functions over a finite volume (uniform distribution). The toolbox is a modular implementation of Neural Importance Sampling (NIS, Muller et al. 2018). The experiments demonstrate the toolbox and its usefulness compared to a popular software (VEGAS) and uniform sampling.\n", "summary_of_the_review": "This is work proposes an organized toolbox and has very limited novelty. Therefore, I’m suggesting a rejection. However, since I’m not a connoisseur of high energy physics, I am not able to exactly gauge the utility of this toolbox to the specific demographic.", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "4: All of the claims and statements are well-supported and correct.", "technical_novelty_and_significance": "1: The contributions are neither significant nor novel.", "empirical_novelty_and_significance": "2: The contributions are only marginally significant or novel."}, "extra_text": {}}]

Author rebuttal:

Meta-review: The paper considers neural importance sampling (that is, importance sampling with a trained flow proposal) and its application to high-energy physics. The two contributions of the paper are: (a) a methodological improvement in the training of the proposal; (b) a description of a software library that implements the framework.

All reviewers were critical of the paper and recommended rejection. The main issue raised was that the methodological contribution was not novel or significant enough, and not sufficiently evaluated. The authors disagreed with the reviewers that the methodological contribution was not significant enough, but they acknowledged that the first version of the paper did not present the contribution clearly; consequently, they submitted a heavily revised second version following the reviewers' feedback.

Although it seems that the second version is an improvement over the first one, it's clear that the paper requires a second round of reviewing to ascertain whether it satisfies the requirements for acceptance. At this stage, the consensus among reviewers remains that the paper should be rejected. For that reason, I cannot recommend acceptance to ICLR. I sincerely hope the reviewers' feedback will be useful to the authors for a future submission to a different venue.
Decision: Reject