A practical guide to unbinned unfolding

TL;DR

This paper provides a guide on using machine learning for unbinned unfolding in particle physics to improve data analysis.

Contribution

The paper offers practical recommendations for implementing unbinned unfolding techniques in real-world particle physics experiments.

Findings

Unbinned unfolding allows for higher-dimensional analyses compared to traditional binned methods.

Machine learning techniques enable more flexible and accurate data unfolding in particle physics.

Researchers have successfully applied these methods to real experimental data.

Abstract

Unfolding, in the context of high-energy particle physics, refers to the process of removing detector distortions in experimental data. The resulting unfolded measurements are straightforward to use for direct comparisons between experiments and a wide variety of theoretical predictions. For decades, popular unfolding strategies were designed to operate on data formatted as one or more binned histograms. In recent years, new strategies have emerged that use machine learning to unfold datasets in an unbinned manner, allowing for higher-dimensional analyses and more flexibility for current and future users of the unfolded data. This guide comprises recommendations and practical considerations from researchers across a number of major particle physics experiments who have recently put these techniques into practice on real data.

Figures

Illustration of the OmniFold method. In Step 1, MC is corrected to match data at the detector level, producing a weighting function $\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$w_1(\vec {x}_{\text {reco}}^{\text {MC}})$$\end{document}$ . In Step 2, a new function $\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{ma

An example of the documentation from the ATLAS $\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$Z+$$\end{document}$ jets measurement [[22](#CR22)] shows how to use the unbinned data to configure plots of not only the unfolded result, but also the uncertainty breakdown (left) and correlation matrix (right) using bins of one’s choice for a given observable such as leading lepton $\documentclass[12pt]{minimal} \usepac

Full text

Introduction

Science is driven by the pursuit of deeper truths about our Universe, but scientists must extract these truths using datasets collected by instruments that carry inherent measurement errors and resolution limits. To assess the validity of a theory by directly comparing it with experimental data, scientists must therefore account for the expected changes that will occur when a pristine theoretical prediction is altered by the imperfections of the measurement device. If the physics process to be measured and the model to be tested are both known a priori, the standard approach used by researchers is “forward modeling”: taking a theoretical prediction and then modifying it using a carefully-constructed, high-accuracy detector simulation that requires expert knowledge of the detector geometry and its calibrations in order to directly compare the prediction with real data. Each time one wants to test a new theory or physics process, however, this time- and compute-intensive process must be fully re-run, making it impractical for testing multiple theories, particularly if a new hypothesis is developed long after the measurements have already taken place.

As an alternative to this paradigm, a real-world dataset can instead be corrected to try to remove any detector distortions from the data, thereby resulting in measurements that can be directly compared between experiments and theoretical predictions. Such a dataset could then be tested against any number of theoretical hypotheses without having to repeat the detector simulation each time. The nuanced process of removing detector effects from data is widespread across many areas of science and can be variously referred to as “deconvolution”, “denoising”, “reconstruction”, or “inverse modeling”. In High-Energy Physics (HEP), it has historically been called “unfolding”.

For decades, unfolding in HEP has been performed on binned datasets of no more than a handful of variables at a time. Algorithms such as Iterative Bayesian Unfolding [1] frame the experimental process as a series of transformations including a detector response matrix applied to the “true” dataset, i.e. an observable prior to the interaction with the detector – or, equivalently, what a detector with perfect resolution and detection efficiency would measure. The detector response matrix encodes bin-by-bin measurement effects such as reconstruction errors, energy smearing, and background noise. Unfolding, in this framing, means effectively inverting this detector response matrix so that it can be applied instead to experimental data and thereby produce a truth-level dataset. Approximating this inverted matrix can produce numerical instabilities that can yield unreliable results, particularly as the number of variables – and therefore the number of independent degrees of freedom – increases. As a result, unfolded results from traditional binned strategies are restricted to fixed histogram bins for only a few simultaneous observables measured at once.

Over the past few years, HEP researchers have proposed new unfolding techniques that use machine learning pipelines to perform unfolding in an unbinned manner [2–18], allowing for entire Monte Carlo simulated datasets to be corrected event-by-event such that the entire dataset aligns with the target. In addition to the inherent flexibility of removing a predetermined choice of bins from a measurement, these methods use neural networks that can readily process high-dimensional inputs, meaning that the unfolding can be performed for dozens or more variables at once without significant computational overhead. One such method that has gained recent prominence in the field is called OmniFold [6, 7]. OmniFold is a density reweighting method that determines a set of weights that can transform a simulated dataset event-by-event to match the target, i.e. the data one would see with an idealized detector. These weights are derived from the ratio of the likelihoods of two samples. These ratios can be directly approximated by harnessing neural network classifiers, as these methods are known to naturally learn likelihood ratios in order to optimally discriminate between classes.

Unbinned unfolding has now progressed well beyond simulation-only and proof-of-concept results. Between mid-2021 and mid-2025, a number of public results have made use of unbinned unfolding methods for measurements on datasets from at least five different HEP detectors: ATLAS, CMS, H1, LHCb, and STAR, as well as a study on highly realistic simulation from one accelerator-based neutrino experiment, T2K. This white paper is designed to further facilitate broader adoption of unbinned unfolding in the field by synthesizing the lessons learned from researchers involved in each of these measurements. In Sect. 2, we summarize the core methodology of OmniFold, the unbinned unfolding technique that enabled each of these public measurements. In Sect. 3, we summarize the specific choices made by the various analyses in areas including hyperparameter optimization, ensembling, uncertainty calculation, validation, and presentation of the final results. Finally, in Sect. 4, we conclude and propose some future directions for unbinned unfolding in HEP using machine learning.

Unbinned unfolding

Fig. 1. Illustration of the OmniFold method. In Step 1, MC is corrected to match data at the detector level, producing a weighting function \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$w_1(\vec {x}_{\text {reco}}^{\text {MC}})$$\end{document} . In Step 2, a new function \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$w_2(\vec {x}_{\text {true}}^{\text {MC}})$$\end{document} is learned based only on particle-level quantities. The method proceeds iteratively, refining the function \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$w_2(\vec {x}_{\text {true}}^{\text {MC}})$$\end{document} such that the particle-level MC can be reweighted to give event yields and kinematics that match those observed in the data

Each unbinned analysis summarized in this work makes use of the classifier-based machine learning method called OmniFold [6, 7], which is illustrated in Fig. 1 and briefly reviewed here. OmniFold requires two data samples that provide distributions of three different sets of observables in total:

1. \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\vec {x}^{\text {MC}}_\text {true}$$\end{document} : a Monte Carlo simulation sample containing a list of events with a number of “truth-level” observables, i.e. particle-level properties as they would be captured by an idealized detector;
2. \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\vec {x}^{\text {MC}}_\text {reco}$$\end{document} : a Monte Carlo simulation sample containing that same list of events, but instead with “reco-level” observables, i.e. the properties of those same events as they would be captured by a realistic detector; and
3. \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\vec {x}^{\text { data}}_\text {reco}$$\end{document} : the experimental data that, by definition, only has “reco-level” observables. Every dataset also specifies an initial “weight” for each event. For the Monte Carlo datasets \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\vec {x}^{\text {MC}}_\text {true}$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\vec {x}^{\text {MC}}_\text {reco}$$\end{document} , these weights are generally given by the Monte Carlo generator, but for the real data \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\vec {x}^{\text { data}}_\text {reco}$$\end{document} these weights are all set to 1. An event indexed i in the MC sample can hence typically be written as \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$(w_{\text {init},i},\vec {x}^{\text {MC}}_{\text {true},i},\vec {x}^{\text {MC}}_{\text {reco},i})$$\end{document} , and the task of OmniFold is to use these three different sets of event properties to adjust the initial weight \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$w_{\text {init},i}$$\end{document} for each MC event such that the weighted sample of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\vec {x}^{\text {MC}}_\text {reco}$$\end{document} agrees with data. At that point, the weighted sample of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\vec {x}^{\text {MC}}_\text {true}$$\end{document} will provide an approximation (or measurement) of the target, i.e. \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\vec {x}^{\text { data}}_\text {true}$$\end{document} . OmniFold and other likelihood-based methods do not change the observables themselves in the dataset – only the weights. This is in contrast to unfolding based on generative models that change both the dataset observables and the weights.

OmniFold makes frequent use of a useful feature of neural network binary classifiers: their learned decision functions can be used to construct a smooth reweighting function between the two datasets that participate in the binary classification. This property holds for a broad range of loss functions (see e.g. Ref. [19]), but the binary cross-entropy loss function is a popular choice. With this choice of loss function, a classifier trained to distinguish whether a set of events \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\vec {x}$$\end{document} are drawn from either probability distribution \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$p_A(\vec {x})$$\end{document} or \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$p_B(\vec {x})$$\end{document} will learn a function \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$f(\vec {x})$$\end{document} that minimizes:

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} \mathcal {L}_{\text {BCE}}[f] = - \int d\vec {x} \left( p_A(\vec {x})\log (f(\vec {x})) + p_B(\vec {x}) \log (1-f(\vec {x}))\right) \end{aligned}$$\end{document}

The function that minimizes this expression will satisfy, for any small variation \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\delta f(\vec {x})$$\end{document} :

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} \delta \mathcal {L}_{\text {BCE}}&= \int d\vec {x} \left( -\frac{p_A(\vec {x})}{f(\vec {x})} + \frac{p_B(\vec {x})}{1 - f(\vec {x})} \right) \delta f(\vec {x}) = 0, \end{aligned}$$\end{document}

therefore the integrand must be identically zero:

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} -\frac{p_A(\vec {x})}{f(\vec {x})} + \frac{p_B(\vec {x})}{1 - f(\vec {x})} = 0. \end{aligned}$$\end{document}

