# Bayesian Probabilistic Modelling for Four-Tops at the LHC

Ezequiel Alvarez<EMAIL_ADDRESS>International Center for Advanced Studies
(ICAS) and CONICET, UNSAM, Campus Miguelete, 25 de Mayo y Francia, CP1650, San
Martin, Buenos Aires, Argentina  Barry M. Dillon<EMAIL_ADDRESS>heidelberg.de Institut für Theoretische Physik, Universität Heidelberg,
Germany  Darius A. Faroughy<EMAIL_ADDRESS>Physik-Institut,
Universität Zürich, CH-8057, Switzerland  Jernej F. Kamenik
<EMAIL_ADDRESS>Jožef Stefan Institute, Jamova 39, 1000 Ljubljana,
Slovenia  Faculty of Mathematics and Physics, University of Ljubljana,
Jadranska 19, 1000 Ljubljana, Slovenia  Federico Lamagna
<EMAIL_ADDRESS>Centro Atómico Bariloche, Instituto Balseiro
and CONICET  Manuel Szewc<EMAIL_ADDRESS>International Center for
Advanced Studies (ICAS) and CONICET, UNSAM, Campus Miguelete, 25 de Mayo y
Francia, CP1650, San Martin, Buenos Aires, Argentina  Jožef Stefan Institute,
Jamova 39, 1000 Ljubljana, Slovenia

###### Abstract

Monte Carlo (MC) generators are crucial for analyzing data in particle
collider experiments. However, often even a small mismatch between the MC
simulations and the measurements can undermine the interpretation of the
results. This is particularly important in the context of LHC searches for
rare physics processes within and beyond the standard model (SM). One of the
ultimate rare processes in the SM currently being explored at the LHC, $pp\to
t\bar{t}t\bar{t}$ with its large multi-dimensional phase-space is an ideal
testing ground to explore new ways to reduce the impact of potential MC
mismodelling on experimental results. We propose a novel statistical method
capable of disentangling the 4-top signal from the dominant backgrounds in the
same-sign dilepton channel, while simultaneously correcting for possible MC
imperfections in modelling of the most relevant discriminating observables –
the jet multiplicity distributions. A Bayesian mixture of multinomials is used
to model the light-jet and $b$-jet multiplicities under the assumption of
their conditional independence. The signal and background distributions
generated from a deliberately mistuned MC simulator are used as model priors.
The posterior distributions, as well as the signal and background fractions,
are then learned from the data using Bayesian inference. We demonstrate that
our method can mitigate the effects of large MC mismodellings in the context
of a realistic $t\bar{t}t\bar{t}$ search, leading to corrected posterior
distributions that better approximate the underlying truth-level spectra.

##  I Introduction

In recent years, the large abundance of LHC data on one hand, and the absence
of clear New Physics (NP) signals in theory driven analyses of this data on
the other, have motivated the development of novel, more data driven
approaches to LHC data analysis and NP searches. In particular, the advent of
unsupervised and weakly-supervised Machine Learning (ML) techniques has
allowed for the development of broad model independent NP search and
characterisation strategies [1]. Simultaneously, there have been important
efforts to reduce reliance of LHC measurements on Monte Carlo (MC) simulations
of hadronic processes [2, 3, 4, 5, 6].

The simultaneous production of four top quarks represents an important NP
benchmark (see e.g. Refs. [7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18]), but
also an interesting point of coalescence for several of these developments
[19]. One of the main issues in studying this final state is its tiny cross-
section (12 fb) compared to its main backgrounds ($\sim 600$ fb), which is
compounded by the challenges to correctly model the complex final states
through MC simulations. To address these issues, we have previously studied
the two lepton same sign channel (2LSS${\pm\pm}$) [20] which in the SM may
contain signal and background events up to the same order of magnitude and
furthermore exhibits somewhat reduced complexity of the (multi jet) final
state, compared to the single lepton channel [21, 22]. In the 2LSS++ channel
$t\bar{t}W^{+}$ production represents the main and most challenging background
for the 4-top signal.111Our results and discussion would apply equally well to
other non-negligible backgrounds such as $t\bar{t}h$ and $t\bar{t}Z$. Recent
experimental analyses in this channel [23, 24] have highlighted difficulties
in reliably modelling the signal and background kinematics using state of the
art MC tools. This in turn hinders the sensitivity of this important signature
to possible NP effects in four-top production.

Using the experimental challenge described above as an example and motivation,
in the present paper we describe a novel Bayesian statistical framework to
disentangle in-situ signal and background distributions of categorical data.
Our method can be used to simultaneously identify and correct potential (MC)
mismodelling of discrete distributions as well as extract signal and
background admixtures in the data close to their truth values.

The paper is organized as follows. In Sec. II we introduce our statistical
model of multinomial mixtures with Bayesian inference and demonstrate its use
on a toy example. We apply the model to jet multiplicity distributions in the
2LSS++ channel of 4-top production at the LHC in Sec. III and show how it can
be used to identify and correct MC mismodelling and extract signal and
background fractions. Sec. IV is devoted to a detailed study of the
assumptions and consistency checks of the model when applied to realistic
datasets. Finally, we summarize our findings in Sec. V.

##  II Categorical mixture model for four-tops

Anticipating the application to 4-top production, in the following we
represent an event generation process by a pair of random variables
$(N_{j},N_{b})$ indicating the number of clustered light-jets and $b$-jets,
respectively. Our starting point is that a collection of such events can be
described using a likelihood with a joint probability density $p(j,b)$ where
$j$ ($b$) are the observed number of light-jets (b-jets) in an event. The most
general discrete model for this likelihood is the multinomial
distribution222Along this work we refer to multinomial distribution although
in all cases it consists of a single drawing per event and therefore it is
also a categorical distribution, which is a special case of the former. with
$d_{j}\times d_{b}-1$ parameters, where $d_{j,b}$ are the number of possible
light-jets and $b$-jets to be expected in an event. However, our goal is to
disentangle the contributions to this joint likelihood arising from four-top
events and $t\bar{t}W$ events. To do so we introduce two mixture components,
one for $t\bar{t}W$ and one for four-top. If we simply describe each mixture
with a multinomial distribution $p(j,b|z)$ with $z\in[0,1]$ representing the
mixture label, we would have a mixture model with $2\times(d_{j}\times
d_{b}-1)+1$ parameters. Since each event is independent and consists of just a
single draw from this distribution, each mixture can describe all possible
combinations of $N_{j}$ and $N_{b}$ values in the data and therefore all
correlations by itself. The model would thus over-parameterize the data making
the inclusion of mixtures redundant.333Note that this would not be the case if
each event was generated by several draws from $p(j,b|z)$, since there would
then be additional correlations between the multiple draws per event. This is
the case in s.c. mixed membership models [25, 26, 27] used in jet substructure
analyses where the mixtures describe correlations between the multiple draws
per event.

Therefore the key insight is to instead write down a mixture model in terms of
$p(j|z)$ and $p(b|z)$, such that the correlations between $N_{j}$ and $N_{b}$
in the dataset are parameterized by the class label alone. The number of
parameters in this model is $2\times(d_{j}+d_{b}-2)+1$. To be explicit, we
optimize the model to parameterize the correlations between $N_{j}$ and
$N_{b}$ in terms of a discrete variable $Z$, and interpret this as a class
label for four-top and $t\bar{t}W$ events. We are making the simplifying
assumption that $N_{j}$ and $N_{b}$ are conditionally independent variables,
that all correlations between them in the dataset arise only from assignments
to the two classes. Conditional independence is of course an approximation. In
particular, in a realistic measurement setting, $N_{j}$ and $N_{b}$ are not
strictly conditionally independent due to mis-tagging or other reconstruction
imperfections. The degree to which the method succeeds is limited by this
approximation. Conversely, a failure of the method to converge to a consistent
description of the measured distributions would be a clear sign that the
assumptions of the statistical model are not respected by the dataset. We
return to this important caveat and discuss its mitigation in Sec. IV.444A
systematic study of statistical models which go beyond strict conditional
independence assumptions is in progress and will be presented elsewhere.
However, as we will show, in the case at hand, the method exhibits good
convergence indicating that conditional independence holds sufficiently well
in practice.

Within the limitations described above, the generative process for the dataset
proceeds as follows: for each event ($n$) a class label $z_{n}$ is first drawn
from a binomial probability distribution parametrized by $\pi\in[0,1]$. Then
$j_{n}$ and $b_{n}$ are sampled from separate multinomials corresponding to
the drawn class and parametrized by $\alpha_{z,i}$ and $\beta_{z,k}$,
respectively, where $i$ and $k$ run up to $d_{j}$ and $d_{b}$, respectively.
We assume that the whole dataset $X$, consisting of $n\in N$ pairs of
measurements $x_{n}$ = ($j_{n}$, $b_{n}$) for the 2LSS++ selected events, is
generated through this probabilistic model and we want to infer the values of
its parameters, namely $\pi,\alpha_{0,j},\beta_{0,i},\alpha_{1,j}$ and
$\beta_{1,i}$, which we collectively indicate as $\theta$. Observe that the
described model corresponds to a special case of a mixture of multinomials
[28].

Adopting a Bayesian framework, we consider the model parameters ($\theta$) to
be random variables as well and we want to update our knowledge of these
random variables after measuring $X$. However, it is more convenient in
practice to consider explicitly also the latent variables $Z$ which represent
the class assignments of each event. Graphically, the probabilistic model can
be represented through the plate diagram in Fig. 1 and leads to the posterior:

$p(Z,\pi,\alpha,\beta|X)=\frac{p(X,Z,\pi,\alpha,\beta)}{p(X)}\,,$ (1)

where the joint distribution $p(X,Z,\pi,\alpha,\beta)$ is given explicitly by

$\displaystyle p(X,Z,\pi,\alpha,\beta)$ $\displaystyle=$
$\displaystyle\prod_{n=1}^{N}p(x_{n}|z_{n},\alpha,\beta)p(z_{n}|\pi)$
$\displaystyle
p(\pi|\eta_{\pi})\prod_{k=0}^{1}p(\alpha_{k}|\eta_{\alpha_{k}})p(\beta_{k}|\eta_{\beta_{k}})\,.$

Here $p(x_{n}|z_{n},\alpha,\beta)=\alpha_{z_{n}j_{n}}\beta_{z_{n}b_{n}}$,
$p(z_{n}|\pi)=\pi_{z_{n}}$ and $p(\pi|\eta^{\pi})$,
$p(\alpha_{k}|\eta^{\alpha_{k}})$ and $p(\beta_{k}|\eta^{\beta_{k}})$ are
Dirichlet distributions with the corresponding $\eta^{i}$ set of parameters.

Figure 1: Plate diagram of (Bayesian) 2-mixture model of multinomials for
$(N_{b},N_{j})$ $N$-event dataset. From the Dirichlet prior distributions
(with hyperparameters $\eta_{i}$) the multinomial parameters ($\pi,\
\alpha_{i}$ and $\beta_{i}$) are sampled, then $N$ events are sampled through
a latent variable $Z$ that determines in turn from which multinomial the two
observables in each event ($N_{j}$ and $N_{b}$) are sampled.

The main idea in this expression is that given the dataset $X$, a
probabilistic model that allows us to write down an expression for
$p(X|\theta)$ and a reasonable prior $p(\theta)$, we can in principle
determine the probability density function (pdf) for the parameters
$p(\theta|X)$. This is a powerful result, since it gives us not only the
fraction of signal to background and its uncertainty through $p(\pi|X)$
marginalizing over the other parameters, but it can also give us the $N_{j}$
and $N_{b}$ distributions of both individual classes. If the probabilistic
model describes the data well and the prior is reasonable, then these should
match within uncertainties the true underlying background and signal $N_{j}$
and $N_{b}$ distributions.

There are many known approaches to solving Eq. (1) using Bayesian Inference;
including mean-field techniques such as Variational Inference (VI) [28] and
numerical Markov Chain Monte Carlo methods such as Gibbs Sampling (GS) [28].
Below we focus on the latter numerical approach which turns out to be
preferred to the mean-field methods which approximate the posterior with a
fully factorized model that neglects possible correlations between the
inferred parameters. As we are interested in finding the correlations between
$N_{j}$ and $N_{b}$ through class assignment, VI is challenged by definition
to find the appropriate correlations.

