Jet-Parton Assignment in ttH Events using Deep Learning

TL;DR

This paper explores deep learning methods for accurately assigning jets to partons in simulated ttH events, achieving over 50% correctness, to improve top-Higgs coupling measurements.

Contribution

It introduces a deep learning approach for jet-parton assignment in ttH events, outperforming traditional methods and enhancing reconstruction accuracy.

Findings

Deep learning achieves over 50% correct assignments.

Compared methods show improved accuracy over traditional techniques.

Enhances the measurement of the top-Higgs coupling.

Abstract

The direct measurement of the top quark-Higgs coupling is one of the important questions in understanding the Higgs boson. The coupling can be obtained through measurement of the top quark pair-associated Higgs boson production cross-section. Of the multiple challenges arising in this cross-section measurement, we investigate the reconstruction of the partons originating from the hard scattering process using the measured jets in simulated ttH events. The task corresponds to an assignment challenge of m objects (jets) to n other objects (partons), where m>=n. We compare several methods with emphasis on a concept based on deep learning techniques which yields the best results with more than 50% of correct jet-parton assignments.

Tables

Table 1. Table 1 : Design parameters of the neural network.

Parameter	Value
Hidden layers	8
Units per layer	500
Activation	ELU
$L^2$ factor ($\lambda$)	$4\cdot {10}^{-6}$
Batch size	10,000
Optimizer	ADAM
Learning rate	0.001
Training epochs	200

Table 2. Table 2 : Efficiency ε ​ ( t ​ t ¯ ​ H ) 𝜀 𝑡 ¯ 𝑡 𝐻 \varepsilon(t\bar{t}H) of finding the correct jet-parton assignments in t ​ t ¯ ​ H 𝑡 ¯ 𝑡 𝐻 t\bar{t}H events, and efficiency ε ​ ( H ) 𝜀 𝐻 \varepsilon(H) of finding the correct jets originating from the Higgs boson decay for all methods.

Method	$\epsilon$($t\overline{t}H$) / %	$\epsilon$($H$) / %
${\chi }^2$ method	37	52
Boosted Decision Trees	45	57
Deep Learning	52	63

Equations

$${\mathrm{input}}_i$$

$$J(W,b|\hat{D},D_{\mathrm{DNN}}(W,b,x))=-\left[\hat{D}\cdot \log (D_{\mathrm{DNN}})+\omega \cdot \left(1-\hat{D}\right)\cdot \log (1-D_{\mathrm{DNN}})\right]+\lambda \sum _j^N{W}_j^2$$

$$\Delta R=\sqrt{\Delta {\varphi }^2+\Delta {\eta }^2},$$

$${\chi }^2={\left(\frac{m_H-{\tilde{m}}_H}{{\tilde{\sigma }}_H}\right)}^2+{\left(\frac{m_{t,lep}-{\tilde{m}}_{t,lep}}{{\tilde{\sigma }}_{t,lep}}\right)}^2+{\left(\frac{m_{t,had}-{\tilde{m}}_{t,had}}{{\tilde{\sigma }}_{t,had}}\right)}^2+{\left(\frac{m_{W,had}-{\tilde{m}}_{W,had}}{{\tilde{\sigma }}_{W,had}}\right)}^2$$

$$I(D)\equiv N_p(D>D^{truth})$$

Full text

1 Introduction

A direct measurement of the Higgs boson coupling to top quarks is considered one of the most important consistency tests of the Higgs particle discovered in 2012 within the Standard Model of particle physics [1, 2]. In this context, the mass-dependent coupling of the Higgs to matter particles and gauge bosons is expected to be largest for top quarks and close to unity.

A promising process to measure the strength of the Higgs-top coupling is to determine the cross section of top quark pair-associated Higgs boson production ($t\bar{t}H$) with the Higgs boson decaying into two bottom quarks and the top quarks decaying semi-leptonically. The expected cross-section is sufficiently large to be detected with the LHC experiments using the data recorded thus far. The analysis is challenged by an irreducible background from top quark pair production with radiated gluons that split into pairs of bottom quarks. A combination of advanced analysis methods is required to separate the $t\bar{t}H$ signal from background processes.

In this work we investigate a central processing step of the analysis, namely the assignment of the jets detected in a detector to the partons of the underlying hard scattering process. Correct assignments of the jets to the Higgs boson and to the top quarks enable calculation of suitable observables such as invariant masses and increase the separation power of signal and background. We compare several methods to determine the jet-parton assignment with emphasis on a method based on deep learning techniques.