Solving this equation, we find that we can reconstruct the likelihood ratio:

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ \frac{f(\vec {x})}{1-f(\vec {x})} = \frac{p_A(\vec {x})}{p_B(\vec {x})}, $$\end{document}

This function of the neural network classifier f(x) can then be used to reweight one sample to look like the other by adjusting the event weight of each data point appropriately. In practice, this method will asymptotically approximate the likelihood ratio \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$p_A(\vec {x})/p_B(\vec {x})$$\end{document}  [20, 21].

Like some binned unfolding methods that have been widely used in HEP for decades, OmniFold operates through multiple iterations. In fact, when applied to binned data, OmniFold is mathematically equivalent to Iterative Bayesian Unfolding [1]. Each iteration of OmniFold consists of two steps both of which involve estimating a reweighting function using the procedure described above:

• Step 1: First, a classifier is trained to distinguish between distributions of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\vec {x}^{\text {MC}}_\text {reco}$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\vec {x}^{\text { data}}_\text {reco}$$\end{document} . This classifier is used to construct a reweighting function \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$w_1(\vec {x}^{\text {MC}}_{\text {reco}})$$\end{document} that can be applied to make \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\vec {x}^{\text {MC}}_{\text {reco}}$$\end{document} match \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\vec {x}^{\text { data}}_{\text {reco}}$$\end{document} at the detector level.
• Step 2: Next, a new classifier is trained to distinguish between distributions of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\vec {x}^{\text {MC}}_\text {true}$$\end{document} and the adjusted events from the previous step taken at the particle level: \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\vec {x}_{\text {true}}^{\text {MC}}$$\end{document} weighted by \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$w_1(\vec {x}_{\text {reco}}^{\text {MC}})$$\end{document} . A new reweighting function \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$w_2(\vec {x}^{\text {MC}}_{\text {true}})$$\end{document} is learned that can be applied to make \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\vec {x}^{\text {MC}}_{\text {true}}$$\end{document} statistically match the target \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\vec {x}^{\text { data}}_{\text {true}}$$\end{document} at the particle level. These two steps are then repeated for multiple iterations, where the next iteration instead takes \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\vec {x}_\text {reco}^\text {MC}$$\end{document} weighted by \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$w_2(\vec {x}_{\text {true}}^{\text {MC}})$$\end{document} as its input to compare with data. The updated reweighting function becomes a product of the previous one until a predefined number of iterations are performed and the method stops.