The goal of the GS algorithm is to approximate the posterior through the use
of a finite number of samples. These samples can then be used to obtain any
desired expected values such as the mean of the relevant parameters
$\mathbb{E}[\theta_{i}]$. To obtain samples from the posterior, each iteration
samples an observation of each parameter $\theta_{i}$ from the marginal
distribution conditioned on the remaining parameters
$p(\theta_{i}|\theta_{\setminus i},X)$. When implementing a Gibbs sampler to
approximate Eq. (1), the conditional distributions can be obtained and sampled
from efficiently, being either Dirichlet or Multinomial distributions. Our
algorithm implemented in python is available at GitHub [29].

In practice, subsequently drawn samples are highly correlated. To mitigate
this we drop the first $M$ samples, which constitute what is called the burn-
in phase, and then apply a ‘thinning’ procedure which consists in only keeping
every $l^{\text{th}}$ sample. We also implement different chains, or walkers,
initialized at different randomly chosen starting points. We estimate
sufficient $M$ and $l$ values by computing the integrated autocorrelation time
$\tau$ as defined in Ref. [30] and adapting its implementation in emcee [31]
accounting for the fact that we do not have an ensemble sampler. We find that
with $30$ walkers and $1000$ saved iterations per walker after thinning with
$l=100$ we have $\tau$’s in the range $\tau\in[1,2.5]$. We consider a burn-in
phase of $M=1000$ after which we save the aforementioned $1000$ samples with
thinning.555The GS algorithm and techniques described above are well known in
other disciplines, in particular computer sciences, however they have to our
knowledge not been applied before in the context of (high energy) physics.
Once we have an accurate approximation of $p(Z,\pi,\alpha,\beta|X)$, we can
marginalize over the class assignments by neglecting the sampled values $Z$.

###  II.1 A simple toy example

To demonstrate the efficiency of this approach, as well as the limitations due
to the approximations we make, we will first apply it to inference in a very
simple toy model. We take a sample of ‘events’, each with just two features.
The sample is comprised of two types of events, which for the sake of analogy
we call background and signal. These signal and background events are sampled
from sets of overlapping distributions in the feature-space. The features for
each event are sampled independently, therefore in this simple toy example
these two features are completely uncorrelated from each other. We consider
the case in which the prior distributions for these features are not too far
from the truth. In contrast, we consider a uniform prior distribution for the
$\pi$ parameter giving the fraction of signal and background in the sample.
This indicates no prior knowledge on how much background and signal we can
expect in the dataset and is the most conservative assumption we can make in
this regard. We show the prior distributions as well as the true values of the
parameters in the upper row of Fig. 2

After numerically solving the Bayes Inference problem using GS, we compare the
class-$0$ and class-$1$ inferred distributions for $N_{j}$ and $N_{b}$ to the
truth-level background and signal distributions in $X$. A good summary to
assess the success of the algorithm is the corner-plot which visualizes the
distribution through marginalizing to either two or one parameter dimensions
and the true values. An excerpt is shown in Fig. 3. In each panel we show the
corresponding prior distribution (red), posterior distribution (black) and the
true values (blue). Quantitatively, one can also compare the level of
improvement between the prior and the posterior by computing their Log-
Likelihood Ratio (LLR) with respect to the true value for each parameter. We
display these numbers above the diagonal panels of the corner-plot, and we see
a robust improvement in most of them. To compute the LLR of the posterior and
prior of the complete model one would in principle need to evaluate the joint
density distributions of all pairs of parameters (off-diagonal elements in the
corner-plot) which is beyond the scope of this work. Instead, as a rough
approximation, neglecting the correlations between the parameters, we obtain a
global LLR as a sum of the individual parameter LLRs, LLR $\approx$ 36\. We
display this global sum as well as partial sums grouping different parameters
together in Fig. 3. We also include the partial and global sums of LLR
obtained when approximating the posterior through VI. We observe that although
VI captures the maximum of the posterior accurately, it consistently
underestimates the variance of the distribution yielding a too narrow
approximation to the GS obtained posterior. This is reflected in a lower
improvement over the prior (LLR $\approx$ 15).

Finally, in the bottom row of Fig. 2 we group together the one-dimensional
marginalized posterior distributions for each parameter to obtain the $N_{j}$,
$N_{b}$ distributions of the signal and background, as well as for the $\pi$
parameter, i.e. the fraction of signal in the sample. In the plot the true
value of the parameters is shown in solid blue. Notice that the posterior
exhibits good convergence to the true values as well as a considerable
reduction of the uncertainty, when compared to the prior, which emulates the
imperfect MC. We find an improvement in both the $N_{j}$ and $N_{b}$
distributions for each process as expected from Fig. 3. It is also interesting
to notice in Fig. 2 how from a complete ignorance of the signal and background
fractions in the sample, the algorithm recovers a pdf for $\pi$ in good
agreement with its true value.

Figure 2: $N_{j}$, $N_{b}$ and $\pi_{1}$ distributions: true values (blue),
priors (red) and posterior (black) for the toy model. Shaded regions in first
four plots indicate the $1\sigma$ uncertainty region. Comparing the posteriors
to the priors one can appreciate the improvement in estimating the true
distributions departing from incorrect and uncertain priors using Bayesian
Inference on the data.  Figure 3: Excerpt from the corner-plot for the toy
model. Red indicates the prior distribution, black the posterior distribution
obtained through GS and blue is the true value. We see how the posterior
distribution captures the correlation between $N_{j}$ and $N_{b}$. The titles
of each 1D histogram contain the Log-Likelihood Ratio between the posterior
and the prior using either GS or VI for the posterior estimation, with the
latter shown in parentheses. The table contains the sum of Log-Likelihood
Ratios per parameter block, again considering the posteriors obtained through
GS and through VI. We see that VI is a bad approximation to GS, failing to
improve on the prior for several parameter blocks.

##  III Application to four-top measurements

In the 2LSS++ channel, the final state is usually characterized by at least
$2\ell^{+}$, at least 2 $b$-tagged jets, and at least 4 light jets. Additional
cuts on missing transverse energy and transverse momentum may be invoked to
enhance the signal fraction in the sample. The exact details of the event
selection are however not important for the purposes of this work. From the
decay products at matrix-element level of the signal, one expects a priori
that the $N_{j}$ and $N_{b}$ distributions to be skewed towards higher values
when compared to the background process, thus providing enough separation for
disentangling them using statistical inference.

In our setup we have simulated 4-top and $t\bar{t}W^{\pm}$ events using
Madgraph [32], Pythia [33] and Delphes [34] to account for matrix level
calculations and showering, hadronization and detector simulation,
respectively. We selected $N=500$ events, roughly equivalent to
$\mathcal{L}=800\mathrm{fb}^{-1}$, in the 2LSS++ channel with 70% background
and 30% signal (we also tested for other signal fractions and obtained similar
results). Using this data we created a dataset $X$, represented by $N$ pairs
($j_{n},b_{n}$), $n=1,\ldots,N$, to serve as our benchmark truth-level sample.
The resulting two dimensional distributions are shown in Fig. 4.

Figure 4: $N_{j}$ and $N_{b}$ two dimensional distributions for
$t\bar{t}W^{\pm}$ and 4-top.

We observe from Fig. 4 that the $N_{j}$ and $N_{b}$ distributions do not
appear to be strictly (conditionally) independent. This is evidenced by the
fact that different rows (columns) show different bin hierarchies depending on
the column (row) they are conditioned on. These effects arise from
experimental systematics such as imperfect b-tagging and different (b-)jet
acceptances, as well as from statistical fluctuations due to finite sample
sizes involved: the sample size of the Monte Carlo simulation and the sample
size of the expected events at to the collider luminosity considered. Category
bins with very small event yields are particularly affected by these later
effects. In Sec IV we study in detail how well our model approximates the true
data distributions even when conditional independence is not exact. We find
that for the foreseen (HL)LHC luminosities our model is statistically
indistinguishable from the data while retaining classification power to infer
the 4-top and $t\overline{t}W^{\pm}$ distributions.

On the other hand, regarding the potential MC mismodelling, we would like to
emphasize that our model is aimed to work directly on data and thus address
this very kind of problem. That is, we care that our model recovers the true
underlying distribution with imperfect (i.e. MC based) priors. In this context
we use MC simulations as stand-in mock data for actual (mixed) distribution
measurements and apply our model to this mock data with imperfect knowledge
encoded in the priors. In order to emulate an imperfect MC prior we skewed the
corresponding $N_{j}$ and $N_{b}$ distributions from $X$ to higher values and
incorporated this into our model through the prior hyperparameters. In
general, we can write the hyperparameters $\eta$ of a $V$-dimensional
Dirichlet distribution of a random variable $\theta$ as $\eta_{v}=\Sigma\cdot
p_{v}$, for $v=1,\dots,V$. Here $p$ is a multinomial probability distribution
and $\Sigma$ is a normalization factor. The role of $p_{v}$ and $\Sigma$ can
be understood by looking at the mean and variance of $\theta_{v}$:

$\displaystyle\mathbb{E}[\theta_{v}]$ $\displaystyle=$ $\displaystyle p_{v}$
$\displaystyle\text{Var }[\theta_{v}]$ $\displaystyle=$
$\displaystyle\frac{p_{v}(1-p_{v})}{\Sigma+1}\,.$ (2)

From these equations, we see that $p_{v}$ represents the expected value of
$\theta_{v}$ while $\Sigma$ controls the confidence we have on that
expectation. We fixed the $p_{v}$ values of the priors for $\alpha$ and
$\beta$ in their respective Dirichlets to the normalized $N_{j}$ and $N_{b}$
populations given by the imperfect MC predictions. To reflect our confidence
in this estimate, in this example we chose $\Sigma=10$ for each Dirichlet. See
Fig. 5 upper row, where we plot the central values and $1\sigma$ ranges for
the prior distributions for $\alpha$ and $\beta$. In an actual experimental
analysis, $\Sigma$ could be chosen such that the priors cover all reasonable
ranges of the modeled observables. As an extreme example, for the prior on the
$\pi$ parameter, giving the fraction of signal and background in the sample,
we take a uniform distribution, indicating no prior knowledge on how much
background and signal we can expect in the dataset.

Figure 5: $N_{j}$, $N_{b}$ and $\pi_{1}$ distributions: true values (blue),
priors (red) and posterior (black). Shaded regions in first four plots
indicate the $1\sigma$ uncertainty region. Comparing the posteriors to the
priors one can appreciate the improvement in estimating the true distributions
departing from incorrect and uncertain priors using Bayesian Inference on the
data.

As we do for the toy model, we study the posterior distribution obtained using
GS through the corner-plot, with its LLR partial and global and sums, and
through histograms that condense the class-$0$ and class-$1$ $N_{j}$ and
$N_{b}$ probability distributions and the $\pi$ probability distribution. We
show an excerpt of the corner-plot in Fig. 6. The global sum of the LLRs is
$\approx$ 20, reflecting an improvement over the prior. In comparison, the VI
estimated posterior does not show an improvement over the prior. This is due
to the narrow width of the approximation which excludes the true values of the
parameters to a higher level than the more accurate GS obtained posterior
estimation.

In Fig. 5 we show the results for $N_{j}$ and $N_{b}$ distributions of the
signal and background, as well as for the $\pi$ parameter, i.e. the fraction
of signal in the sample. As in the toy model case, the posterior exhibits good
convergence to the true values as well as a considerable reduction of the
uncertainty when compared to the prior which emulates the imperfect MC.
However, in this case the improvement is different for each feature. The
$N_{j}$ distribution shows a larger improvement, as expected from Fig. 6,
while the $N_{b}$ distribution is harder to reconstruct due to the much larger
fraction of events populating the first bin. Similar results are obtained for
other cases which differ in signal-to-background ratio and number of events.
It is also interesting to notice in Fig. 5 how again from a complete ignorance
of the signal and background fractions in the 2LSS++ sample the algorithm
recovers a pdf for $\pi$ in good agreement with its true value. We also
checked that this agreement holds for other truth values of $\pi$, and that
the matching only worsens as the value of $\pi$ approaches the boundaries of
$[0,1]$.