The application of deep learning methods in various areas of fundamental research is receiving increasing attention. This is motivated by the impressive successes of the deep learning ansatz in handwriting recognition [3, 4], speech recognition [5], challenges with humans regarding image identification [6, 7], and in playing games [8]. For a review, see [9].

Deep neural networks have recently been investigated for challenges in particle physics. Various network designs have been successfully applied to extract a simulated, new exotic particle or Higgs boson signal from a background-dominated data sample [10, 11, 12], to identify the underlying parton flavor of a jet or to measure jet substructure [13, 14], and to reconstruct the neutrino flavor in neutrino-nucleus interactions [15].

As a new application category of deep neural networks, we investigate the capability to select the single correct assignment from a number of possibilities of associating $m$ objects with $n$ other objects ($m\geq n$).

As we investigate the Higgs boson decaying to bottom quark pairs and top quark pairs in the semi-leptonic decay channel, the desired final state consists of the lepton and neutrino of one of the $W$ boson decays, two light quark jets of the other $W$ decay, and at least four bottom quark jets ($b$-jets). A typical event situation is shown in Figure 1.

Note that the $b$-jets of the top and Higgs decays have similar transverse momenta such that simple kinematic cuts do not yield promising results.

The challenge thus consists of assigning $6$ out of $m\geq 6$ jets to $n=6$ partons. Higher-order radiation processes cause the number of jets to exceed $m=6$. We reduce the number of possible jet-parton assignments (permutations) by ignoring the interchange of the two jets associated with the Higgs decay or the $W$ decay, respectively. We also take advantage of bottom quark identification algorithms ($b$-tag) with high efficiency for correct $b$-tags and only a small probability of light quark jets being incorrectly $b$-tagged.

For events with $6$ jets where exactly $4$ jets are $b$-tagged, we assign the two untagged jets to the $W$ boson and obtain $4!/2=12$ permutations for the $b$-tagged jets. Only one of the $12$ permutations corresponds to the correct jet-parton assignment. Owing to higher-order radiation effects with correspondingly more jets observed in the detector, the majority of events allows for a larger number of possible permutations which typically reaches several hundred.

We will evaluate the efficiency of the deep neural network to find the correct permutations and compare them with the results of a boosted decision tree and a $\chi^{2}$ measure derived from the hypothetical masses of a permutation.

This work is structured as follows: Initially, the deep network design is explained. The simulated dataset is then introduced and input observables are presented in the following sections. We evaluate the efficiencies of finding the correct assignment with different methods. We also demonstrate the effect on reconstructed observables which are typically used in $t\bar{t}H$ analyses. Finally, we present our conclusions.

2 Neural network design

For the network architecture we use various technical concepts developed within deep learning research. As we work with a fixed number of features, i.e. input observables, our basic concept is a fully connected network consisting of $8$ hidden layers with $500$ nodes in each layer (Figure 2(a), Table 1).

In addition, we apply the so-called residual network concept [7] where the input to the layer $i$ not only consists of the output of the previous layer $i-1$, but also from another previous layer $i-j$, where we choose $j=3$:

[TABLE]

This residual concept is found to accelerate the training of the deep neural network. For modeling non-linearities, we apply an exponential linear (ELU) activation function on the output of each node.

The network is trained to find the best match between the partons of the fundamental hard scattering process and the final state particles as measured in a collider detector. In Figure 2(b) we exemplarily show input variables and illustrate the training principle. The partons $b_{\rm lep},b_{\rm had},...$ and the high level variables calculated from them define the input layer of the network. A correct jet-parton assignment is represented by ${\rm jet}_{1}$ being the $b$-jet truly originating from the bottom quark $b_{\rm lep}$ of the leptonic top quark decay, ${\rm jet}_{2}$ truly originating from the bottom quark $b_{\rm had}$ of the hadronic top quark decay, etc. For the incorrect jet-parton assignments, the order of the jets on input is changed accordingly.

The cost function to be minimized in the training process is described by the binomial cross-entropy:

[TABLE]