Practical considerations

Experimental overview

Table 1 summarizes the eleven analyses considered.Table 1. An overview of the recently-published experimental results that use unbinned unfolding methods. For most results, the unfolded dimensionality is the same at reconstruction-level and at truth-level, but * indicates the analyses for which the unfolded dimensionality at reconstruction-level was the full phase space, but the truth-level unfolded result is 8-dimensional. Full details of the phase space for each measurement, including \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\eta $$\end{document} selections, are listed in the individual papersExperimentPaper linkDimensionsFinal stateMomentum selectionATLAS [22]2405.2004124Z+jets \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$p_{\text {T}}^{\ell \ell } > 200$$\end{document} GeVATLAS [23]2502.020626Dijets \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$p_{\text {T}}^{j_1} > 240$$\end{document} GeV & \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$p_{\text {T}}^{j_1} < 1.5\ p_{\text {T}}^{j_2}$$\end{document} CMS [24]2505.178508Minimum bias \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$> 2$$\end{document} charged particles with \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$p_{\text {T}} > 0.5$$\end{document} GeVH1 [25]2108.123768High \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$Q^2$$\end{document} DIS \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$Q^2 > 150$$\end{document} GeV \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^2$$\end{document} H1 [26]2303.1362010High \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$Q^2$$\end{document} DIS \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$Q^2 > 150$$\end{document} GeV \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^2$$\end{document} H1 [27]2412.140928High \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$Q^2$$\end{document} DIS \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$Q^2 > 150$$\end{document} GeV \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^2$$\end{document} H1 [28]H1prelim-25-031VariableHigh \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$Q^2$$\end{document} DIS \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$Q^2 > 150$$\end{document} GeV \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^2$$\end{document} LHCb [29]2208.116914Z+hadrons in jets \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$20< p_{\text {T}}^{j} < 100$$\end{document} GeV and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$p_{\text {T}}^{h} > 0.25$$\end{document} GeVSTAR [30]2307.077186Jets \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$20< p_{\text {T}}^{j} < 50$$\end{document} GeVSTAR [31]2403.139217Jets in heavy ions \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$20< p_{\text {T}}^{j} < 45$$\end{document} GeVT2K [32]2504.068576Muon + Proton \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$p^p > 450$$\end{document} MeV for single transverse variables

Hyperparameter optimization

The unbinned unfolding process includes a number of hyperparameters that affect either the core unfolding methodology (e.g. number of unfolding iterations) or the neural network training (e.g. batch size). The OmniFold codebase [33] includes a suite of recommended default hyperparameters that were generally sufficient for some analyses. Other analyses chose to optimize hyperparameters based on the particular needs of their experiments. Those who did choose to perform a hyperparameter optimization generally tracked the performance of the model in unfolding with two different datasets of MC simulations.