In summary, we find that the algorithm successfully infers the $N_{j}$ and
$N_{b}$ distributions as well as the signal/background fractions. Notably, the
best inference occurs for the $N_{j}$ distribution, which is usually the
hardest to predict correctly through MC simulations based on perturbative QCD
calculations matched to parton shower algorithms.

Figure 6: Excerpt from the corner-plot. Red indicates the prior distribution,
black the posterior distribution obtained through GS and blue is the true
value. We see how the posterior distribution captures the correlation between
$N_{j}$ and $N_{b}$. The titles of each 1D histogram contain the Log-
Likelihood Ratio between the posterior and the prior using either GS or VI for
the posterior estimation, with the latter shown in parentheses. The table
contains the sum of Log-Likelihood Ratios per parameter block, again
considering the posteriors obtained through GS and through VI. We see that VI
is a bad approximation to GS, failing to improve on the prior for several
parameter blocks.

##  IV Testing model validity

Our method hinges on the validity of the underlying statistical (generative)
model. Thus it is imperative to understand how well our model that assumes
conditional independence, approximates the true data distributions even when
their conditional independence is not exact. To quantify the agreement between
the data and our model we consider the mutual information (MI)
$I(N_{j},N_{b})$ between $N_{j}$ and $N_{b}$,

$\displaystyle I(N_{j},N_{b})$ $\displaystyle=$ $\displaystyle
D_{\text{KL}}(p({j},{b})||p({j})p({b}))$ $\displaystyle=$
$\displaystyle\sum_{{j}=4}^{9}\sum_{{b}=2}^{5}p({j},{b})\text{ Ln
}\frac{p({j},{b})}{p({j})p({b})}\,.$

The MI encodes how much information is lost by approximating the full
distribution with the product of the two marginal distributions. We can also
condition the MI on the class label and obtain the MI for each process
$I(N_{j},N_{b}|z)$. By combining the per process MI, we build the conditional
MI $I(N_{j},N_{b}|Z)=\sum_{z}p(z)I(N_{j},N_{b}|z)$ which encodes our exact
model hypothesis: the data follows a probability distribution which can be
written as a combination of two processes, each of which presents a factorized
probability distribution. We should note that $I(N_{j},N_{b}|z)$ and
$I(N_{j},N_{b}|Z)$ depend explicitly on the availability of labelled data and
thus are not computable purely from measured distributions. However, because
we expect the simulations to be qualitatively reasonable approximations to
measurements, studying the validity of the modelling hypothesis using MC
simulations is justified.

Using our finite 4-top and $t\bar{t}W^{\pm}$ dataset, we can estimate the
relevant probability distributions and obtain finite sample estimations of the
relevant MIs. In the large statistics limit, the estimator follows compact
asymptotic distributions [35]. However, we are dealing with finite event
samples where some category bins are scarcely populated. Thus, in order to
quantify the compatibility of our model with the data, we do a series of
pseudo-experiments according to the following procedure:

  1. 1.

We take the expected event rates obtained from the Madgraph+Pythia+Delphes
pipeline and their uncertainties to generate 2500 pseudo-datasets. For each
pseudo-dataset, we sample the expected event rate for each bin according to a
Gaussian centered in the MC central value and with the appropriate
uncertainty. Then, we sample the observed events for that bin through a
Poisson distribution.

  2. 2.

For each of these pseudo-datasets, we compute the two-dimensional probability
distribution and the marginals for each process and for the full dataset. With
these, we obtain the estimators of all four relevant MIs
$\hat{I}(N_{j},N_{b}|z)$, with $z=t\overline{t}W^{\pm},\text{4-top}$,
$\hat{I}(N_{j},N_{b})$ and $\hat{I}(N_{j},N_{b}|Z)$.

  3. 3.

We use these estimators to study the validity of approximating the joint
probability distribution with a certain modelling hypothesis. To this end we
construct the probability distribution of the estimator by generating another
batch of 2500 pseudo-datasets. This time, each pseudo-dataset is generated
using the relevant approximation: for $I(N_{j},N_{b}|z)$, we generate the
pseudo-datasets with $p({j}|z)p({b}|z)$; for $I(N_{j},N_{b})$, we generate the
pseudo-datasets with $p({j})p({b})$; and for $I(N_{j},N_{b}|Z)$, we generate
the pseudo-datasets with $\sum_{z}p(z)p({j}|z)p({b}|z)$. The hypothesis that
the obtained estimators $\hat{I}$ are sampled according to the model is the
null hypothesis $H_{0}$.

  4. 4.

Having obtained the probability distribution of each estimator conditioned on
its null hypothesis $H_{0}$ using these additional pseudo-datasets, we compute
the one-sided p-value for the ”measured” estimator which allows us to discard
the null hypothesis with a certain confidence level666Although not explicit,
there is an assumed alternative hypothesis $H_{1}$: the saturated model. For a
given pseudo-dataset of $N$ events sampled from a multinomial distribution,
its MI is nothing more than $2N$ times its saturated log-likelihood [36].. The
p-value can be computed as

$\text{p-value}=\int_{\hat{I}}^{\infty}p(I|H_{0})dI$ (4)

where $I$ can be any of the four metrics considered and $H_{0}$ its associated
null hypothesis. In the large statistics limit, this one-sided test
asymptotically converges to the compact formulae considered in Ref. [35].

We show the results of this procedure in Fig. 7 for four types of pseudo-
datasets. In solid black line we show the pseudo-dataset generated with the
expected events as obtained from the Madgraph+Pythia+Delphes pipeline. In
dashed black line we consider the event rates we obtain when considering
perfect b-tagging in Delphes. We do this to verify whether the introduction of
imperfect b-tagging, and the resulting correlations between the number of
light- and b-jets, spoil conditional independence. In green solid and dashed
lines we modify the sampled expected event rates to ensure conditional
independence for realistic and perfect b-tagging. These two pseudo-datasets
thus agree with our modelling hypothesis and provide a self-consistency check.
One should note that the Poisson sampling with relatively small event rates
induces a slight violation of conditional independence as it is done in a bin
by bin basis.

Figure 7: Top row: We show in solid (dashed) black lines the MI between
$N_{j}$ and $N_{b}$ with realistic (perfect) b-tagging. In solid (dashed)
green we show the MI distribution for the expected event rates which respect
conditional independence with realistic (perfect) b-tagging. Bottom row: We
show with the same color and line conventions the p-value of the null
hypothesis distribution of each estimator. We show in red the $p=0.05$
conventional exclusion value. We can see that for the considered luminosity,
$N_{j}$ and $N_{b}$ cannot be ruled out to be conditionally mutually
independent. See text for details.

In Fig. 7, we observe that for the considered luminosity $\mathcal{L}\simeq
800$ fb-1, the data and our model are not statistically distinguishable from
each other. This can be seen from the first, second and fourth columns, where
the null hypothesis coincides with the green curves. The p-value distributions
in the first and second column imply that 4-top and $t\overline{t}W^{\pm}$
cannot be ruled out to have factorized $(N_{j},N_{b})$ distributions while the
fourth column implies the same for the full data and the model which assumes
conditional independence. For the third column, both the data and the model
are different from the null hypothesis that considers full independence
between $N_{j}$ and $N_{b}$. In that case, both the model and the data show
slight disagreements with the null hypothesis although they remain compatible
with it. We observe how the p-value distribution is more tilted towards the
discarded region for the MI compared to the Conditional MI for the full data
distribution, specially for perfect b-tagging. Conditional independence is
thus a reasonable modelling hypothesis that yields qualitatively different
behavior than assuming a single process with a factorized ($N_{j}$,$N_{b}$)
distribution. Because conditional independence assumes that correlations
between light- and b-jets are induced by marginalizing over the labels, the
model acquires classification power for the underlying processes (that we can
match to 4-top and $t\overline{t}W^{\pm}$) by learning the induced
correlations to achieve explanatory power over the full data distribution.

Figure 8: Same as Fig. 7 but for projected High-Luminosity expected event
rates.

The different hypotheses become better distinguishable at larger luminosities.
This is seen in Fig. 8 where we show the results for High-Luminosity LHC
projections with $\mathcal{L}=3000\text{ fb}^{-1}$. We observe that 4-top
exhibits larger deviations from independence than $t\overline{t}W^{\pm}$. In
particular $N_{j}$ and $N_{b}$ independence can be ruled out for the 4-top
distribution with realistic b-tagging. This in turn causes the full data
distribution to be tilted towards lower p-values for the conditionally
independent null hypothesis. The $t\overline{t}W^{\pm}$ does not exhibit the
same behavior. We verify that the MI of both processes decreases considerably
in the case of (near) perfect b-tagging. In particular, joint (black) 4-top
distribution is much closer to its marginalized (green) counterpart which is
also reflected in the full data conditional MI distribution This implies that
imperfect b-tagging is indeed an important factor behind observed deviations
from the conditional independence hypothesis although it is not the only one.
Because we are considering a probabilistic model for the data, a feasible
sophistication of this model that includes b-tagging efficiencies as a random
variable could restore conditional independence while keeping the number of
parameters under control. Such incorporation of the b-tagging efficiency would
be the analogue to the introduction of an associated nuisance parameter in
traditional statistical analyses. Another key feature at HL-LHC luminosity is
that for all four pseudo-datasets full independence between $N_{j}$ and
$N_{b}$ can be ruled out, as evidenced by the third column. For perfect
b-tagging, we can conclude that conditional independence is a valid
approximation which yields learnable distributions with discriminatory power
between processes. If imperfect b-tagging is taken into account in the
generative model then conditional independence remains a valid modelling
hypothesis with explanatory power for the full range of luminosities expected
at the LHC.

##  V Conclusions

In summary, we have proposed a new technique to extract signal and background
features and fractions relevant for measurements of four-top production at the
LHC using Bayesian Inference on the $N_{j}$ and $N_{b}$ jet multiplicity
distributions. It relies on the assumption of conditional (upon signal and
background class) independence of the inferred distributions and harnesses the
resulting correlations between $N_{j}$ and $N_{b}$ within each class. The
algorithm is weakly-supervised since, in addition to data (in the signal
region), it only relies on imperfect a priori knowledge how the signal and
background differ in their $N_{j}$ and $N_{b}$ distributions. Using these
results we have proposed a novel approach to test or tune MC predictions in
the signal region. Alternatively, it could allow to measure four-top
production cross-section and/or test for NP effects in a novel way that
alleviates the dependence on MC simulations altogether, as also proposed in
Ref. [19]. One could for instance tune the MC in the signal region using the
class-$0$ (background) $N_{j}$ and $N_{b}$ distributions and then simulate the
signal using the tuned MC to check whether its predicted fraction in the
2LSS++ sample agrees with the predictions in $p(\pi|X)$. Moreover, one can
also check whether the MC signal $N_{j}$ and $N_{b}$ distributions match the
$p(\alpha_{1,i}|X)$ and $p(\beta_{1,i}|X)$ inferred by the algorithm. Using
these ideas one would be able effectively to compute acceptances with a MC
tuned in-situ in the signal region, while simultaneously measure the four-top
cross-section, or study potential NP contributions to the signal or the
backgrounds.

Certainly, our method as presented in Sec. II is general and applicable also
to other particle physics scenarios beside four-top production and potentially
opens new venues of searches for NP at colliders. Certainly however, much
further work is needed to implement these techniques into feasible
experimental analyses.

Acknowledgements JFK acknowledges the financial support from the Slovenian
Research Agency (grant No. J1-3013 and research core funding No. P1-0035). BD
acknowledges funding from BMBF. DAF has received funding from the European
Research Council (ERC) under the European Union’s Horizon 2020 research and
innovation program under grant agreement 833280 (FLAY), and by the Swiss
National Science Foundation (SNF) under contract 200021-175940. We thank the
Referee for his/her report, which has considerably improved the content of the
article.