Here, $D_{\rm DNN}$ is the output of the neural network given the parameters $W$ and $b$ for input observables $x$, which contain all information on the final state particles. $\hat{D}$ is the training target whose value is $1$ for correct jet-parton assignments and [math] otherwise. As the number of incorrect permutations by far exceeds the number of correct ones, incorrect permutations are scaled by $\omega$ to give the same weight in the training. Below we will apply a preselection of at most $50$ permutations such that the weights for incorrect permutations vary within $\omega\in[1/12,1/50]$.

During training, the weights $W$ of the network and the biases $b$ are adjusted such that the network returns the result of a logistic function as the discriminator $D_{\rm DNN}\,\epsilon\,[0,1]$ which scores a permutation given the input $x$. For all possible assignments within one event, the trained discriminator of a perfectly trained network should give the maximum discriminator $D_{\rm DNN}=D_{max}$ for the correct jet-parton assignment, and a smaller value $D_{\rm DNN}<D_{max}$ for all other permutations. A set of $10^{4}$ permutations $(x_{i},\hat{D}_{i})$ forms a training batch, randomly selected from all simulated collision events.

In order to ensure generalization and reduce prediction errors when applied to a new, independent dataset, we introduce regularization measures. We use $L^{2}$ regularization (last term in eq.(2.2)) which penalizes large weights during training and prevents the network from learning features unique to the training set.

We initialize the weights $W$ according to a Gaussian distribution with squared width $\sigma^{2}=2/(l_{in}+l_{out})$ where $l_{in}$ and $l_{out}$ are the number of ingoing and outgoing connections of a layer [16]. We also pre-process the input observables using feature scaling where we normalize $x$ to mean [math] and standard deviation $1$. For the optimization we use the concept of the stochastic gradient descent with an adaptive learning rate (ADAM) [17].

We use the TensorFlow software package to set up the network [18]. The network training is performed on a cluster with GPU cards of GeForce GTX 1080 type with up to 64 GB of system memory. A typical training time amounts to $6$ h.

3 Simulated dataset

To simulate $t\bar{t}H$ events we use the Pythia $8.2.19$ program package [19]. The matrix elements for $t\bar{t}H$ processes include angular correlations of the decay products from heavy resonances. The beam conditions correspond to LHC proton-proton collisions at $\sqrt{s}=13$ TeV. We use only the dominant gluon-gluon process and the Higgs boson decay into a bottom quark pair. Hadronization is performed with the Lund string fragmentation concept.

In order to analyze a typical final state as observed in a LHC detector we use the DELPHES package to simulate the CMS detector [20]. The DELPHES project provides a modular framework that simulates a multipurpose detector in a parameterized way. It includes all major effects such as pile-up, charged particle deflections in magnetic fields, electromagnetic and hadronic calorimetry, and muon detection systems. The simulated output consists of muons, electrons, photons, jets and missing transverse momentum from a particle flow algorithm as well as identifiers for bottom jets and tau leptons.

For the jet finding the anti-$k_{t}$ algorithm is used with the size parameter $R=0.5$. For the $b$-tag algorithm both the efficiency $\varepsilon_{b}$ of correct $b$-tags and the probability $\rho_{b}$ of jets from light quarks being incorrectly $b$-tagged have a dependency on the jet transverse momentum $p^{jet}_{t}$ with the values $\varepsilon_{b}=0.71$ and $\rho_{b}=0.014$ at $p^{jet}_{t}=100$ GeV.

In our analysis of the simulated data, we select events with lepton transverse momenta above $p^{l}_{t}=20$ GeV (electrons or muons) and their pseudorapidities within $|\eta^{l}|<2.1$. The jet transverse momenta are above $p^{jet}_{t}=25$ GeV, with their pseudorapidities being within $|\eta^{jet}|<2.5$. We require the number of jets to be within $6\leq n^{jet}\leq 10$; of these at least $n^{jet}=2$ must be $b$-tagged.

The correct assignment is obtained by matching partons with jets using generator information. A match must fulfill $\Delta R_{jet,parton}<0.3$ with

[TABLE]

where $\phi$ is the azimuthal angle in the detector and $\eta$ denotes the pseudorapidity. Ambiguities are resolved by minimizing the sum $\sum{\Delta R}$ of all viable matches. In order to define valid training targets, we consider only events for which all partons could be matched.