Unfolding hyperparameters

Number of iterations: OmniFold [6], like its binned analog IBU [1], is an iterative method. However, the number of iterations needed to see good performance can vary. To estimate this number, most analyses used dedicated studies with a known target to measure the convergence of the unfolded dataset with the target over a large number of iterations, sometimes reaching up to 70 or 100. In practice, far fewer iterations were used for the actual analysis result, typically around 5, due to limitations in the detector resolution. One notable exception is the T2K study [32], which preferred as many as 20–40 iterations. It also explored using per-event weight changes after each iteration to deal with the still open question of how to choose the appropriate number of iterations in an unbinned manner. The default choice in the OmniFold software is set to 3 iterations. Analyzers can also report the unfolding results for different numbers of iterations to assess the impact of this choice, as in Refs. [34, 35].

Network initialization: Most analyses chose to train the networks for each unfolding step from scratch with the onset of each new iteration, which is also the default choice in the OmniFold software. However, using a pre-trained network as a starting point for training subsequent steps can potentially reduce the overall training time and improve the stability of the result, particularly in cases for which there is a large imbalance in statistics between simulation and data [28], or where the available data statistics are just generally low [32]. For example, a classifier with the same architecture as the models used in steps 1 and 2 of OmniFold might be pre-trained to distinguish between two different MC samples, as in Ref. [36]. Alternatively, a pre-trained particle physics foundation model (see e.g. Refs. [37–39] and Ref. [40] for a recent review) could be fine-tuned at each iteration.

Step 2 task: Step 2 of the unfolding procedure is designed to learn a function that maps the truth-level MC to the target of unfolded data. At each iteration, however, this task can be framed either as performing this mapping in full or in part. “In part”, in this case, means learning to map from the previous iteration’s reweighted MC to the target and then composing each of these mappings together after training to learn the full mapping from initial MC simulation to the target. The default behavior in OmniFold is to learn this mapping in full. In the CMS analysis [24], researchers found that learning the mapping from the previous iteration’s reweighted MC sample yielded worse performance and larger variance in the OmniFold weights, though this behavior has not yet been tested extensively.

Neural network hyperparameters

Network architecture: The classifiers across the various analyses took the form of dense neural networks with a small number of of hidden layers – typically three, although some analyses also explored architectures with two or four layers. The T2K study [32] found measurably worse performance with smaller networks. Hidden layers mostly used the ReLU activation function and the final output layer used a sigmoid activation. Hidden layers generally had a size of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mathcal {O}(100)$$\end{document} , e.g. 100 or 200 nodes. Batch size: Larger batch sizes of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mathcal {O}(10^3)$$\end{document} events are preferable in order to maximize GPU usage and improve stability of the results. The default recommendation is 128, but in practice, different analyses selected batch sizes ranging from 1028 [24] to 50,000 [30].

Number of training epochs: Instead of using a fixed number of epochs, most analyses implemented early stopping using either a pre-determined patience of e.g. 10 epochs or a minimum \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\varDelta $$\end{document} , i.e. a lower bound on the change in the loss that counts as an improvement to the model. Early stopping was the most common form of regularization used, but a few analyses also employed dropout and batch normalization.

Preprocessing

Data preprocessing can have major effects on the ultimate performance of a neural network training. Each analysis had to determine not only how to represent the features used as inputs, but also how to standardize the inputs and weights.

Input features: Many analyses used jet features represented as four-vectors, e.g. \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$p_T$$\end{document} , \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\eta $$\end{document} , and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\phi $$\end{document} . In some cases, these features were represented relative to another reference frame such as the recoiling electron [28]. A few analyses chose to apply Graph Neural Networks (GNNs) in order to process the particles as 3D point clouds [26, 28]. The ATLAS \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$Z+$$\end{document} jets analysis [22] chose to represent \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\phi $$\end{document} coordinates as a combination of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\text {sin}(\phi )$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\text {cos}(\phi )$$\end{document} after a dedicated study that suggested this representation was more effective at reducing discontinuities in the final results. In addition to raw input features, summary statistics such as moments and cumulants can also be successfully unfolded, even higher-order ones with relatively small values (e.g. \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mathcal {O}(10^{-5}$$\end{document} )) [23]. Inputs were generally standardized using their z-scores, i.e. \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$z = (x - \mu ) / \sigma $$\end{document} , before using them in the training process.

The T2K study [32] identified two outgoing particles for each event, the muon and proton, and used the angles and forward momentum of these particles as inputs. However, in events for which there was no outgoing proton – whether because it was not reconstructed or because no proton was kicked out at truth level – placeholder values were used for the proton kinematics. Additionally, this analysis separated events into different samples in detector space based on the locations within the detector that the particles were reconstructed, and separated them into different reaction topologies in truth space. This qualitative information was provided as one-hot encoded input to the neural network.