## References

  * Kasieczka  _et al._ [2021] G. Kasieczka _et al._ , The LHC Olympics 2020: A Community Challenge for Anomaly Detection in High Energy Physics, (2021), arXiv:2101.08320 [hep-ph] .
  * Kasieczka  _et al._ [2020] G. Kasieczka, B. Nachman, M. D. Schwartz, and D. Shih, ABCDisCo: Automating the ABCD Method with Machine Learning, (2020), arXiv:2007.14400 [hep-ph] .
  * Ghosh  _et al._ [2021] A. Ghosh, B. Nachman, and D. Whiteson, Uncertainty Aware Learning for High Energy Physics, (2021), arXiv:2105.08742 [physics.data-an] .
  * Benkendorfer  _et al._ [2020] K. Benkendorfer, L. L. Pottier, and B. Nachman, Simulation-Assisted Decorrelation for Resonant Anomaly Detection, (2020), arXiv:2009.02205 [hep-ph] .
  * Choi  _et al._ [2020] S. Choi, J. Lim, and H. Oh, Data-driven Estimation of Background Distribution through Neural Autoregressive Flows, (2020), arXiv:2008.03636 [hep-ph] .
  * Flesher  _et al._ [2020] F. Flesher, K. Fraser, C. Hutchison, B. Ostdiek, and M. D. Schwartz, Parameter Inference from Event Ensembles and the Top-Quark Mass, (2020), arXiv:2011.04666 [hep-ph] .
  * Lillie  _et al._ [2008] B. Lillie, J. Shu, and T. M. P. Tait, Top Compositeness at the Tevatron and LHC, JHEP 04, 087, arXiv:0712.3057 [hep-ph] .
  * Kumar  _et al._ [2009] K. Kumar, T. M. P. Tait, and R. Vega-Morales, Manifestations of Top Compositeness at Colliders, JHEP 05, 022, arXiv:0901.3808 [hep-ph] .
  * Acharya  _et al._ [2009] B. S. Acharya, P. Grajek, G. L. Kane, E. Kuflik, K. Suruliz, and L.-T. Wang, Identifying Multi-Top Events from Gluino Decay at the LHC, (2009), arXiv:0901.3367 [hep-ph] .
  * Kim  _et al._ [2016] J. H. Kim, K. Kong, S. J. Lee, and G. Mohlabeng, Probing TeV scale Top-Philic Resonances with Boosted Top-Tagging at the High Luminosity LHC, Phys. Rev. D94, 035023 (2016), arXiv:1604.07421 [hep-ph] .
  * Liu and Mahbubani [2016] D. Liu and R. Mahbubani, Probing top-antitop resonances with $t\bar{t}$ scattering at LHC14, JHEP 04, 116, arXiv:1511.09452 [hep-ph] .
  * Aguilar-Saavedra and Santiago [2012] J. A. Aguilar-Saavedra and J. Santiago, Four tops and the $t\bar{t}$ forward-backward asymmetry, Phys. Rev. D85, 034021 (2012), arXiv:1112.3778 [hep-ph] .
  * Camargo-Molina  _et al._ [2018] J. E. Camargo-Molina, A. Celis, and D. A. Faroughy, Anomalies in Bottom from new physics in Top, Phys. Lett. B 784, 284 (2018), arXiv:1805.04917 [hep-ph] .
  * Alvarez  _et al._ [2019] E. Alvarez, A. Juste, and R. M. S. Seoane, Four-top as probe of light top-philic New Physics, JHEP 12, 080, arXiv:1910.09581 [hep-ph] .
  * Darmé  _et al._ [2021] L. Darmé, B. Fuks, and F. Maltoni, Top-philic heavy resonances in four-top final states and their EFT interpretation, (2021), arXiv:2104.09512 [hep-ph] .
  * Khatibi and Khanpour [2021] S. Khatibi and H. Khanpour, Probing four-fermion operators in the triple top production at future hadron colliders, Nucl. Phys. B 967, 115432 (2021), arXiv:2011.15060 [hep-ph] .
  * Banelli  _et al._ [2021] G. Banelli, E. Salvioni, J. Serra, T. Theil, and A. Weiler, The Present and Future of Four Top Operators, JHEP 02, 043, arXiv:2010.05915 [hep-ph] .
  * Cao  _et al._ [2021] Q.-H. Cao, J.-N. Fu, Y. Liu, X.-H. Wang, and R. Zhang, Probing Top-philic New Physics via Four-Top-Quark Production, (2021), arXiv:2105.03372 [hep-ph] .
  * Alvarez  _et al._ [2020] E. Alvarez, F. Lamagna, and M. Szewc, Topic Model for four-top at the LHC, JHEP 01, 049, [JHEP20,049(2020)], arXiv:1911.09699 [hep-ph] .
  * Alvarez  _et al._ [2017] E. Alvarez, D. A. Faroughy, J. F. Kamenik, R. Morales, and A. Szynkman, Four Tops for LHC, Nucl. Phys. B 915, 19 (2017), arXiv:1611.05032 [hep-ph] .
  * Aad  _et al._ [2021] G. Aad _et al._ (ATLAS), Measurement of the $t\bar{t}t\bar{t}$ production cross section in $pp$ collisions at $\sqrt{s}$=13 TeV with the ATLAS detector, (2021), arXiv:2106.11683 [hep-ex] .
  * Sirunyan  _et al._ [2019] A. M. Sirunyan _et al._ (CMS), Search for the production of four top quarks in the single-lepton and opposite-sign dilepton final states in proton-proton collisions at $\sqrt{s}$ = 13 TeV, JHEP 11, 082, arXiv:1906.02805 [hep-ex] .
  * Sirunyan  _et al._ [2020] A. M. Sirunyan _et al._ (CMS), Search for production of four top quarks in final states with same-sign or multiple leptons in proton-proton collisions at $\sqrt{s}=$ 13 TeV, Eur. Phys. J. C 80, 75 (2020), arXiv:1908.06463 [hep-ex] .
  * Aad  _et al._ [2020] G. Aad _et al._ (ATLAS), Evidence for $t\bar{t}t\bar{t}$ production in the multilepton final state in proton–proton collisions at $\sqrt{s}=13$ TeV with the ATLAS detector, Eur. Phys. J. C 80, 1085 (2020), arXiv:2007.14858 [hep-ex] .
  * Dillon  _et al._ [2019] B. M. Dillon, D. A. Faroughy, and J. F. Kamenik, Uncovering latent jet substructure, Phys. Rev. D100, 056002 (2019), arXiv:1904.04200 [hep-ph] .
  * Dillon  _et al._ [2020] B. M. Dillon, D. A. Faroughy, J. F. Kamenik, and M. Szewc, Learning the latent structure of collider events, JHEP 10, 206, arXiv:2005.12319 [hep-ph] .
  * Dillon  _et al._ [2021] B. M. Dillon, T. Plehn, C. Sauer, and P. Sorrenson, Better Latent Spaces for Better Autoencoders, (2021), arXiv:2104.08291 [hep-ph] .
  * Bishop [2006] C. M. Bishop,  _Pattern recognition and machine learning_, Information science and statistics (Springer, New York, NY, 2006) softcover published in 2016\.
  * cod [2021] Bayesian inference for four tops at the lhc, https://github.com/ManuelSzewc/bayes-4tops (2021).
  * Sokal [1996] A. Sokal, Monte carlo methods in statistical mechanics: Foundations and new algorithms note to the reader (1996).
  * Foreman-Mackey  _et al._ [2013] D. Foreman-Mackey, D. W. Hogg, D. Lang, and J. Goodman, emcee: The mcmc hammer, Publications of the Astronomical Society of the Pacific 125, 306–312 (2013).
  * Alwall  _et al._ [2014] J. Alwall, R. Frederix, S. Frixione, V. Hirschi, F. Maltoni, O. Mattelaer, H. S. Shao, T. Stelzer, P. Torrielli, and M. Zaro, The automated computation of tree-level and next-to-leading order differential cross sections, and their matching to parton shower simulations, JHEP 07, 079, arXiv:1405.0301 [hep-ph] .
  * Sjöstrand  _et al._ [2015] T. Sjöstrand, S. Ask, J. R. Christiansen, R. Corke, N. Desai, P. Ilten, S. Mrenna, S. Prestel, C. O. Rasmussen, and P. Z. Skands, An Introduction to PYTHIA 8.2, Comput. Phys. Commun. 191, 159 (2015), arXiv:arXiv:1410.3012 [hep-ph] .
  * de Favereau  _et al._ [2014] J. de Favereau, C. Delaere, P. Demin, A. Giammanco, V. Lemaître, A. Mertens, and M. Selvaggi (DELPHES 3), DELPHES 3, A modular framework for fast simulation of a generic collider experiment, JHEP 02, 057, arXiv:1307.6346 [hep-ex] .
  * Goebel  _et al._ [2005] B. Goebel, Z. Dawy, J. Hagenauer, and J. Mueller, An approximation to the distribution of finite sample size mutual information estimates (2005) pp. 1102 – 1106 Vol. 2.
  * Baker and Cousins [1984] S. Baker and R. D. Cousins, Clarification of the use of chi-square and likelihood functions in fits to histograms, Nuclear Instruments and Methods in Physics Research 221, 437 (1984).


# Spectrum of Large N Glueballs: Holography vs Lattice

Anatoly Dymarskya and Dmitry Melnikovb

###### Abstract

Recently there has been a notable progress in the study of glueball states in
lattice gauge theories, in particular extrapolating their spectrum to the
limit of large number of colors $N$. In this note we compare the large $N$
lattice results with the holographic predictions, focusing on the Klebanov-
Strassler model. We note that glueball spectrum demonstrates approximate
universality across a range of gauge theory models. Because of this
universality the holographic models can give reliable predictions for the
spectrum of pure $SU(N)$ Yang-Mills theories with and without supersymmetry.
This is especially important for the supersymmetric theories, for which no
firm lattice predictions exist yet, and the holographic models remain the most
tractable approach. For non-supersymmetric pure $SU(N)$ theories with large
$N$ we find an agreement within 5-8% between the lattice and holographic
predictions for the mass ratios of the lightest states in various sectors. In
particular both lattice and holography give predictions for the $2^{++}$ and
$1^{--}$ mass ratio, consistent with the known constraints on the pomeron and
odderon Regge trajectories.

aDepartment of Physics and Astronomy,
University of Kentucky, Lexington, KY, USA 40506
bInternational Institute of Physics, Federal University of Rio Grande do
Norte,
Campus Universitário, Lagoa Nova, Natal-RN 59078-970, Brazil

In [1] the D0 and TOTEM collaborations announced a $3.4\sigma$ divergence of
the $pp$ and $p\bar{p}$ cross sections, compatible with a $t$-channel exchange
of a colorless $C$-odd particle, an “odderon” [2]. Further work is under way
to improve the statistics of this result. The principle candidate for the
odderon is a $C$-odd three-gluon bound state classified as $1^{--}$ in the
$J^{PC}$ notations. The experimental evidence for odderon highlights
importance of pure glue states, the glueballs. Such states are important for
understanding of the detailed structure of strong interactions [3] as well as
possible extensions of physics beyond the Standard Model, see e.g. [4, 5, 6].

Existence of pure glue bound states was conjectured long ago, but their
experimental detection, and elucidating their physical properties, turned out
to be a complicated task. The main theoretical challenge is the heavy mixing
of such nonquarkyonic short-lived states with heavy excited states of ordinary
mesons. As a result most of the presently available data concerning glueballs
come from the lattice simulations. Mixing with quarks is a challenge for the
lattice as well, so most of the time the studies assume the quenched
approximation, effectively studying the pure glue (or pure Yang-Mills) sector
of QCD.

The standard reference for the spectrum of glueballs is the work of
Morningstar and Peardon [7], which found 12 lightest states with spins varying
from $J=0$ to $J=3$ in the pure glue $SU(3)$ Yang-Mills theory.111See [8, 9,
10] for yet earlier work on the lattice spectrum. The ensuing lattice studies
improved the original predictions of the glueball masses [11, 12, 13],
observing additional excited states, studied the effects of introducing quarks
[14, 15, 16, 17, 18, 19, 20, 21], investigated dependence on the gauge group
and its rank [22, 23, 24, 25, 26, 27, 28] and attempted to estimate the decay
constants of a few lightest states [11], cf. [29, 30, 31, 32].

The glueball spectrum exhibits relatively mild dependence on the number of
colors, in the ’t Hooft limit222In the $Sp(N)$ case the leading correction to
the $N\to\infty$ value is $O(1/N)$.

$\begin{array}[]{c}\displaystyle m\ \simeq\
m_{\infty}+\frac{c}{N^{2}}+\dots\end{array}$ (1)