From the total of $10^{8}$ generated events, we could identify all partons of the hard scattering process in $44\%$ of the events. After applying the selection criteria for reconstructed jets and leptons $5.6\%$ of the generated events remained. The additional requirement of the above-mentioned jet-parton match was then fulfilled by $0.7\%$ of the events. Therefore, we perform the investigations with $738,270$ events, using $40\%$ for training purposes and $10\%$ for immediate validation. To determine the correct assignment efficiency we use $50\%$ of the events.

4 Detector observables

We use two categories of variables as the input observables for solving the jet-parton assignment in $t\bar{t}H$ events:

Basic particle-related variables:

1. (a)

For all jets and leptons: four-vector values such as the momenta in both spherical and Cartesian coordinates. 2. (b)

For the charged lepton: isolation variables and sum of transverse momenta $\Sigma p_{T}$ of charged and neutral particles within $\Delta R\leq 0.5$ relative to the lepton momentum. 3. (c)

For jets: b- and $\tau$-tag values, jet area, the number of constituents, and electromagnetic and hadronic calorimeter contributions. 2. 2.

Observables derived from combinations of several objects:

1. (a)

The masses of the Higgs boson, the top quarks, the hadronically decaying $W$ boson, and the individual $\chi^{2}$ values compared to the expected masses. Also the combined mass of the $t\bar{t}$ system. 2. (b)

The transverse momenta of the $t\bar{t}$ system and the $t\bar{t}H$ system, and the ratio of the jet scalar transverse momentum sum to that of all final state particles including the neutrino. 3. (c)

The distance $\Delta R$ between the two jets originating from the Higgs boson decay, and their differences in pseudorapidity and in azimuthal angles. The distance $\Delta R$ between the two jets originating from the $W$ boson decay, and their differences in pseudorapidity and in azimuthal angles. The difference between the pseudorapidities of the lepton and the $b$-jet of the leptonically decaying top quark. 4. (d)

The spatial angle $\theta^{*}$ between the charged lepton in the $W$ boson rest frame and the $W$ boson direction when boosted into the rest frame of its corresponding top quark. The spatial angle $\theta^{*}$ between the softer jet of the $W$ boson decay in the $W$ boson rest frame and the $W$ boson direction in the rest frame of its corresponding top quark.

The neutrino was reconstructed using the charged lepton, and the missing transverse energy according to  [21].

As an example of how different the distributions are when selecting the correct jet-parton assignment or permutations, respectively, in Figure 3(a) we show the transverse momentum distribution of the $p_{t}$ leading light quark jet resulting from the $W$ decay. The full red curve represents the distribution from correct assignments of the jet to a quark from the $W$ decay. The dashed blue curve illustrates the distribution resulting from incorrectly assigning a jet to the $W$ decay. The combinatorics of wrong assignments by far exceeds the single true assignment solution.

In Figure 3 we show similar distributions for b) the pseudorapidity of the jet induced by the bottom anti-quark from the Higgs boson decay, c) the distance $\Delta R$ between the jets originating from the Higgs boson, d) the azimuthal angular distance of the two jets from the $W$ boson decay, e) the reconstructed mass of the Higgs boson, and f) the cosine of the spatial angle $\theta^{*}$ of the leptonic branch as defined in item $2(d)$ above.

In total, we use $140$ basic kinematic observables of the above first category, and $22$ combined observables of the second category.

5 Jet-parton assignments

The efficiencies of finding the correct jet-parton assignment using three different methods are described in the following. As a benchmark, drawing the correct jet-parton assignments by chance is largest for events with $6$ jets where the four $b$-jets were correctly tagged: $1/12\approx 8\%$. As the experimental $b$-tagging efficiency is $\varepsilon<1$ and the probability of successfully tagging $k=4$ $b$-jets is small, $\varepsilon^{k}\ll 1$, this is a rare case within our simulations. The majority of $t\bar{t}H$ events has more than $6$ jets which results in a negligible rate of correct assignments by chance.

5.1 $\chi^{2}$ method

Initially, we use the reconstructed masses of the Higgs boson, both top quarks, and the $W$ boson of the hadronically decaying top quark in the context of a $\chi^{2}$ measure:

[TABLE]

Here, the index $H$ refers to the Higgs boson, $({t,lep})$ to the top quark of the leptonic branch, and $({t,had})$ and $({W,had})$ to the top quark and $W$ boson of the hadronic branch. The tilde symbol indicates the reconstructed average masses $\tilde{m}$ using the correct jet-parton assignments, and the corresponding widths $\tilde{\sigma}$ of the mass distributions.