Negative weights: Simulated samples of high-energy physics data can sometimes be associated with MC weights representing wide magnitude ranges, even including negative values, which can introduce difficulties when using these weighted samples in training neural networks. In the ATLAS Z+jets measurement [22], several of the input MC samples suffered from a significant fraction of negative event weights and large spread of weights. For these samples, the original MC weight was replaced by a new weight following a procedure similar to the method described in Ref. [41]. The resulting samples had positive weights with a reduced spread (standard deviation), yielding datasets that are statistically compatible to those produced using the original weights.

Calculating the differential cross-section: Most analyses performed the unfolding task on normalized samples such that the neural networks learned the relative shape differences between the datasets. This choice allows for the use of standardization in the weights themselves as well as the input features. The measured differential cross-sections can then estimated by scaling the weights based on a calculation using the efficiency and fiducial factor from the nominal MC sample as well as the fiducial cross-section and number of data events:

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ \frac{\mathcal {L}\cdot \sigma _{\text {fiducial}}\cdot \epsilon }{f_{\text {fiducial}}} = n_{\text {data}}. $$\end{document}

These corrections are implicitly accounted for if the unfolding is done without normalizing the samples, as was done in the T2K study [32]. Alternatively, if the goal is to measure the shape of an observable, the measurement can be normalized again after unfolding such that the weighted area sums to 1.

Simultaneously unfolding event- and particle-level features: Unfolded features can potentially include both aggregate and individual characteristics of event constituents, e.g. a jet and its constituent particles. In the LHCb measurements [29], to account for this, an initial weight of unity for each particle was scaled down by the particle multiplicity within the jet to allow the network learn the substructure of jets and consistently normalize the synthetic data by the number of jets.

Background and acceptance effects

These analyses broadly benefited from a lack of major background sources that would require a detailed background subtraction process. However, their methodologies can in some cases be extended for analyses with nontrivial backgrounds. We note that “backgrounds” can comprise irreducible backgrounds, i.e. physics processes with the same final states as the signal, as well as reducible backgrounds, i.e. events that have been misclassified as signal events due to detector artifacts or other errors during particle reconstruction and identification.

Effects from irreducible backgrounds: Backgrounds were estimated to be small or negligible for the analyses covered in this report, and were generally either ignored or were treated by assigning an overall uncertainty. For example, the ATLAS \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$Z+$$\end{document} jets measurement [22] had a primary background of top events, but these were small enough ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$<0.25\%$$\end{document} of total events) to simply subtract them from a set of pseudodata, run the full unfolding procedure, and then take the difference between this and the nominal result as an uncertainty. The CMS measurement [24] targeted minimum bias proton-proton collisions, and the fiducial space for the measurement is as inclusive as possible. In the STAR heavy ions jets analysis [31], backgrounds were modeled using an embedding simulation to remove residual effects. To mitigate any discrepancies between data and MC due to choice of embedding strategy, an extra weight was then multiplied to the reco-level jets to make their multiplicity and luminosity distributions look similar to real minimum-bias data. For future users dealing with nontrivial irreducible backgrounds, we recommend incorporating these backgrounds into the initial MC dataset as MC events with negative weights such that in total, the initial MC dataset corresponds to data with background subtracted.

Acceptance and reducible background effects: These effects emerge due to events that pass the truth-level selections but not reconstruction-level selections, or vice-versa. Generally, the former is the dominant acceptance effect due to tracking or trigger inefficiencies. The default OmniFold [6] procedure handles these events by assigning them the average weight in their region of phase space, which has been found to yield better convergence than directly including them in the unfolding procedure [32]. Unfolding a larger fiducial region than the one reported in the final measurement can also help mitigate acceptance effects – for instance, the ATLAS \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$Z+$$\end{document} jets measurement [22] measures the phase space \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$p_T^{\ell \ell } > 190$$\end{document} GeV, but only reports results with \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$p_T^{\ell \ell } > 200$$\end{document} GeV in order to mitigate acceptance effects due to migrations across this threshold.

In the LHCb measurements of charged particles in jets [29], the track reconstruction efficiencies along with particle identification corrections were applied to each reconstructed particle in data being unfolded based on a look-up table prepared prior to unfolding.

In the STAR measurement of jet substructure variables [30], this correction is instead done after the unfolding, and involves applying an efficiency correction by dividing the normalized unfolded dataset over the efficiency function bin-by-bin. Backgrounds due to misidentified jets were estimated using simulations and were used to construct initial weights for the data by subtracting the fake rates [30].

The minor background from the minimum-bias CMS measurement [24] came from fake tracks induced by detector noise or misidentification of the algorithms, which resulted in events that passed the detector-level selection but failed the truth-level selection. Acceptance effects came from tracking inefficiency corresponding to events that passed the truth-level selection but failed the detector-level selection. Both effects were taken into account by adding two steps1 to each iteration of the unfolding, as described in Ref. [7].