expected from a large $N$ expansion and corroborated by the lattice analysis
[33, 34, 35]. In table 1 we show lattice results for $SU(N)$ [23] and $Sp(N)$
[27] for the mass of the lightest $0^{++}$ glueball and its ratio with the
mass of the lightest excitation, the $2^{++}$ state. For $SU(N)$ with $N$
ranging from $N=2$ to $N=\infty$ the variation of masses is within 14-16%,
while the variation of the ratios is even smaller, about 10%. In the $Sp(N)$
case the variation of the ratios is even smaller, about 5%.

$G$ | $SU(2)$ | $SU(3)$ | $SU(4)$ | $SU(6)$ | $SU(8)$ | $SU(\infty)$
---|---|---|---|---|---|---
$m_{0^{++}}$ | 3.78 | 3.55 | 3.36 | 3.25 | 3.55 | 3.31
$m_{2^{++}}/m_{0^{++}}$ | 1.44 | 1.35 | 1.45 | 1.46 | 1.32 | 1.46
$G$ | $Sp(1)$ | $Sp(2)$ | $Sp(3)$ | $Sp(4)$ | – | $Sp(\infty)$
$m_{2^{++}}/m_{0^{++}}$ | 1.41 | 1.41 | 1.48 | 1.41 | – | 1.47
Table 1: Lattice masses (in QCD string tension units [27]) and mass ratios of
the $0^{++}$ and $2^{++}$ glueballs in the pure glue $SU(N)$ [23] and $Sp(N)$
theories [27].

Another important observation is a reasonably small effect of quark mixing in
the unquenched version of the Yang-Mills theory. The phenomenological OZI rule
[36, 37, 38] does not favor quark mixing with purely gluonic states. Lattice
simulations confirm this effect to certain extent, as can be seen from the
data in table 2. The match between quenched and unquenched cases is better for
the lighter states (5% for $0^{++}$), and the ratio $m_{2^{++}}/m_{0^{++}}\sim
1.46$ obtained for QCD with three flavors is compatible with the pure glue
$SU(N)$ and $Sp(N)$ results (approximate universality of this ratio was also
acknowledged in [26]). Mass ratios appears to be even more stable with respect
to model variation, which might indicate that not only $m_{\infty}$, but also
$c$ in equation (1) are universal.

$J^{PC}$ | $0^{++}$ | $2^{++}$ | $1^{+-}$ | $0^{-+}$ | $m_{2^{++}}/m_{0^{++}}$ | $m_{1^{+-}}/m_{2^{++}}$ | $m_{0^{-+}}/m_{2^{++}}$
---|---|---|---|---|---|---|---
YM | 1710 | 2390 | 2980 | 3640 | 1.40 | 1.25 | 1.52
QCD3 | 1795 | 2620 | 3270 | 4490 | 1.46 | 1.25 | 1.71
Table 2: Comparison of the lattice glueball masses (measured in MeV) in the
quenched Yang-Mills and QCD with three flavors [18]. The ratio of the masses
of $m_{2^{++}}$ and $m_{0^{++}}$ for the two cases are $1.40$ and $1.46$
respectively.

Approximate universality of the glueball spectrum, in the sense of weak
$N$-dependence, makes it natural to study it using holographic models.
Classical gravity limit of holography [39, 40, 41], the so-called
gauge/gravity duality, can only capture the limit of large number of colors
$N\to\infty$ and large ’t Hooft coupling constant
$\lambda=g_{YM}^{2}N\to\infty$ of Yang-Mills theory. Without taking into
account quantum corrections to gravity one can only hope to access some
universal large $N$ properties of hadrons. Lattice simulations indicate that
the glueball spectrum, which shows approximate $N$ independence, is one of
such quantities which can be accessed holographically.

Sparing the details, the holographic calculation in the classical gravity
limit boils to choosing a background – a solution of the appropriate
(super)gravity equations – that is dual to the gauge theory of interest, and
calculating the spectrum of linearized excitations by expanding equations of
motion to linear order above the background solution. The spectrum of gravity
excitations is then equivalent to the spectrum of gauge theory. Higher order
correlation functions are also available though require calculations beyond
the linear order.

Early take on glueballs in holography [42, 43] was in the case of the so-
called Witten’s model [44]. The background of this model is the near-horizon
limit of the world-volume theory of D4 branes compactified on a circle. This
background preserves no supersymmetry and can be used to qualitatively
describe gauge theory spectrum. As a drawback, it exhibits features
incompatible with the Yang-Mills theory, namely mass degeneracy of different
glueball families, which points to additional symmetries beyond those of the
Yang-Mills theory. It also exhibits larger spacing of masses, as demonstrated
in table 3. In Witten’s model the ratio $m_{2^{++}}/m_{0^{++}}\sim 1.74$,
significantly larger than in the pure $SU(N)$ case. The discrepancy is even
bigger for heavier states.

$J^{\,PC}$ | $SU(\infty)$ | $Sp(\infty)$ | BMT | BBCD | BBCN
---|---|---|---|---|---
${0^{++}}$ | 1 | 1 | 1 | 1 | 1
${2^{++}}$ | 1.49 | 1.47 | 1.74 | 1.48 | 1.56
${0^{-+}}$ | 1.53 | 1.54 | 2.09 | – | –
$1^{+-}$ | 1.88 | – | 2.70 | – | –
$1^{--}$ | 2.32 | – | 3.37 | – | –
$0^{+-}$ | 3.01 | – | – | – | –
$0^{++\ast}$ | 1.89 | 1.94 | 2.53 | 1.63 | 1.83
$2^{++\ast}$ | 2.11 | – | 2.76 | 2.15 | 2.49
Table 3: Spectra of the lightest glueballs in Witten’s (BMT [43]) and hard
wall models (BBC [45]) compared to the state of the art lattice extrapolations
of $m_{\infty}$ for $SU(N)$ [28] and $Sp(N)$ [27]. The masses are normalized
to the mass of the $0^{++}$ state. Models BBCN and BBCD used different
boundary conditions (Neumann or Dirichlet) in the calculation of the spectrum.

A straightforward approach to model hadron physics holographically can be
taken via bottom-up models, the so-called hard wall model would be the
simplest case [46, 47, 48]. For the purpose of this note, the Light-Front
Holography models [49] can also be attributed to this class. In the hard wall
model the background is five-dimensional anti de Sitter space ($AdS_{5}$),
whose group of symmetries is precisely the conformal group in 3+1 dimensions,
and a cutoff is introduced to break this symmetry explicitly. The dual theory
to such a background at high energies would qualitatively be a
(3+1)-dimensional conformal gauge theory, with an IR scale defining the
masses, as in pure Yang-Mills or QCD. One then considers wave equations for
matter fields of various spin in the $AdS_{5}$ background with the boundary
conditions at the cut off radius, i.e. at the hard wall. This mathematical
problem leads to the spectrum of light states, as calculated by [45] and shown
in table 3.

Remarkably, the plain hard wall model predicts the value
$m_{2^{++}}/m_{0^{++}}\sim 1.48$, that is very close to the extrapolated ratio
of $m_{\infty}$ calculated on the lattice. This “magic number” can be
expressed as a ratio of the first non-trivial zeroes $x_{2,1}$ and $x_{4,1}$
of Bessel functions $J_{2}(x)$ and $J_{4}(x)$,

$\begin{array}[]{c}\displaystyle\text{hard
wall:}\quad\frac{m_{2^{++}}}{m_{0^{++}}}\ =\ \frac{x_{4,1}}{x_{2,1}}\ \simeq\
1.47759\,.\end{array}$ (2)

The hard-wall model also has its drawbacks. In particular, it does not include
states with non-trivial parity and charge conjugation, since the model does
not have a natural way to implement corresponding symmetries. In contrast,
top-down holographic constructions inherit these symmetries from string
theory, as in the case of Witten’s model.

The discussion above, which loosely follows historic developments, suggests
that a more nuanced top-down model would be necessary to accurately describe
large $N$ physics in pure Yang-Mills and QCD-like theories. The Klebanov-
Strassler (KS) model [50] was proposed with this goal in mind. It is a top-
down holographic model, which possesses rich IR physics. It is a low-energy
limit of D3 and D5 branes in type IIB string theory compactified on a six-
dimensional cone (conifold). The compactification gives rise to a type IIB
supergravity background on a space which is a warped product of $AdS_{5}$ and
a pair of spheres $S^{3}\times S^{2}$ [51, 52, 53]. Smoothing out the
singularity of the cone (deformation of the conifold) introduces an IR scale
that breaks (3+1)-dimensional conformal invariance and provides an interesting
example of a holographic RG flow with a logarithmic running of the coupling
constants. The background has $SU(2)\times SU(2)\times U(1)_{\rm B}$ global
symmetries, so its spectrum is organized in the irreducible representations of
this group.

The Klebanov-Strassler background preserves ${\cal N}=1$ supersymmetry, so the
dual gauge theory is a non-conformal ${\cal N}=1$ supersymmetric Yang-Mills
theory with additional matter fields transforming under the global symmetry
group, in the presence of a particular superpotential. In the IR the theory
flows to a strongly coupled fixed point (the corresponding RG flow is known as
a “cascade” of Seiberg dualities [50, 54]) which was initially thought to be
in the universality class of the pure supersymmetric $SU(N)$ Yang-Mills
theory. It was then understood that the IR limit of the Klebanov-Strassler
theory has additional massless states due to spontaneous breaking of the
baryon $U(1)_{\rm B}$ symmetry by the baryonic operators [55, 56, 57].

As a close relative of the pure ${\cal N}=1$ Yang-Mills, the Klebanov-
Strassler gauge theory shares with it many core properties. Apart from the
mentioned logarithmic running of the coupling constants, the holographic
theory exhibits the same mechanism of breaking of the $U(1)$ R-symmetry [58].
The singlet sector mostly contains the super Yang-Mills operators [59, 60]. In
view of the approximate universality of the spectrum, one can hope that
$SU(2)\times SU(2)$-singlet states reproduce the appropriate limit of the low-
energy spectrum of the supersymmetric, or even pure bosonic $SU(N)$ Yang-
Mills. With this logic in mind, we summarize below the results for the singlet
states masses in the Klebanov-Strassler theory.

The spectrum of singlet fluctuations of the Klebanov-Strassler background was
studied in a series of papers including [61, 62, 63, 64, 65, 66, 67, 68, 69,
70]. (See also [71, 72, 73, 74, 75, 76] for the studies of non-singlet states
or deformations of the Klebanov-Strassler theory.) The structure of the
spectrum in the supergravity limit was explained in [68, 69, 70]: it contains
two massless scalar supermultiplets and thirteen massive supermultiples
including one graviton, two gravitino, four vector and six scalar
supermultiplets. The latter are massive ${\cal N}=1$ supermultiplets labeled
by the highest spin component, except for the scalar multiplet, which is
conventionally labeled by its lowest spin component.

Our goal is to compare this spectrum of Klebanov-Strassler theory with the
pure gauge ${\cal N}=1$ Yang-Mills. Since the KS theory exhibits additional
global symmetries only certain states (multiplets) are the counterparts of the
pure YM. We elaborate this point below. The massless states of the KS theory
are not part of the Yang-Mills spectrum. One should also exclude multiplet
coming from the $U(1)_{\rm B}$ sector: a scalar multiplet as well as a vector
multiplet containing $0^{+-}$ and $1^{+-}$ [66, 67]. The remaining eleven
massive multiplets exist in the spectrum of the pure supersymmetric Yang-
Mills. They are listed in table 4 and their spectrum is illustrated by figure
1.

| $J^{\,PC}$ | Multiplet | $m^{\rm GS}$ | $m^{\ast}$ | Ref.
---|---|---|---|---|---
1 | $1^{++},2^{++}$ | graviton | 1 | 1.51 | [63]
2 | $1^{+-},1^{--}$ | gravitino | 1.30 | 1.85 | [67]
3 | $1^{+-},1^{--}$ | gravitino | 1.64 | 2.15 | [67]
4 | $0^{--},1^{--}$ | vector | 1.47 | 2.00 | [66]
5 | $0^{+-},1^{+-}$ | vector | 2.01 | 2.55 | [66]
6 | $0^{++},1^{++}$ | vector | 1.99 | 2.53 | [68]
7 | $0^{++},0^{-+}$ | scalar | 0.421 | 0.894 | [64]
8 | $0^{++},0^{-+}$ | scalar | 0.640 | 1.25 | [64]
9 | $0^{++},0^{-+}$ | scalar | 1.11 | 1.58 | [64]
10 | $0^{++},0^{-+}$ | scalar | 1.36 | 1.84 | [64]
11 | $0^{++},0^{-+}$ | scalar | 2.32 | 2.87 | [64]