In Figure 4(a) we show the distribution of our discriminator $D_{\chi^{2}}\equiv\exp{(-\chi^{2})}$ for permutations of an exemplary event. We use this definition of $D_{\chi^{2}}$ for consistency with the subsequent methods, which return better results at higher scores. Correspondingly, the permutation with the highest $D_{\chi^{2}}$ value (smallest $\chi^{2}$) is considered the selected jet-parton assignment and eventually denotes the full event reconstruction. The correct jet-parton assignment is typically contained within the $50$ permutations with the largest $D_{\chi^{2}}$ discriminator values. Therefore we investigate only these $50$ permutations.

In the same figure, we show the discriminator $D_{\chi^{2}}^{truth}$ resulting from the correct jet-parton assignment as a vertical line which coincides in this event with the reconstructed assignment. Also shown in Figure 4(b) are the $D_{\chi^{2}}$ distributions for all events with correct assignment (full red curve) and incorrect assignments (dashed blue curve).

Using the number of events with correct jet-parton assignments in comparison to the total number of events, we find an efficiency of $\varepsilon(t\bar{t}H)=37\%$ to correctly assign all jets to the $t\bar{t}H$ parton final state. If only the correct reconstruction of the Higgs decay is of interest, the efficiency of the correctly assigned $b$-jets is greater and amounts to $\varepsilon(H)=52\%$ (Table 2).

In order to obtain more information on the quality of the selected jet-parton assignment, per event we count the permutations with $D_{\chi^{2}}$ value larger than that of the correct assignment $D_{\chi^{2}}^{truth}$ (see Figure 4(a)):

[TABLE]

The distribution of $I(D_{\chi^{2}})$ can be used to assess the quality of a discriminator. For correctly selected assignments, $I(D_{\chi^{2}})$ vanishes (full curve in Figure 4(c)). However, in the case of an incorrectly selected assignment, $I(D_{\chi^{2}})>0$, it should be as small as possible for good discriminators. For the $\chi^{2}$ method, the resulting distribution of $I(D_{\chi^{2}})$ is shown for all events by the dashed blue curve in Figure 4(c).

5.2 Boosted decision trees

As an alternative method we use boosted decision trees (BDT) as implemented in [22]. We apply the AdaBoost algorithm for boosting. With respect to the default hyperparameters we found improved results when using an increased number of $400$ decision trees, but no substantial change when varying the maximum depth of a tree, or the learning rate of the boosting algorithm. The discriminator $D_{\rm BDT}$ of the BDT is the weighted average of all decision trees each classifying a given jet-parton assignment as signal or background alternatively.

For the BDT we use the same observables as described in section 4. We apply the same training principle as described in section 2 and illustrated in Figure 2(b). For training the BDT, a signal indicator is given in the case that the order of the observables follows the correct jet-parton assignment. Correspondingly, a background flag is given for other permutations.

On evaluation of the jet-parton assignment in a reconstructed event, first a preselection is performed using the above-mentioned $\chi^{2}$ method. Only the $50$ permutations with the largest $D_{\chi^{2}}$ discriminator (smallest $\chi^{2}$) values are considered for evaluation by the BDT method. The permutation with the largest discriminator $D_{\rm BDT}$ is considered as the selected jet-parton assignment. The efficiency of finding the correct parton-jet assignment is $\varepsilon(t\bar{t}H))=45\%$ which is better compared to the $\chi^{2}$ method. Correspondingly, also the efficiency of providing the correct $b$-jets of the Higgs decay is improved: $\varepsilon(H)=57\%$ (see Table 2).

As for the $\chi^{2}$ method, we assess the performance of the BDT discriminator $D_{\rm BDT}$ by examining the number of permutations yielding a higher value than the correct one, $I(D_{\rm BDT})$ (5.2). This distribution is shown in Figure 4(c) by the dashed green curve. The efficiency of the BDT method is better than that of the $\chi^{2}$ method as can be seen at $I(D_{\rm BDT})=0$. Also, the distribution flattens out faster for $I(D_{\rm BDT})>0$ compared to the $\chi^{2}$ method, indicating a better overall performance of the BDT method.

5.3 Neural network

In our third method, we investigate the performance of the neural network. Network architecture, training and evaluation procedures are described in section 2, and the observables in section 4.