The T2K study [32] dealt with both significant inefficiencies from acceptance effects and significant backgrounds, which are common features of measurements in neutrino detectors due to their large volumes and relatively low fidelity. The backgrounds and effects from non-uniform detector acceptance were dealt with by splitting the detector-level events into distinct samples based on the topology with which they were reconstructed in the detector, including background-enhanced samples to constrain background contributions to the overall signal. This qualitative sample information was included as additional input to the neural network used in the unfolding. Additionally, results on single transverse variables were only reported for events with proton momentum \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$>450$$\end{document} MeV, but the unfolding was performed over the entire available phase space without this restriction.

Ensembling

Most analyses found that ensembling improved the stability of the unfolded result. “Ensembling” refers to training multiple independent instances of the same model and then using the mean or median weight. In principle, this ensembling could occur at the level of binned observables, but ensembling the weights themselves facilitates publishing a fully unbinned dataset that does not depend on any particular choice of binning. The standard error on this quantity across the trained models is then included as a source of uncertainty. This uncertainty can be thought of as measuring the impact of the stochastic nature of the neural network training – i.e. choice of random seed – on the final unfolded result. The magnitude of this uncertainty is generally small ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$<2\%$$\end{document} ). Most analyses used some form of ensembling, with numbers of ensembles typically ranging from 4 to 10. Both the ATLAS \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$Z+$$\end{document} jets measurement [22] and the ATLAS jet track functions measurement [23] used an ensemble of 100 models, which helped to mitigate the compounding effects of these statistical fluctuations over a large number of uncertainties. Ensembling was usually performed by re-running the entire result, but as an alternative, the ensembling could instead be done independently for Step 1 and Step 2 of each unfolding iteration [36].

Uncertainties

Despite the introduction of neural networks underlying the unfolding procedures in these analyses, the procedure of estimating uncertainties is, with a few exceptions, relatively standard for particle physics measurements. These uncertainty estimation strategies include bootstrapping for estimating statistical uncertainties and calculating the effects of varying the MC sample for estimating systematic uncertainties. One new uncertainty included in these measurements was an uncertainty associated with the neural network initialization, i.e. the inherent variation in repeating the neural network training with a different random seed but otherwise an identical configuration. These variations are typically small, but they affect every other reported uncertainty, so they are important to constrain. Additionally, while OmniFold will asymptotically (i.e. with enough training data and a flexible enough network architecture and training paradigm) approach the maximum likelihood estimate of the true particle-level distribution and should therefore not be dependent on mismodeling of the true prior (see Appendix A of Ref. [6]), in practice the finite dataset size, model size, and number of iterations can still yield effects due to the mismodeled particle-level prior in the MC sample. These effects are measured as part of the overall unfolding uncertainty.

ATLAS unfolding uncertainty

In the ATLAS \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$Z+$$\end{document} jets measurement, the unfolding uncertainty consisted of two parts: an uncertainty due to the truth-level prior as well as an uncertainty due to the potential mismodeling of the detector effects [22]. These effects were estimated separately and reported as two separate sources of unfolding uncertainty. Together, the unfolding uncertainty represented the dominant uncertainty across most of the measured phase space.

The first component of the ATLAS unfolding uncertainty, which measures the sensitivity of the result to the modeling of the particle-level sample, is estimated via a data-driven process: the MC is reweighted at truth-level such that all of its reco-level observables agree with data, and the result of unfolding using this dataset is compared with the nominal unfolding result. Crucially, the reweighting function used here is a sequence of one-dimensional Gaussian-kernel functions, not the OmniFold procedure that is itself being probed. The second component of the ATLAS unfolding uncertainty, which measures the sensitivity of the result to potential mismodeling of the detector response for features not included in the unfolding, is sometimes called a “hidden variable” uncertainty. To estimate this component, a MC sample from a different generator is used in the unfolding, but it is first reweighted at truth-level to match the MC sample from the original generator as closely as possible.

CMS unfolding uncertainty

In the CMS measurement [24], uncertainty templates are implemented as weighting functions applied to nominal MC events. These functions are learned by classifiers trained to distinguish nominal and alternative MC samples, either at truth- or reco-level, depending on the source of the uncertainty. The templates for the CMS unfolding bias, which originate from the mismodeled truth-level prior on the target variables to be unfolded, are constructed as weighting functions \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$w_\alpha (\vec {x}_\text {true}^\text {MC})$$\end{document} at the truth level from the nominal MC sample to the MC sample from each alternative model \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\alpha $$\end{document} . When applied to the nominal MC events, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$w_\alpha (\vec {x}_\text {true}^\text {MC})$$\end{document} makes their truth-level samples match the ones from the alternative model. Similarly, the templates for mismodeled detector responses, which come from the discrepancies between the reconstruction in the simulation and truth, are constructed as weighting functions \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$w_\text {det}(\vec {x}_\text {reco}^\text {MC})$$\end{document} at the detector level, mapping the nominal MC sample to the MC sample from the same truth-level model but alternative detector simulation.