$J^{\,PC}$ | $m^{\rm GS}$ | $m^{\ast}$
---|---|---
$2^{++}$ | 1 | 1.43
$1^{+-}$ | 1.25 | 1.62
$1^{--}$ | 1.58 | $\geq 1.88$
|  |
$0^{+-}$ | $\geq 2.01$ | –
|  |
|  |
$0^{++}$ | 0.668 | 1.27
$0^{-+}$ | 1.02 | 1.53
|  |
$0^{++}$ | – | –

Table 4: Spectrum of the lightest supermultiplets of the supersymmetric Yang-
Mills sector of the Klebanov-Strassler theory (left table). The bosonic
members of the multiplets are indicated. The masses of the ground and first
excited states, in units of the mass of the ground state of $2^{++}$, are
shown, as well as the references to the works, from which the values are
extracted. For the scalar sector, the assignment was made based on the
analysis of [70, 77]. These masses are compared with the masses recently
obtained for the $SU(\infty)$ bosonic Yang-Mills on the lattice [28] (right
table). Lower bounds indicate states with known masses, but with a problem of
confirming the $J^{PC}$ numbers in the continuum limit.

Figure 1: Full light $SU(2)\times SU(2)$-singlet spectrum of the bosonic
states in the Klebanov-Strassler theory (left). The ground states of each
supermultiplet are shown in full color, while the excited states are faded.
All the mass values are computed in the units of the ground state of $2^{++}$.
The discrepancy of the $0^{++}$ [64] and $0^{-+}$ [70] sectors is visible. The
right plot compares the masses of the ground states in the supersymmetric
Yang-Mills sector (filled rectangles) and the ground states in $SU(N)$ bosonic
Yang-Mills, extrapolated to $N=\infty$ (empty rectangles). In the
$0^{++}/0^{-+}$ sector the results of [64] are shown in full color, while that
of [70] are faded. States that do not have pairs correspond to hybrid
glueballs with no analog in the bosonic theory, except for the heaviest
$0^{++}$ that has not been approached on the lattice.

In table 4 we measure all masses in the units of the mass of the ground state
of the graviton multiplet, the $2^{++}$ glueball. This state is particularly
easy to analyze on the supergravity side, since the $2^{++}$ fluctuation above
the background decouples from all the other modes [63, 78]. It thus serves a
natural reference point for the remaining states. The charge conjugation-odd
($C$-odd) sector is also relatively simple to analyze holographically. The
masses of the $1^{+-}$, $1^{--}$ and $0^{+-}$ glueballs can be compared with
the $SU(N)$ bosonic Yang-Mills extrapolated to $N=\infty$ [28]. We note the
holographic masses have the same hierarchy as on the lattice. Even more, they
match the $SU(\infty)$ ground state estimates within 5%. For the excited
states the correspondence is not so good: the holographic supersymmetric
masses of excited states grow more rapidly with the excitation number than
their lattice bosonic counterparts.

Masses of some of the states calculated on the lattice in [28] and shown in
table 4, appear as lower bounds, as for the $1^{--\ast}$ and $0^{+-}$ states.
This is because the unambiguous identification of the $J^{PC}$ quantum numbers
of these states is not yet possible on the lattice, and one can only classify
them tentatively as $1^{--\ast}$ and $0^{+-}$, based on general mass hierarchy
arguments. We note that in this case the holographic results support this
identification.

Not all the multiplets of the supersymmetric Yang-Mills theory survive in the
quenched approximation. Some of the multiplets only contain composite states
of fermion operators. The two examples are vector multiplets 4 and 6 in table
4. For such states there are no lattice results yet to compare.

The most complex is the scalar sector, originally studied in [62, 64]. There
are six pairs of scalar (and also six pseudoscalar) multiplets which mix
together. The scalars contain the purely bosonic $0^{++}$ and $0^{-+}$ modes
and scalars associated with the composite gluino operators. The scalar
multiplet containing the gluino bilinear $\lambda\lambda$ state is expected to
be the lightest in the spectrum of the supersymmetric Yang-Mills. Apart from
${\rm tr}\,\lambda\lambda$ and the scalars due to ${\rm
tr}\,F_{\mu\mu}F^{\mu\nu}$ and ${\rm tr}\,F_{\mu\nu}\tilde{F}^{\mu\nu}$
operators, the holographic model includes other $0^{++}$ and $0^{-+}$ gluon-
gluino states associated with the composite operators of dimension
$\Delta=3,4,5,6,7,8$ (see [59, 60, 69, 70] for the classification).

The holographic spectrum of $0^{++}$ states was originally calculated in [64].
The result of that calculation is a list of lowest masses, which on physical
grounds are expected to belong to independent towers associated with seven
individual particles and their excitations. In the conformal case these would
be related to seven composite operators discussed in [59, 60]. It is not clear
how to split excitations into towers from the first principles. That is except
for one tower associated with a particle from the vector supermultiplet. It
can be separated by comparison with its superpartner vector mode, which
decouples from all other modes on gravity side, see [68]. For all other
states, below we use that masses approximately follow “Regge trajectories”
associated with the radial quantum number $n$, and identify ground states for
each tower.

In [70] the spectrum of $0^{-+}$, the superpartners of $0^{++}$, was
calculated as a consistency check. It was proposed that the following ordering
of states is consistent with [64, 70] and with the lattice results [77],

$\begin{array}[]{c}\displaystyle
m_{\lambda\lambda}<m_{0^{++}}<m_{\lambda\lambda}^{*}<m_{0^{-+}}<m_{0^{++}}^{*}<\cdots\,,\end{array}$
(3)

where $m_{\lambda\lambda}$ is the mass of the ground state of the gluino
bilinear, $m_{\lambda\lambda}^{\ast}$ is its first excited state, $m_{0^{++}}$
and $m_{0^{-+}}$ are masses of the states in the spectrum of the two gluon
operators ${\rm tr}\,F_{\mu\nu}F^{\mu\nu}$ and ${\rm
tr}\,F_{\mu\nu}\tilde{F}^{\mu\nu}$. Unlike gluino operators, which pose a
challenge for lattice methods, the latter operators are commonly studied on
the lattice.

In the present note we extend (3) and complete the hierarchy of masses by
adding two remaining multiplets. One of these correspond to four-gluon
operator ${\rm tr}\,(F_{\mu\nu}F^{\mu\nu})^{2}$, which is in principle
accessible on the lattice. This is done by applying quadratic fit to the
values of masses squared. This approach worked very well in the past for other
sectors of the KS theory [63, 65, 66, 67, 69]. The results of the fitting of
the pseudoscalar sector [70] are shown in figure 2. The following fits work
very well for the heavy part of the spectrum, while for the light states we
end up with a few noticeable deviations,

$\displaystyle m_{\lambda\lambda}^{2}\ =\ 0.223+0.499n+0.257n^{2}\,,$
$\displaystyle m^{\rm GS}\ =\ 0.472\,,$ (4) $\displaystyle m_{0^{++}}^{2}\ =\
0.480+0.854n+0.260n^{2}\,,$ $\displaystyle m^{\rm GS}\ =\ 0.693\,,$ (5)
$\displaystyle m_{0{-+}}^{2}\ =\ 1.22+1.26n+0.257n^{2}\,,$ $\displaystyle
m^{\rm GS}\ =\ 1.11\,,$ (6) $\displaystyle m_{10}^{2}\ =\
1.97+1.15n+0.267n^{2}\,,$ $\displaystyle m^{\rm GS}\ =\ 1.40\,,$ (7)
$\displaystyle m_{AB}^{2}\ =\ 4.24+2.14n+0.261n^{2}\,,$ $\displaystyle m^{\rm
GS}\ =\ 2.06\,,$ (8) $\displaystyle m_{11}^{2}\ =\ 5.92+2.49n+0.260n^{2}\,,$
$\displaystyle m^{\rm GS}\ =\ 2.43\,.$ (9)

Here all the values are given in units of the mass of the ground state of
$2^{++}$ and $m^{\rm GS}$ stands for $n=0$, i.e. extrapolated ground state
value of the corresponding tower. Here $m_{10}$ and $m_{11}$ are pseudoscalar
predictions for the masses of the ground states in the entries 10 and 11 of
table 4. Masses $m_{AB}$ cannot be compared with the spectrum of the pure
supersymmetric Yang-Mills, as the corresponding operators come from the
$U(1)_{\rm B}$ sector and would be absent in pure glue theory. We note that
using $m^{\rm GS}=m(n=0)$ rather than actual lowest mass slightly improves the
convergence to the $SU(\infty)$ lattice values.

Figure 2: Quadratic fits of the $0^{-+}$ spectrum of $m^{2}$ found in [70].
The black dots (as well as the grids) indicate the masses and the curves are
fits (4)-(9). $n$ is the “radial” quantum number labeling excited states. The
left panel shows the behavior of the fits for $m^{2}\leq 10$. The intercepts
correspond to the ground states. The right panel illustrates the quality of
the fits for $m^{2}\geq 10$. All mass values are shown in the units of the
ground state mass of $2^{++}$.

It is necessary to mention that the spectrum of pseudoscalars computed in
[70], although similar with the spectrum of $0^{++}$ from [64], and features
the same hierarchy and large $n$ asymptotic, exhibits a consistent divergence
of the masses of two out of six towers. This divergence is particularly
pronounced for low masses. For example, in [64] the masses of the three
lightest scalar ground states are $m_{\lambda\lambda}\simeq
0.511\,m_{2^{++}}$, $m_{0^{++}}\simeq 0.701\,m_{2^{++}}$ and $m_{0^{-+}}\simeq
1.15\,m_{2^{++}}$, cf. (4)-(9), to be compared with the entries 7, 8 and 9 in
table 4. See also figure 1. The spectrum of [64] was independently confirmed
in [74], so the status of the pseudoscalar calculation remains unclear.

We conclude by discussing relevance of holographic results in the context of
lattice simulations. Recent progress on the lattice largely confirms the
expectation of universality – approximate $N$ independence in the sense of (1)
– of the glueball masses. This allows to extrapolate the spectrum to the limit
of infinite number of colors, making possible comparison of the lattice
results for $SU(\infty)$ with the Yang-Mills subsector of the Klebanov-
Strassler theory. Before proceeding, let us stress a subtlety important for
any comparison of this kind.

The models conventionally studied on the lattice are asymptotically free and
have no tunable parameters except for dimensionful $\Lambda_{QCD}$, or
equivalently the confining string tension $\sigma$. Hence, the spectrum on the
lattice can be unambiguously expressed in units of $\sqrt{\sigma}$. Meanwhile,
in the KS theory, as well as other models of gauge/(super)gravity duality,
there is a free parameter, an exactly marginal coupling associated with the
value of the dilaton. On the dual field theory side this is the ’t Hooft
coupling $\lambda$.333In the KS theory, in addition to an exactly marginal
$\lambda$ there is also a running coupling. Accordingly in holographic models
the ratio $m_{J^{PC}}/\sqrt{\sigma}$ is $\lambda$-dependent, making direct
comparison with the lattice problematic. Furthermore in the strong coupling
limit $\lambda\gg 1$, when gravity approximation is valid,
$m_{J^{PC}}/\sqrt{\sigma}$ is suppressed by a positive power of $\lambda$, see
[50, 79]. Strictly speaking, the mass-tension ratios vanish in the
supergravity limit. However, for glueballs of small spin the leading scaling
of $m_{J^{PC}}/\sqrt{\sigma}$ with $\lambda$ is universal, rendering glueball
mass ratios $\lambda$-independent. They can be compared directly with the
lattice results, as discussed below.

For the higher spin states with $J\geq 2$, except for the $2^{++}$ which
behaves as $J<2$ states, in the strong coupling limit, the
$\lambda$-dependence is different. In fact these masses diverge in
supergravity limit, as they belong to massive string sector and their masses
are of order $1/\sqrt{\alpha^{\prime}}$. In other words these states are
beyond (super)gravity approximation, making comparison with the lattice
presently impossible.

