Obtaining the distribution of a physical quantity is a frequent objective in experimental physics. In cases where the distribution of the relevant quantity cannot be accessed experimentally, it has to be reconstructed from distributions of correlated quantities that are measured, instead. This reconstruction is called deconvolution. While several approaches to deconvolution exist in experimental physics, the problem is rather unknown in data science. Both fields can benefit from each other, but notational differences tend to prevent this. In this work, we outline the unification of existing approaches in order to pave the way for new joint research directions. Cherenkov astronomy is employed as a use case which illustrates the deconvolution problem.
Obtaining the distribution of a physical quantity is a frequent objective in experimental physics. An accurate and reliable estimate of the sought-after distribution, e.g. the energy spectrum of an astrophysical particle source, is crucial to understand the underlying physical principles. It is also important for the verification—or exclusion—of certain theoretical models. In cases where the distribution of the relevant quantity cannot be accessed experimentally, it has to be reconstructed from distributions of correlated quantities that are measured, instead. This reconstruction is called deconvolution (or unfolding). It states a prominent problem faced in experimental physics.
While several approaches to deconvolution exist in physics [
Cherenkov astronomy studies the energy distribution of cosmic γ radiation to reason about the characteristics of celestial objects emitting such radiation. One prominent tool in this field are Imaging Air Cherenkov Telescopes (IACTs), which are placed at high altitudes on the surface of our planet. Since it is not viable to detect high-energy γ particles directly there, IACTs record Cherenkov light emitted by a cascade of secondary particles, instead. Such cascades, called air showers, are induced by γ particles interacting with Earth’s atmosphere. Since the γ radiation is not directly measured by the telescopes, deconvolution is applied to reconstruct the γ energy distribution from the related Cherenkov light recorded by IACTs.
γ Particle
The Cherenkov light emitted by each shower is recorded as a video sequence, from which features are extracted. These features represent each air shower on a higher level of abstraction, e.g. by its geometrical shape and size. Deconvolution adopts each feature as one observable quantity, from which the energy distribution of observed γ particles is reconstructed. 2
Formally, the task of deconvolution is to estimate a density function f : Y → R of
a random variable Y which represents a physical quantity, where Y denotes the
state space of Y . In Cherenkov astronomy, the most prominent Y is the energy
of γ particles of which the Cherenkov light is observed by an IACT. The task is
aggravated by the following deficiencies inherent to the experimental setup [
A reliable estimate of the relevant density f is only obtained, if these deficiencies are rectified, i.e. if f is appropriately reconstructed from the measured density g. The name of deconvolution is motivated by g being modeled as a convolution of f with a detector response function R : X × Y → R.
Z
Y g(x) =
R(x | y) · f (y) dy R provides the link between X and Y , incorporating background knowledge about their relation. Specifically, R(x | y) represents the conditional probability of measuring some x ∈ X when the actual value of the relevant quantity is y ∈ Y. To obtain f , this model of g has to be inverted—it has to be “deconvolved”. This is done by fitting f to the model of g, when g and R are given1.
Classical deconvolution algorithms solve a discrete variant of the continuous deconvolution problem from Equation 1. This discrete variant presented in Equation 2 is a linear system of equations where f , g, and R from the continuous problem are replaced by discretized versions. The relation between the continuous case (Eq. 1) and the discrete case (Eq. 2) is clarified by Equation 3, where we denote the discretized problem for each discretized state of the observed data g = (g1, . . . gJ ) ∈ RJ . (1) (2) (3) g = R f ⇔ gj = PI i=1 Rij fi 1 ≤ j ≤ J Here, each dimension of the vector g corresponds to a specific discrete value discg(x), and each dimension of the vector f corresponds to a specific discrete value discf (y). Since the domain X of g is usually multi-dimensional, clustering is employed to map each vector-valued observation x ∈ X to the index j of the cluster which the observation belongs to, effectively aggregating the multidimensional feature space into a single dimension. The conditional probabilities in R are estimated from a set of training observations. Two deconvolution methods are considered state-of-the-art, even though they have not yet been bench-marked on IACT data: 1 Readers familiar with deconvolution will notice that our nomenclature is slightly different from the notation that is conventional in particle physics. Specifically, the roles of the letters x and y are exchanged here for compliance with the notation that is common in machine learning: Here, y refers to the value of the target variable and x is the observed value. This notation clarifies the relation between the deconvolution problem and classification—a relation that is important to our unifying view.
Regularized Unfolding performs a maximum likelihood fit of Equation 2,
assuming that the absolute number of events for each discretized state of g(x)
is Poisson-distributed conditioned on f [
The novel Dortmund spectrum estimation algorithm (DSEA) [
Pˆ(Y ≡ i) =
X Pˆ(Y ≡ i|X = x) · Pˆ(X = x) x∈X ˆfi = 1 N
N X cM(i | xn) n=1 (4) (5) The deconvolution result is then improved by iterating the reconstruction, updating the distribution of the training set by weighting the examples according to the latest estimate of the target distribution. We suggest to scale the update step between iterations to ensure convergence of the algorithm. In fact, the original DSEA diverges from the true f after having found a suitable estimate in some iteration. Our suggestion is inspired by a common algorithmic pattern in numerical optimization, where the basic algorithm selects a search direction and a sub-routine called line search selects an appropriate step size for this direction. α(k) = √1 k
We have illustrated the problem of deconvolution, employing Cherenkov astronomy as a use case. Our presentation provides a basis for unifying the view on current state-of-the-art algorithms from particle physics in the context of machine learning. Against the background of many possible synergies arising from collaborative research on deconvolution, we have given one example of an algorithmic improvement inspired by numerical optimization.
Even though deconvolution methods have been surveyed elsewhere [