A two-stage reweighting scheme handles mismodeled priors in kinematic variables not directly unfolded. First, a truth-level weighting \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$v^1_\alpha (\vec {x}_\text {true}^\text {MC})$$\end{document} is applied to MC events from the alternative model \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\alpha $$\end{document} to match their truth-level samples to the ones from the nominal model. Then, a joint weighting of both truth and detector level, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$v_\alpha (\vec {x}_\text {true}^\text {MC},\vec {x}_\text {reco}^\text {MC})$$\end{document} , is applied to the nominal MC events to match their samples at both levels to the ones from the alternative samples weighted by \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$v^1_\alpha (\vec {x}_\text {true}^\text {MC})$$\end{document} . The resulting weighting function keeps the truth-level of the nominal sample while mapping its detector response to that of the alternative model. These reweighted templates allow smooth interpolation between models via nuisance parameters and enable uncertainty estimation via bootstrapping.

Validation

Many analyses were performed “blinded”, i.e. not considering the real data while configuring the analysis parameters. This meant using a different MC sample in place of the real data while developing the analysis, and then only using the real data after the analysis passed a number of predetermined validation criteria. Even if the analysis itself was not formally blinded, all analyses validated their unfolding procedure using alternate MC samples.

Pseudodata: The ATLAS \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$Z+$$\end{document} jets analysis [22] constructed a realistic set of “pseudodata” in place of the actual data to validate their procedure before obtaining results using the data. The pseudodata has a known underlying dataset \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\vec {x}_\text {true}$$\end{document} , so the measured OmniFold result of any quantity can be compared to the desired “target” value. To produce the pseudodata, alternate signal and background MC samples were reweighted to match data using a dedicated reweighting algorithm. The full analysis was performed using this pseudodataset, and the unfolded results were evaluated by analyzing the p-value distributions between the unfolded result and the target truth pseudodata. The agreement was checked using not only the 24 input observables, but also on some select composite observables (i.e. observables that can be calculated using the 24 but are not included in the unfolding), for a series of 2D distributions and for a series of phase space cuts.

The CMS analysis [24] used a few sets of “pseudodata” from detector simulations of alternative models that are not used in unfolding or uncertainty estimation. The pseudodata sets are generated from distinct generators, parton distribution functions, or underlying-event tunes. The unfolding results of pseudodata sets are compared with their corresponding truth-level samples. The agreements are used for hyperparameter optimization and have been confirmed for all the variables to be unfolded and their multi-dimensional distributions.

Validation: Technical closure tests and stress tests served to validate the unfolding methodology itself using different sets of MC samples. In the ATLAS Z+jets measurement, the stress tests were performed using a single MC sample to unfold itself, but the half used for data had additional stress weights applied. Two tests were performed: one where the weights were based on a deterministic function of the observables, and another where the weights were based on a stochastic function of the observables. In the T2K study [32], the test was performed using a set of fake data with reweightings that were functions of non-observable parameters. This resulted in the unfolding needing to learn a reweighting that is not a simple function of the quantities made available to it – a type of stress test that is frequently done for neutrino cross-section measurements.

In the CMS result, bias and coverage tests were performed to evaluate the bias of the unfolded results compared to the truth and to validate the coverage of the unfolding uncertainties. The tests were based on unfolding pseudodata using a frequentist framework. Multiple toy experiments of unfolding were performed, in which both the MC and pseudodata samples were randomly deviated from the original ones according to the considered uncertainties. The cumulative behaviors of these toy experiments indicate the fluctuations if the measurement is repeated several times. The bias was evaluated as the deviations between the central values from the toy experiments and the pseudodata truth, and the coverage was estimated as the frequencies that the toy experiment results are compatible with the pseudodata truth, considering the uncertainties.

Furthermore, the CMS analysis performed a bottom-line test to compare the differences between the unfolded results and the truth-level nominal MC sample to the ones between the (pseudo)data and the detector-level nominal MC sample. Since the unfolded results cannot be more discriminative to the models than the detector-level data, the former should always be smaller than the latter ones. A failure of the bottom-line test could come from the algorithm attributing the data-MC discrepancies from statistical fluctuations to the truth-level mismodeling and overfitting the given samples.

Presentation of results

Some community-sourced guidelines have already proposed how to best present unbinned measurements [3]. Across all of the results presented in this paper, however, though the unfolding procedure itself is performed unbinned, the unfolded results are presented in a binned histogram format for all but the ATLAS Z+jets measurement, which was presented in an unbinned manner.

Binned format: Measurements are presented as individual histograms corresponding to each unfolded dimension. The specifications and bin counts of these histograms can then be uploaded to HEPData.