With this said we proceed by comparing the holographic mass ratios with those
on the lattice. For the five lightest states with $J<2$ and for the $2^{++}$
state of the bosonic Yang-Mills theory we observed the same hierarchy of the
spectrum in both the Klebanov-Strassler theory and on the lattice.
Furthermore, for the five mass ratios of these six lightest states we observe
a better than ten percent consistency of the numerical values, as summarized
in table 5 for both $SU(\infty)$ and $SU(3)$ theories. These results provide
an additional evidence for the universality of the spectrum, and an
independent check for the consistency of the lattice approach.

$J^{\,PC}$ | $2^{++}$ | $0^{++}$ | $0^{-+}$ | $1^{+-}$ | $1^{--}$ | $0^{+-}$
---|---|---|---|---|---|---
Holography/$SU(\infty)$ | 1 | 0.959 | 1.081 | 1.037 | 1.038 | 0.999
Holography/$SU(3)$ | 1 | 0.920 | 1.027 | 1.048 | 0.967 | 1.057
Table 5: Ratios of glueball masses measured in the units of the $2^{++}$ mass,
$({m_{\rm hol}^{J^{PC}}/m_{\rm hol}^{2^{++}})/(m_{\rm latt}^{J^{PC}}}/m_{\rm
latt}^{2^{++}})$, for holographic (Klebanov-Strassler) and lattice predictions
for the lightest glueball ground states in the $SU(\infty)$ and $SU(3)$ Yang-
Mills theories [28].

Provided the observed universality is not coincidental, the holographic
glueball spectrum provides a rare opportunity to test non supersymmetry-
protected sectors of holographic correspondence. It shows that despite the
peculiarities of the holographic limit, first principle lattice calculations
can be used to access physical observables in this regime.

Note that the Klebanov-Strassler theory is not expected to give precisely the
same results as the $SU(\infty)$ bosonic Yang-Mills. The results summarized in
figure 1 is rather the spectrum of the supersymmetric $SU(\infty)$ Yang-Mills,
deformed by the presence of additional matter fields. The numerical
consistency of the results in the KS theory and in the bosonic Yang-Mills
suggests that the KS theory would give a similar, or even better match with
the spectrum of the supersymmetric pure Yang-Mills. The supersymmetric case
would give us a larger base for comparison with the multitude of additional
states, but unfortunately it is a much more challenging case to consider on
the lattice. We hope that such a comparison will be possible in the future.

As a separate question we would like to discuss implications of holographic
results for the odderon physics. First subtlety is that glueball Regge
trajectories are expected to have flatter slopes as compared to the universal
slope of other hadrons. This is attributed to the structure of glueballs,
which unlike mesons or baryons, are comprised of particles connected by a flux
tube in the adjoint representation. Alternatively glueballs can be seen as
closed fundamental flux tubes. Two different pictures indicate there is no
universal glueball Regge trajectory slope, but in both cases the value is
smaller than for other hadrons. These theoretical expectations are supported
by the lattice studies of pure glue theories and more broadly by experimental
studies of the $pp$ and $p\bar{p}$ cross sections. Lattice predicts the slope
of the leading pomeron trajectory, of the $2^{++}$ glueball, to be
$\alpha^{\prime}_{\rm P}\sim 0.28\alpha^{\prime}$ [80], consistent with the
experimentally observed value [81], where $\alpha^{\prime}$ is the slope of
the meson trajectory. The subleading trajectory of the lightest $0^{++}$
glueball state on the lattice exhibits a larger slope than
$\alpha^{\prime}_{\rm P}$ but still smaller than $\alpha^{\prime}$.

Now, if one assumes the leading odderon trajectory to have the same slope with
the leading pomeron trajectory, than the intercepts of the two are related by

$\begin{array}[]{c}\displaystyle j^{-}\ =\
1-(2-j^{+})\left(\frac{m_{1^{--}}^{2}}{m_{2^{++}}^{2}}\right).\end{array}$
(10)

The intercept of the leading pomeron trajectory is known experimentally to be
slightly above unity, $j^{+}\sim 1.08$ [82]. One could expect that the
robustness of the mass ratios, through (10), can then provide a good estimate
for $j^{-}$. However, this relation predicts a rather low value for the
intercept of the $1^{--}$ trajectory $j^{-}$, using lattice result for
${m_{1^{--}}^{2}}/{m_{2^{++}}^{2}}$ we find $j^{-}\sim-1.3$ (for holographic
models we find $j^{-}\sim-1.5$, which is consistent with the value obtained in
early lattice estimates [83, 84]). Above we used the mass ratios,

$\begin{array}[]{c}\displaystyle\frac{m_{1^{--}}}{m_{2^{++}}}\ \simeq\
1.58~{}\text{(lattice~{}\cite[cite]{[\@@bibref{}{Athenodorou:2021qvs}{}{}]})}\qquad\text{or}\qquad\frac{m_{1^{--}}}{m_{2^{++}}}\
\simeq\
1.64~{}\text{(holography~{}\cite[cite]{[\@@bibref{}{Dymarsky:2008wd}{}{}]})}\,.\end{array}$
(11)

Similar picture was outlined in [80], which used lattice to study Regge
trajectories. If one assumes the $1^{--}$ state to be on the same trajectory
with $3^{--}$, then this trajectory’s slope is numerically close to
$\alpha^{\prime}_{\rm P}$. At the same time, a perturbative analysis (hard
odderon) predicts the intercept $j^{-}\sim 1$ and a slope smaller than
$\alpha^{\prime}_{P}$ [85]. Interpolation between strong and weak coupling was
studied in [86] in ${\cal N}=4$ supersymmetric Yang-Mills theory. It shows
$j^{-}$ and $\alpha^{\prime}_{\rm P}$ can vary substantially between weak and
strong coupling, thus providing consistency to the picture above.

A common experimental view, however, is that the intercept of the odderon
trajectory is slightly below unity. Should $1^{--}$ belong to this trajectory,
the slope of the latter would be very small compared to that of the pomeron.
This is in sharp contrast to the picture above. Hence, this may indicate that
$1^{--}$ belongs to a subleading (rather than leading) trajectory. A possible
scenario is that the leading trajectory could be that of the $3^{--}$ state,
as pointed out in [80]. This scenario though is beyond the scope of the
holographic approximation discussed in this note. In other words, the
structure of the leading odderon trajectory is an interesting open question,
which will hopefully be addressed experimentally.

## Acknowledgements.

The authors thank I. Klebanov, C. Rodrigues Filho and C. Royon for
discussions. AD is supported by the NSF under grant PHY-2013812. The work of
DM was done with support of the Simons Foundation, award #884966, via the
Association International Institute of Physics (AIIF) and of the grant of the
agency CNPq of the Brazilian Ministry of Science, Technology and Innovation
#433935/2018-9.

This note is an extended version of a contribution to proceedings of the Low X
2021 conference on Elba, based on a talk given at that event.