For evaluation of the jet-parton assignment in a reconstructed event, we again preselect only the $50$ permutations with the largest $D_{\chi^{2}}$ values of the $\chi^{2}$ method. On these permutations, the network delivers a corresponding discriminator $D_{\rm DNN}\,\epsilon\,[0,1]$ where the maximum discriminator value leads to the selected jet-parton assignment. In Figure 4(d) the distributions of the discriminator $D_{\rm DNN}$ are shown for all events with correct assignment (full red curve) and incorrect assignments (dashed blue curve).

For the majority of the events, the neural network provides the correct jet-parton assignment. The corresponding efficiency is $\varepsilon(t\bar{t}H))=52\%$, which is better than the BDT method (see Table 2). In $\varepsilon(H)=63\%$ of the events, the network provides the correct jets of the Higgs boson decay.

We also compare the performance of the network by counting permutations with the discriminator $D_{\rm DNN}$ value exceeding that of the correct jet-parton assignment. The resulting distribution of $I(D_{\rm DNN})$ (5.2) is shown in Figure 4(c) for all events by the dashed red curve. Again, the high efficiency of finding the correct jet-parton assignments is visible at $I(D_{\rm DNN})=0$. The distribution of the non-zero counts $I(D_{\rm DNN})$ is consistently below the distributions of the two other methods.

6 Impact on $t\bar{t}H$ analyses

Finding the correct jet-parton assignment in $t\bar{t}H$ events can have a major impact on the sensitivity of a corresponding cross-section analysis. As a quantitative measure of the improvement in analysis sensitivity is beyond the scope of this work, we will show distributions of key observables typically used in such analyses.

In Figure 5(a) we show the distribution of the distance $\Delta R$ of the jets originating from the Higgs boson decay for the correct jet-parton assignment (black curve) in comparison to the three reconstruction methods investigated above. The result of the neural network is shown by the dashed red curve and consistently depicts the best distribution of correct jet-parton assignments. The BDT method is illustrated by the dashed green curve which on average yields $\Delta R$ values that are too large. The $\chi^{2}$ method (dashed blue curve) on average returns $\Delta R$ values that are too small.

Figure 5(b) shows distributions of the azimuthal angular distance of the two jets of the $W$ boson decay with a similar quality in the description for the BDT and neural network methods but worse for the $\chi^{2}$ method. The reconstructed mass of the Higgs boson can be seen in Figure 5(c). Both the $\chi^{2}$ and BDT methods have a bias towards the center of the distribution while the best possible reconstruction using the correct jet-parton assignment is less pronounced. The latter distribution is best described by the neural network.

In Figure 5(d) we also show the angular characteristics of the leptonic $W$ decay. This is encoded in the distribution of the angle $\theta^{*}$ between the leptonically decaying $W$ in the rest frame of its top quark, and the charged lepton which is boosted into the $W$ rest frame. The neural network yields the best description among the three tested methods.

Overall, the neural network provides the best description of the distributions obtained with the correct jet-parton assignment (truth).

7 Conclusions

In this work, we investigated the jet-parton assigments in simulated $t\bar{t}H$ events using deep learning techniques. With a high fraction of correct assignments, high-level variables such as reconstructed masses or decay angular characteristics can be determined more precisely, which is advantageous when separating signal from background processes.

Our study was based on simulated $t\bar{t}H$ events generated with a parameterized detector simulation. Our investigations were carried out using events which contain all information for reconstructing the $t\bar{t}H$ parton final state.

The fully connected architecture of the neural network contained 8 hidden layers with 500 nodes each and additional connections following the approach of residual networks. For training, the order of the input variables was according to the correct jet-parton assignment (signal) or a permutation of that order (background), respectively.

During evaluation, a discriminator was calculated for each permutation of the jet-parton assignments, and the assignment with the largest discriminator value was considered the selected assignment which is equivalent to a full event reconstruction. When compared to two other commonly used methods, the neural network approach returned the best results by far, correctly reconstructing the entire $t\bar{t}H$ event in $52\%$, and only the Higgs boson in $63\%$ of the cases.

Acknowledgments

This work is supported by the Ministry of Innovation, Science and Research of the State of North Rhine-Westphalia, and the Federal Ministry of Education and Research (BMBF). We wish to thank David Walz for his valuable comments on the manuscript.