Unbinned format: In the ATLAS \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$Z+$$\end{document} jets measurement [22], the final result is presented as a set of Pandas DataFrames in which each DataFrame row corresponds with a truth-level Monte Carlo event, inspired by the example set by Ref. [35] using CMS Open Data. More than one DataFrame is needed for the full ATLAS result in order to account for all sources of systematic uncertainty. Columns correspond with the 24 unfolded dimensions, i.e. truth-level particle quantities, along with hundreds of weights used to construct the nominal measurement and calculate uncertainties. These datasets are published on Zenodo [42] and are paired with a detailed codebase illustrating how to use them to replicate the results in the paper as well as how to produce new measurements (see Fig. 2). This presentation format was closely informed by the recent white paper that discusses the presentation of unbinned results [3]. The CMS results will also be released in a similar unbinned format.Fig. 2. An example of the documentation from the ATLAS \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$Z+$$\end{document} jets measurement [22] shows how to use the unbinned data to configure plots of not only the unfolded result, but also the uncertainty breakdown (left) and correlation matrix (right) using bins of one’s choice for a given observable such as leading lepton \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$p_T$$\end{document}

Computing requirements

Across each of these measurements, an individual unfolding procedure could, in most cases, be run on a single NVIDIA A100 GPU with 40 GB or 80 GB of RAM in \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mathcal {O}(1$$\end{document} hour), generally between 1 and 4 h. The main factors influencing this training time include dimensionality of the unfolding, number of trainable model parameters, and number of unfolding iterations. If ensembling is used, identical copies of the unfolding setup are run up to \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mathcal {O}(100)$$\end{document} times in order to reduce the uncertainty due to the stochastic nature of the training process. These ensembling runs are able to be fully trained in parallel, if compute availability permits.

A full physics measurement, however, includes the calculation of a number of sources of uncertainty, and the number of uncertainties included can directly scale the computational resources required. Computational complexity for systematic uncertainties depends on the uncertainty sources and estimation strategies, requiring \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mathcal {O}$$\end{document} (100) unfolding runs if bootstrapping is used. Statistical uncertainties estimated through bootstrapping can also require about 50–100 individual bootstrap runs.

In total, most analyses reporting significant computational needs for the unfolding portion of their experiments estimated that they used between approximately 500 and 10,000 GPU hours in total to perform their measurement once the full procedure had been fixed. Necessary computational resources depended heavily on the number of models used for ensembling and ranged from the STAR heavy ions analysis [31], which did not use any ensembling and was able to run on a single GPU, to the 24-dimensional ATLAS \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$Z+$$\end{document} jets measurement [22], which used an ensemble of 100 models to calculate each weight and took approximately 25,000 GPU hours in total.

Conclusions and outlook

The recent emergence of more than 10 different high-energy physics publications using unbinned unfolding methods underscores that these methodologies are undoubtedly publication-ready. These results cover a variety of final states, six different experimental detectors, a range of unfolded dimensions (from four to 24 to the full phase space). Importantly, these methods have also been tested and validated in a number of different scenarios and have yielded unbinned cross-section measurements using large-scale physics datasets and full sets of uncertainties. As this field continues to expand, we expect that many more analyses will opt to both perform and publish their measurements in an unbinned fashion using the frameworks these first publications have established as a guide.

Future directions

While unbinned unfolding has now been proven effective for a range of particle physics applications, there remain several interesting directions for the development of future work in this area, including:

• Pre-trained models: Would fine-tuning pre-trained models instead of training from scratch each time help improve the stability of the result as well as save computational resources, and should this become the new default methodology? (See Ref. [36] for some additional discussion.)
• Statistical uncertainties: Are there viable alternatives to bootstrapping, to reduce the computational requirements needed in calculating statistical uncertainties? How can we best capture the statistical uncertainty not only in terms of the limitation of the total number of data events, but also in terms of the resulting instability in the neural network training?
• Unfolding after many iterations: While unfolding with many iterations has not been widely studied across multiple detector settings, some analyzers observed that unfolding performance suddenly plummeted after \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mathcal {O}(100)$$\end{document} iterations. Though this behavior did not affect the quality of the results ultimately presented here, which all relied on small numbers of iterations (5 or fewer), this behavior could merit future study.
• Unfolding the whole phase space: While the methodology for unfolding variable-length events is well-established in theory, we are only recently starting to see examples of this put into practice (e.g. Ref. [28]). To what extent do these recommendations apply to unfolding the whole phase space? Some researchers observed that when unfolding the whole phase space with datasets that were too similar, the unfolding performance could oscillate with each iteration – what is responsible for this behavior, and how could it be mitigated?
• Goodness-of-fit tests for unbinned data: How can we best measure closure with a known target in an unbinned manner, therefore bypassing the typical \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\chi ^2$$\end{document} and p-value binned statistical tests? Existing metrics include the Wasserstein distance [43] and permutation tests [44], among others, but many open questions still remain in this area.
• Generative methods: Each of these measurements has employed a version of OmniFold [6, 7]. Could other ML-based unbinned unfolding methods, including generative methods, yield competitive results?
• Nontrivial background subtraction: While there are proposed methods for putting nontrivial background estimation into practice, there is a strong need for researchers to put these methods to the test and explore which are most effective.