## References

  * [1] TOTEM, D0 collaboration, V. M. Abazov et al., _Odderon Exchange from Elastic Scattering Differences between $pp$ and $p\bar{p}$ Data at 1.96 TeV and from pp Forward Scattering Measurements_, _Phys. Rev. Lett._ 127 (2021) 062003 [2012.03981].
  * [2] L. Lukaszuk and B. Nicolescu, _A Possible interpretation of p p rising total cross-sections_ , _Lett. Nuovo Cim._ 8 (1973) 405.
  * [3] Particle Data Group collaboration, P. A. Zyla et al., _Review of Particle Physics_ , _PTEP_ 2020 (2020) 083C01.
  * [4] J. Kang and M. A. Luty, _Macroscopic Strings and ’Quirks’ at Colliders_ , _JHEP_ 11 (2009) 065 [0805.4642].
  * [5] J. E. Juknevich, D. Melnikov and M. J. Strassler, _A Pure-Glue Hidden Valley I. States and Decays_ , _JHEP_ 07 (2009) 055 [0903.0883].
  * [6] J. E. Juknevich, _Pure-glue hidden valleys through the Higgs portal_ , _JHEP_ 08 (2010) 121 [0911.5616].
  * [7] C. J. Morningstar and M. J. Peardon, _The Glueball spectrum from an anisotropic lattice study_ , _Phys. Rev. D_ 60 (1999) 034509 [hep-lat/9901004].
  * [8] C. Michael and M. Teper, _The Glueball Spectrum in SU(3)_ , _Nucl. Phys. B_ 314 (1989) 347.
  * [9] C. Michael, _The Glueball Spectrum From Lattice Gauge Theory_ , _Acta Phys. Polon. B_ 21 (1990) 119.
  * [10] UKQCD collaboration, G. S. Bali, K. Schilling, A. Hulsebos, A. C. Irving, C. Michael and P. W. Stephenson, _A Comprehensive lattice study of SU(3) glueballs_ , _Phys. Lett. B_ 309 (1993) 378 [hep-lat/9304012].
  * [11] Y. Chen et al., _Glueball spectrum and matrix elements on anisotropic lattices_ , _Phys. Rev._ D73 (2006) 014516 [hep-lat/0510074].
  * [12] H. B. Meyer, _Glueball regge trajectories_ , Oxford PhD thesis, 2004.
  * [13] A. Athenodorou and M. Teper, _The glueball spectrum of SU(3) gauge theory in 3 + 1 dimensions_ , _JHEP_ 11 (2020) 172 [2007.06422].
  * [14] TXL, T(X)L collaboration, G. S. Bali, B. Bolder, N. Eicker, T. Lippert, B. Orth, P. Ueberholz et al., _Static potentials and glueball masses from QCD simulations with Wilson sea quarks_ , _Phys. Rev. D_ 62 (2000) 054503 [hep-lat/0003012].
  * [15] UKQCD collaboration, A. Hart and M. Teper, _On the glueball spectrum in O(a) improved lattice QCD_ , _Phys. Rev. D_ 65 (2002) 034502 [hep-lat/0108022].
  * [16] UKQCD collaboration, A. Hart, C. McNeile, C. Michael and J. Pickavance, _A Lattice study of the masses of singlet 0++ mesons_ , _Phys. Rev. D_ 74 (2006) 114504 [hep-lat/0608026].
  * [17] UKQCD collaboration, C. M. Richards, A. C. Irving, E. B. Gregory and C. McNeile, _Glueball mass measurements from improved staggered fermion simulations_ , _Phys. Rev. D_ 82 (2010) 034501 [1005.2473].
  * [18] E. Gregory, A. Irving, B. Lucini, C. McNeile, A. Rago, C. Richards et al., _Towards the glueball spectrum from unquenched lattice QCD_ , _JHEP_ 10 (2012) 170 [1208.1858].
  * [19] W. Sun, L.-C. Gui, Y. Chen, M. Gong, C. Liu, Y.-B. Liu et al., _Glueball spectrum from $N_{f}=2$ lattice QCD study on anisotropic lattices_, _Chin. Phys. C_ 42 (2018) 093103 [1702.08174].
  * [20] R. Brett, J. Bulava, D. Darvish, J. Fallica, A. Hanlon, B. Hörz et al., _Spectroscopy From The Lattice: The Scalar Glueball_ , _AIP Conf. Proc._ 2249 (2020) 030032 [1909.07306].
  * [21] F. Chen, X. Jiang, Y. Chen, K.-F. Liu, W. Sun and Y.-B. Yang, _Glueballs at Physical Pion Mass_ , 2111.11929.
  * [22] M. J. Teper, _Glueball masses and other physical properties of SU(N) gauge theories in D = (3+1): A Review of lattice results for theorists_ , hep-th/9812187.
  * [23] B. Lucini, M. Teper and U. Wenger, _Glueballs and k-strings in SU(N) gauge theories: Calculations with improved operators_ , _JHEP_ 06 (2004) 012 [hep-lat/0404008].
  * [24] B. Lucini, A. Rago and E. Rinaldi, _Glueball masses in the large N limit_ , _JHEP_ 08 (2010) 119 [1007.3879].
  * [25] J. Holligan, E. Bennett, D. K. Hong, J.-W. Lee, C. J. D. Lin, B. Lucini et al., _$Sp(2N)$ Yang-Mills towards large $N$_, _PoS_ LATTICE2019 (2019) 177 [1912.09788].
  * [26] E. Bennett, J. Holligan, D. K. Hong, J.-W. Lee, C. J. D. Lin, B. Lucini et al., _Color dependence of tensor and scalar glueball masses in Yang-Mills theories_ , _Phys. Rev. D_ 102 (2020) 011501 [2004.11063].
  * [27] E. Bennett, J. Holligan, D. K. Hong, J.-W. Lee, C. J. D. Lin, B. Lucini et al., _Glueballs and strings in $Sp(2N)$ Yang-Mills theories_, _Phys. Rev. D_ 103 (2021) 054509 [2010.15781].
  * [28] A. Athenodorou and M. Teper, _SU(N) gauge theories in 3+1 dimensions: glueball spectrum, string tensions and topology_ , _JHEP_ 12 (2021) 082 [2106.00364].
  * [29] J. Sexton, A. Vaccarino and D. Weingarten, _Scalar glueball decay_ , _Nucl. Phys. B Proc. Suppl._ 42 (1995) 279 [hep-lat/9412042].
  * [30] J. Sexton, A. Vaccarino and D. Weingarten, _Coupling constants for scalar glueball decay_ , _Nucl. Phys. B Proc. Suppl._ 47 (1996) 128 [hep-lat/9602022].
  * [31] N. Yamanaka, H. Iida, A. Nakamura and M. Wakayama, _Glueball scattering cross section in lattice SU(2) Yang-Mills theory_ , _Phys. Rev. D_ 102 (2020) 054507 [1910.07756].
  * [32] F. J. Llanes-Estrada, _Glueballs as the Ithaca of meson spectroscopy_ , 2101.05366.
  * [33] B. Lucini and M. Teper, _SU(N) gauge theories in four-dimensions: Exploring the approach to N = infinity_ , _JHEP_ 06 (2001) 050 [hep-lat/0103027].
  * [34] B. Lucini and M. Panero, _SU(N) gauge theories at large N_ , _Phys. Rept._ 526 (2013) 93 [1210.4997].
  * [35] B. Lucini and M. Panero, _Introductory lectures to large- $\scriptsize{N}$ QCD phenomenology and lattice results_, _Prog. Part. Nucl. Phys._ 75 (2014) 1 [1309.3638].
  * [36] S. Okubo, _Phi meson and unitary symmetry model_ , _Phys. Lett._ 5 (1963) 165.
  * [37] G. Zweig, _An SU(3) model for strong interaction symmetry and its breaking. Version 2_.  2, 1964.
  * [38] J. Iizuka, _Systematics and phenomenology of meson family_ , _Prog. Theor. Phys. Suppl._ 37 (1966) 21.
  * [39] J. M. Maldacena, _The Large N limit of superconformal field theories and supergravity_ , _Int. J. Theor. Phys._ 38 (1999) 1113 [hep-th/9711200].
  * [40] S. S. Gubser, I. R. Klebanov and A. M. Polyakov, _Gauge theory correlators from noncritical string theory_ , _Phys. Lett._ B428 (1998) 105 [hep-th/9802109].
  * [41] E. Witten, _Anti-de Sitter space and holography_ , _Adv. Theor. Math. Phys._ 2 (1998) 253 [hep-th/9802150].
  * [42] C. Csaki, H. Ooguri, Y. Oz and J. Terning, _Glueball mass spectrum from supergravity_ , _JHEP_ 01 (1999) 017 [hep-th/9806021].
  * [43] R. C. Brower, S. D. Mathur and C.-I. Tan, _Glueball spectrum for QCD from AdS supergravity duality_ , _Nucl. Phys._ B587 (2000) 249 [hep-th/0003115].
  * [44] E. Witten, _Anti-de Sitter space, thermal phase transition, and confinement in gauge theories_ , _Adv. Theor. Math. Phys._ 2 (1998) 505 [hep-th/9803131].
  * [45] H. Boschi-Filho, N. R. F. Braga and H. L. Carrion, _Glueball Regge trajectories from gauge/string duality and the Pomeron_ , _Phys. Rev. D_ 73 (2006) 047901 [hep-th/0507063].
  * [46] J. Erlich, E. Katz, D. T. Son and M. A. Stephanov, _QCD and a holographic model of hadrons_ , _Phys. Rev. Lett._ 95 (2005) 261602 [hep-ph/0501128].
  * [47] J. Polchinski and M. J. Strassler, _Hard scattering and gauge / string duality_ , _Phys. Rev. Lett._ 88 (2002) 031601 [hep-th/0109174].
  * [48] J. Polchinski and M. J. Strassler, _Deep inelastic scattering and gauge / string duality_ , _JHEP_ 05 (2003) 012 [hep-th/0209211].
  * [49] S. J. Brodsky and G. F. de Teramond, _Hadronic spectra and light-front wavefunctions in holographic QCD_ , _Phys. Rev. Lett._ 96 (2006) 201601 [hep-ph/0602252].
  * [50] I. R. Klebanov and M. J. Strassler, _Supergravity and a confining gauge theory: Duality cascades and chi SB resolution of naked singularities_ , _JHEP_ 08 (2000) 052 [hep-th/0007191].
  * [51] I. R. Klebanov and E. Witten, _Superconformal field theory on three-branes at a Calabi-Yau singularity_ , _Nucl. Phys. B_ 536 (1998) 199 [hep-th/9807080].
  * [52] I. R. Klebanov and N. A. Nekrasov, _Gravity duals of fractional branes and logarithmic RG flow_ , _Nucl. Phys. B_ 574 (2000) 263 [hep-th/9911096].
  * [53] I. R. Klebanov and A. A. Tseytlin, _Gravity duals of supersymmetric SU(N) x SU(N+M) gauge theories_ , _Nucl. Phys. B_ 578 (2000) 123 [hep-th/0002159].
  * [54] M. J. Strassler, _The Duality cascade_ , in _Theoretical Advanced Study Institute in Elementary Particle Physics (TASI 2003): Recent Trends in String Theory_ , pp. 419–510, 5, 2005, hep-th/0505153, DOI.
  * [55] S. S. Gubser, C. P. Herzog and I. R. Klebanov, _Symmetry breaking and axionic strings in the warped deformed conifold_ , _JHEP_ 09 (2004) 036 [hep-th/0405282].
  * [56] S. S. Gubser, C. P. Herzog and I. R. Klebanov, _Variations on the warped deformed conifold_ , _Comptes Rendus Physique_ 5 (2004) 1031 [hep-th/0409186].
  * [57] M. K. Benna, A. Dymarsky and I. R. Klebanov, _Baryonic Condensates on the Conifold_ , _JHEP_ 08 (2007) 034 [hep-th/0612136].
  * [58] I. R. Klebanov, P. Ouyang and E. Witten, _A Gravity dual of the chiral anomaly_ , _Phys. Rev. D_ 65 (2002) 105007 [hep-th/0202056].
  * [59] R. Apreda, _Nonsupersymmetric regular solutions from wrapped and fractional branes_ , hep-th/0301118.
  * [60] D. Cassani and A. F. Faedo, _A Supersymmetric consistent truncation for conifold solutions_ , _Nucl. Phys. B_ 843 (2011) 455 [1008.0883].
  * [61] X. Amador and E. Caceres, _Spin two glueball mass and glueball regge trajectory from supergravity_ , _JHEP_ 11 (2004) 022 [hep-th/0402061].
  * [62] M. Berg, M. Haack and W. Mueck, _Bulk dynamics in confining gauge theories_ , _Nucl. Phys. B_ 736 (2006) 82 [hep-th/0507285].
  * [63] A. Y. Dymarsky and D. G. Melnikov, _On the glueball spectrum in the Klebanov-Strassler model_ , _JETP Lett._ 84 (2006) 368.
  * [64] M. Berg, M. Haack and W. Mueck, _Glueballs vs. Gluinoballs: Fluctuation Spectra in Non-AdS/Non-CFT_ , _Nucl. Phys. B_ 789 (2008) 1 [hep-th/0612224].
  * [65] A. Dymarsky and D. Melnikov, _Gravity Multiplet on KS and BB Backgrounds_ , _JHEP_ 05 (2008) 035 [0710.4517].
  * [66] M. K. Benna, A. Dymarsky, I. R. Klebanov and A. Solovyov, _On Normal Modes of a Warped Throat_ , _JHEP_ 06 (2008) 070 [0712.4404].
  * [67] A. Dymarsky, D. Melnikov and A. Solovyov, _I-odd sector of the Klebanov-Strassler theory_ , _JHEP_ 05 (2009) 105 [0810.5666].
  * [68] I. Gordeli and D. Melnikov, _On I-even Singlet Glueballs in the Klebanov-Strassler Theory_ , _JHEP_ 08 (2011) 082 [0912.5517].
  * [69] I. Gordeli and D. Melnikov, _Calculation of glueball spectra in supersymmetric theories via holography_ , in _International Workshop on Low X Physics_ , 11, 2013, 1311.6537.
  * [70] D. Melnikov and C. R. Filho, _Pseudoscalar glueballs in the Klebanov-Strassler Theory_ , _JHEP_ 01 (2021) 024 [2007.11952].
  * [71] E. Caceres and C. Nunez, _Glueballs of super Yang-Mills from wrapped branes_ , _JHEP_ 09 (2005) 027 [hep-th/0506051].
  * [72] D. Elander, _Glueball Spectra of SQCD-like Theories_ , _JHEP_ 03 (2010) 114 [0912.1600].
  * [73] M. Bianchi and W. de Paula, _On Exact Symmetries and Massless Vectors in Holographic Flows and other Flux Vacua_ , _JHEP_ 04 (2010) 113 [1003.2536].
  * [74] D. Elander, _Light scalar from deformations of the Klebanov-Strassler background_ , _Phys. Rev. D_ 91 (2015) 126012 [1401.3412].
  * [75] D. Elander and M. Piai, _Calculable mass hierarchies and a light dilaton from gravity duals_ , _Phys. Lett. B_ 772 (2017) 110 [1703.09205].
  * [76] D. Elander and M. Piai, _Glueballs on the Baryonic Branch of Klebanov-Strassler: dimensional deconstruction and a light scalar particle_ , _JHEP_ 06 (2017) 003 [1703.10158].
  * [77] C. Rodrigues Filho, _Glueballs in the Klebanov–Strassler Theory: Pseudoscalars vs Scalars_ , _Braz. J. Phys._ 51 (2021) 788 [2011.12689].
  * [78] C. Bachas and J. Estes, _Spin-2 spectrum of defect theories_ , _JHEP_ 06 (2011) 005 [1103.2800].
  * [79] A. Dymarsky, I. R. Klebanov and N. Seiberg, _On the moduli space of the cascading SU(M+p) x SU(p) gauge theory_ , _JHEP_ 01 (2006) 155 [hep-th/0511254].
  * [80] H. B. Meyer and M. J. Teper, _Glueball Regge trajectories and the pomeron: A Lattice study_ , _Phys. Lett. B_ 605 (2005) 344 [hep-ph/0409183].
  * [81] G. A. Jaroszkiewicz and P. V. Landshoff, _Model for diffraction excitation_ , _Phys. Rev. D_ 10 (1974) 170.
  * [82] A. Donnachie and P. V. Landshoff, _Total cross-sections_ , _Phys. Lett. B_ 296 (1992) 227 [hep-ph/9209205].
  * [83] A. B. Kaidalov and Y. A. Simonov, _Glueball spectrum and the Pomeron in the Wilson loop approach_ , _Phys. Atom. Nucl._ 63 (2000) 1428 [hep-ph/9911291].
  * [84] A. B. Kaidalov and Y. A. Simonov, _Glueball masses and Pomeron trajectory in nonperturbative QCD approach_ , _Phys. Lett. B_ 477 (2000) 163 [hep-ph/9912434].
  * [85] C. Ewerz, _The Odderon in quantum chromodynamics_ , hep-ph/0306137.
  * [86] R. C. Brower, M. S. Costa, M. Djurić, T. Raben and C.-I. Tan, _Strong Coupling Expansion for the Conformal Pomeron/Odderon Trajectories_ , _JHEP_ 02 (2015) 104 [1409.2730].