# Link Prediction with Mutual Attention for Text-Attributed Networks

Robin Brochier  Université de Lyon, Lyon 2, ERIC EA3083 robin.brochier@univ-
lyon2.fr ,  Adrien Guille  Université de Lyon, Lyon 2, ERIC EA3083
<EMAIL_ADDRESS>and  Julien Velcin  Université de Lyon, Lyon 2,
ERIC EA3083<EMAIL_ADDRESS>
(2019)

###### Abstract.

In this extended abstract, we present an algorithm that learns a similarity
measure between documents from the network topology of a structured corpus. We
leverage the Scaled Dot-Product Attention, a recently proposed attention
mechanism, to design a mutual attention mechanism between pairs of documents.
To train its parameters, we use the network links as supervision. We provide
preliminary experiment results with a citation dataset on two prediction
tasks, demonstrating the capacity of our model to learn a meaningful textual
similarity.

representation learning; link prediction; attributed network; natural language
processing

††journalyear: 2019††copyright: iw3c2w3††conference: Companion Proceedings of
the 2019 World Wide Web Conference; May 13–17, 2019; San Francisco, CA,
USA††booktitle: Companion Proceedings of the 2019 World Wide Web Conference
(WWW ’19 Companion), May 13–17, 2019, San Francisco, CA, USA††doi:
10.1145/3308560.3316587††isbn: 978-1-4503-6675-5/19/05††ccs: Computing
methodologies Unsupervised learning††ccs: Computing methodologies Artificial
intelligence††ccs: Information systems Data mining

##  1\. Related Works

In this section, we relate recent works in the fields of network embedding
(NE) and attention mechanism for natural language processing (NLP).

###  1.1. Attributed Network Embedding

DeepWalk (Perozzi et al., 2014) first proposed to derive the word embedding
algorithm Word2vec (Mikolov et al., 2013) by generating paths of nodes, akin
to sentences, with truncated random walks. DeepWalk and other variants are
generalized into a common matrix factorization framework in NetMF (Qiu et al.,
2018). To extend DeepWalk for text-attributed networks, TADW (Yang et al.,
2015) expresses this latter as a matrix factorization problem and incorporates
a matrix of textual features $T$, produced by latent semantic indexing (LSI),
into the factorization so that the vertex similarity matrix can be
reconstructed as the product of three matrices $V^{T}$, $H$ and $T$.

###  1.2. Attention Mechanisms for NLP

The Transformer (Vaswani et al., 2017) is a novel neural architecture that
outperforms state-of-the-art methods in neural machine translation (NMT)
without the use of convolution nor recurrent units. The Scaled Dot-Product
Attention (SDPA) is the main constituting part of the Transformer that
actually performs attention over a set of words. It takes as input a query
vector $q$ and a set of key vectors $K$ of dimensions $d_{k}$ and value
vectors $V$ of dimensions $d_{v}$. One weight for a value is generated by a
compatibility function with its corresponding key and the query. Formally, the
attention vectors are generated in parallel for multiples queries $Q$,
following the formula:
$\text{Attention}(Q,K,V)=\text{softmax}(\frac{QK^{T}}{\sqrt{d_{k}}})V$. The
result is a set of $L$ attention vectors ($L$ being the number of queries) of
dimension $d_{v}$. The matrices $Q$, $K$ and $V$ are produced by projection of
initial words representations $W$ with three matrices $P^{Q}$, $P^{K}$ and
$P^{V}$ whose parameters are meant to be learned.

Several works (Devlin et al., 2018; Radford et al., 2018) adapted the
Transformer architecture beyond the task of NMT. Their main idea is to train
the Transformer in an unsupervised fashion over large corpora of texts and
further refine its parameters on a wide variety of supervised tasks. Motivated
by these recent works, we present a model, MATAN (Mutual Attention for Text-
Attributed Networks), that derives from the the SDPA to address the task of
link prediction in a network of documents.

##  2\. Proposed Model

We propose an algorithm for link prediction in text-attributed network. Our
model is trained under a NE procedure, presented in Section 2.1. The
optimization of the reconstruction error is performed via dot-product between
contextual document representations $e^{v}_{u}$ and $e^{u}_{v}$. These
embeddings are generated with a mutual attention mechanism over their textual
contents only, described in Section 2.2.

###  2.1. Overall Optimization

The model takes as input a network of documents $G=(V,E,T)$, $T$ being the
textual content of the documents. We precompute word embeddings $W$ and a
normalized similarity measure between nodes $M$ designed from the adjacency
matrix $A$ of the network. Each document $t_{u}$ is associated with a bag of
word embeddings $W^{t_{u}}$ matrix. For any pair of node $(u,v)\in V^{2}$,
mutual embeddings are generated with an asymmetric mutual attention function
$f^{A}_{\Theta}$ for both documents given their bags of word embeddings
$e^{v}_{u}=f^{A}_{\Theta}(W^{t_{u}},W^{t_{v}})$ and
$e^{u}_{v}=f^{A}_{\Theta}(W^{t_{v}},W^{t_{u}})$. We define the unormalized
similarity between the two nodes as the dot product of their mutual embeddings
$e^{v}_{u}.e^{u}_{v}$. We aim at learning the parameters $\Theta$ by
minimizing the KL divergence from the similarity distributions $M$ (from the
graph) to that of the normalized distribution of the dot products between the
mutual embeddings (Tsitsulin et al., 2018) (text associated to the nodes). We
achieve this by employing noise-contrastive estimation (Tsitsulin et al.,
2018), minimizing the following objective function: $J=-\sum_{(u,v)\in
C}\Big{(}\log\sigma(e^{v}_{u}\cdot e^{u}_{v})+\sum_{i=1}^{k}\mathbb{E}_{z\sim
q}\big{[}\log\sigma(-e^{z}_{u}\cdot e^{u}_{z})\big{]}\Big{)}$, where $\sigma$
is the sigmoid function. $C$ is a corpus of pairs of nodes generated by
drawing uniformly existing links from the empirical distribution of links $M$.
k negative nodes are uniformly drawn for each positive pair. To minimize this
objective, we employ stochastic gradient descent using ADAM (Kingma and Ba,
2014).

###  2.2. Mutual Attention Mechanism

The role of $f^{A}_{\Theta}(W^{t_{u}},W^{t_{v}})$ is to generate a contextual
representation of $t_{u}$ given $t_{v}$. The parameters $\Theta$ we aim to
learn are composed of three matrices $\Theta=\\{P^{Q},P^{K},P^{V}\\}$ of
dimension $D\times D$ each. For all words of the target document $t_{u}$, we
create queries $Q_{u}=W^{t_{u}}P^{Q}$. We similarly create keys and values
from the contextual document $t_{v}$, such that $K_{v}=W^{t_{v}}P^{Q}$ and
$V_{v}=W^{t_{v}}P^{V}$. Attention representations for each target word are
then computed, following the SDPA formula:
$\text{SDPA}_{\Theta}(W^{t_{u}},W^{t_{v}})=\text{softmax}(\frac{Q_{u}K_{v}^{T}}{\sqrt{D}})V_{v}$.
Note that $\text{SDPA}_{\Theta}(W^{t_{u}},W^{t_{v}})$ has dimension $L\times
D$, that is, we have a mutual attention representation of each word of
document $t_{u}$ given $t_{v}$. Finally, the representation for document
$t_{u}$ is obtained by averaging its word mutual attention vectors:
$e^{v}_{u}=f^{A}_{\Theta}(W^{t_{u}},W^{t_{v}})=\sum_{i=0}^{L}\text{SDPA}_{\Theta}(W^{t_{u}},W^{t_{v}})_{i}$.
Similarly, $e^{u}_{v}$ is generated by flipping indices $u$ and $v$. The
intuition behind this model is that the matrices $P^{Q}$ and $P^{K}$ learn to
project pairs of words that explain links in the network such that their dot-
products produce large weights. $P^{V}$ is then meant to project the word
vectors such that their average produces similar representations for nodes
that are close in the network and dissimilar for nodes that are far in the
network.

##  3\. Experiments

To assess the quality of our model, we perform two tasks of link prediction on
a dataset of citation links between scientific abstracts: Cora 111Get the
data: https://linqs.soe.ucsc.edu/data. The first prediction evaluation, called
edges-hidden, consists in hiding a percentage of the links given a network of
documents and measuring the ability of the model to predict higher scores to
hidden links than to non-existing ones by computing the ROC AUC. The second
evaluation, called nodes-hidden, consists in splitting the network into two
unconnected networks, keeping a percentage of the nodes in the training
network.

We precompute on the full corpus word embeddings using GloVe (Pennington et
al., 2014) of dimension 256 with a co-occurrence threshold
$x_{\text{max}}=10$, a window size $w=5$ and 50 epochs. We precompute LSI
(Deerwester et al., 1990) vectors of dimension 128. For the edge-hidden
prediction task, we provide results performed by NetMF with $k=10$ negative
samples. TADW is run with $20$ epochs and MATAN is performed with $k=1$
negative sample and $10^{5}$ sampled pairs of documents. The empirical
similarity matrix between the nodes we chose is the normalized adjacency
matrix. All produced representations are of dimension 256.

###  3.1. Results

Table 1. Edges-hidden link prediction ROC AUC % of training data | 10% | 20% | 30% | 40% | 50%
---|---|---|---|---|---
NeMF | 59.0 | 67.2 | 77.5 | 83.2 | 87.2
TADW | 68.0 | 82.0 | 87.1 | 93.2 | 94.5
MATAN | 82.3 | 87.1 | 88.6 | 90.9 | 91.0
Table 2. Nodes-hidden link prediction ROC AUC % of training data | 10% | 20% | 30% | 40% | 50%
---|---|---|---|---|---
TADW | 64.2 | 75.8 | 80.3 | 81.9 | 82.3
MATAN | 69.4 | 73.0 | 75.4 | 77.9 | 78.6

Tables 1 and 2 show the results of our experiments. MATAN shows promising
results for learning on a small percentage of training data on both
evaluations. TADW has better scores for nodes-hidden predictions which might
be explained by the capacity of LSI to learn discriminant features on a small
dataset unlike GloVe. In future work we would like to deal with bigger
datasets from which word embedding methods might capture richer semantic
information.

## References

  * (1)
  * Deerwester et al. (1990) Scott Deerwester, Susan T. Dumais, George W. Furnas, Thomas K. Landauer, and Richard Harshman. 1990\.  Indexing by latent semantic analysis.  _JOURNAL OF THE AMERICAN SOCIETY FOR INFORMATION SCIENCE_ 41, 6 (1990), 391–407.
  * Devlin et al. (2018) Jacob Devlin, Ming-Wei Chang, Kenton Lee, and Kristina Toutanova. 2018.  Bert: Pre-training of deep bidirectional transformers for language understanding.  _arXiv preprint arXiv:1810.04805_ (2018).
  * Kingma and Ba (2014) Diederik P Kingma and Jimmy Ba. 2014.  Adam: A method for stochastic optimization.  _arXiv preprint arXiv:1412.6980_ (2014).
  * Mikolov et al. (2013) Tomas Mikolov, Kai Chen, Greg Corrado, and Jeffrey Dean. 2013\.  Efficient estimation of word representations in vector space.  _arXiv preprint arXiv:1301.3781_ (2013).
  * Pennington et al. (2014) Jeffrey Pennington, Richard Socher, and Christopher Manning. 2014.  Glove: Global vectors for word representation. In _Proceedings of the 2014 conference on empirical methods in natural language processing (EMNLP)_. 1532–1543.
  * Perozzi et al. (2014) Bryan Perozzi, Rami Al-Rfou, and Steven Skiena. 2014\.  Deepwalk: Online learning of social representations. In _Proceedings of the 20th ACM SIGKDD international conference on Knowledge discovery and data mining_. ACM, 701–710.
  * Qiu et al. (2018) Jiezhong Qiu, Yuxiao Dong, Hao Ma, Jian Li, Kuansan Wang, and Jie Tang. 2018\.  Network embedding as matrix factorization: Unifying deepwalk, line, pte, and node2vec. In _Proceedings of the Eleventh ACM International Conference on Web Search and Data Mining_. ACM, 459–467.
  * Radford et al. (2018) Alec Radford, Karthik Narasimhan, Tim Salimans, and Ilya Sutskever. 2018.  Improving language understanding by generative pre-training.  _URL https://s3-us-west-2. amazonaws. com/openai-assets/research-covers/language-unsupervised/language_ understanding_paper. pdf_ (2018).
  * Tsitsulin et al. (2018) Anton Tsitsulin, Davide Mottin, Panagiotis Karras, and Emmanuel Müller. 2018.  VERSE: Versatile Graph Embeddings from Similarity Measures. In _Proceedings of the 2018 World Wide Web Conference on World Wide Web_. International World Wide Web Conferences Steering Committee, 539–548.
  * Vaswani et al. (2017) Ashish Vaswani, Noam Shazeer, Niki Parmar, Jakob Uszkoreit, Llion Jones, Aidan N Gomez, Łukasz Kaiser, and Illia Polosukhin. 2017.  Attention is all you need. In _Advances in Neural Information Processing Systems_. 5998–6008.
  * Yang et al. (2015) Cheng Yang, Zhiyuan Liu, Deli Zhao, Maosong Sun, and Edward Y Chang. 2015.  Network representation learning with rich text information.. In _IJCAI_. 2111–2117.


Version 3 as of

Primary authors: MIT, NTU, INFN Bari, GSU, CU Boulder, CERN

MITHIG-MOD-21-001

# Jet energy spectrum and substructure in $e^{+}e^{-}$ collisions at 91.2 GeV
with ALEPH Archived Data

Yi Chen<EMAIL_ADDRESS>Massachusetts Institute of Technology, Cambridge,
Massachusetts, USA  Anthony Badea  Harvard University, Cambridge,
Massachusetts, USA  Austin Baty  Rice University, Houston, Texas, USA  Paoti
Chang  National Taiwan University, Taipei, Taiwan  Yang-Ting Chien  Georgia
State University, Atlanta, Georgia, USA  Gian Michele Innocenti  CERN, Geneva,
Switzerland  Marcello Maggi  INFN Sezione di Bari, Bari, Italy  Christopher
McGinn  University of Colorado, Boulder, USA  Dennis V. Perepelitsa
University of Colorado, Boulder, USA  Michael Peters  Massachusetts Institute
of Technology, Cambridge, Massachusetts, USA  Tzu-An Sheng  Massachusetts
Institute of Technology, Cambridge, Massachusetts, USA  Jesse Thaler
Massachusetts Institute of Technology, Cambridge, Massachusetts, USA  Yen-Jie
Lee<EMAIL_ADDRESS>Massachusetts Institute of Technology, Cambridge,
Massachusetts, USA

###### Abstract

The first measurements of energy spectra and substructure of anti-$k_{T}$ jets
in hadronic $Z^{0}$ decays in $e^{+}e^{-}$ collisions are presented. The
archived $e^{+}e^{-}$ annihilation data at a center-of-mass energy of 91.2 GeV
were collected with the ALEPH detector at LEP in 1994. In addition to
inclusive jet and leading dijet energy spectra, various jet substructure
observables are analyzed as a function of jet energy which includes groomed
and ungroomed jet mass to jet energy ratios, groomed momentum sharing, and
groomed jet radius. The results are compared with perturbative QCD
calculations and predictions from the Sherpa, Herwig v7.1.5, Pythia6, Pythia8
and Pyquen event generators. The jet energy spectra agree with perturbative
QCD calculations which include the treatment of logarithms of the jet radius
and threshold logarithms. None of the event generators give a fully
satisfactory description of the data.

##  I Introduction

Jets, sprays of particles originating from the fragmentation and hadronization
of a fast-moving quark or gluon, are some of the most useful tools for
studying the Quantum Chromodynamics (QCD) and searching for new physics beyond
the Standard Model (SM) in high energy colliders. Since the end of the Large
Electron Positron Collider (LEP) operation in 2000, significant progress has
been made on jet definitions and jet clustering algorithms. Information about
the parton shower is not only carried by the four-momentum of jets but also by
their internal structure. Therefore, active developments on jet substructure
observables have been carried out in order to fully exploit jets Abdesselam et
al. (2011); Salam (2010); Altheimer et al. (2012); Larkoski et al. (2020).
However, those jet reconstruction and jet substructure algorithms, which are
widely used in the data analyses of proton-proton Kogler et al. (2019) and
heavy-ion collisions Connors et al. (2018); Cao and Wang (2021); Cunqueiro and
Sickles (2022), have not been used in the most elementary electron-positron
$e^{+}e^{-}$ collision system. Phenomenological models, such as Pythia
Sjöstrand et al. (2001), Sherpa Gleisberg et al. (2009), and Herwig Reichelt
et al. (2017), are tuned with hadron level or hadronic event shape observables
from LEP experiments and are employed to predict or to describe the jet
spectra and jet substructure in more complicated hadron-hadron collision
environments.

Studies of jets, using identical clustering algorithms as those used in high
energy hadron colliders such as the Large Hadron Collider (LHC) and
Relativistic Heavy-Ion Collider (RHIC), are of great interest. Unlike hadron-
hadron collisions, $e^{+}e^{-}$ annihilations do not have beam remnants,
gluonic initial state radiation, or the complications of parton distribution
functions. The $e^{+}e^{-}$ data provide the cleanest test for analytical QCD
calculations and phenomenological models that are tuned with hadronic event
shapes. Unlike the smooth and steeply falling jet transverse momentum spectra
at RHIC and LHC, the inclusive jet and leading jet energy spectra in
$e^{+}e^{-}$ annihilation are peaked near the half of the collision energy.
Therefore, leading and inclusive jet spectra are sensitive to the jet energy
loss from wide-angle radiation that is not clustered into the jet. These
spectra provide new tests of the calculation of leading-jet cross-sections,
which is sensitive to the treatment of logarithms of the jet radius and
threshold logarithms Neill et al. (2021). Moreover, fully reconstructed jets
provide us with an opportunity to inspect quark and gluon fragmentation in
great detail on a shower-by-shower basis. Finally, studies of jet substructure
observables and their comparison to modern event generators are of great
interest since jet substructure observables are novel tools for jet flavor
identification, electroweak boson and top tagging, and studies of the Quark-
Gluon Plasma properties at hadron colliders Andrews et al. (2020). Those
results can also serve as references for understanding proton-proton collision
data and future electron-ion collisions (EIC) Abdul Khalek et al. (2021).

In this paper, the first measurement of inclusive jet momentum spectra, jet
splitting functions, (groomed) jet mass, and groomed jet radius distribution
in hadronic $Z^{0}$ boson decays are presented, with two types of jet
selections. The inclusive observables include all jets above an energy
threshold and inside a defined acceptance. They are sensitive to higher-order
corrections of jet spectra in perturbative QCD. They include low momentum jets
in events that are typically associated with soft gluon radiation or random
combinatorial jets with hadrons from different partons. On the other hand, the
leading dijet observables consider the leading and subleading jet in the
event. This type of observable focuses more on the dominant energy flow and is
less sensitive to soft radiation.

##  II ALEPH Detector

The ALEPH detector is described in ref. Decamp et al. (1990). The central part
of the detector is designed for the efficient reconstruction of charged
particles. The trajectories of them are measured by a two-layer silicon strip
vertex detector, a cylindrical drift chamber, and a large time projection
chamber (TPC). Those tracking detectors are inside a 1.5 T axial magnetic
field generated by a superconducting solenoidal coil. The charged particle
transverse momenta are reconstructed with a resolution of $\delta p_{\rm
T}/p_{\rm T}=6\times 10^{-4}p_{\rm T}\oplus 0.005$ ($p_{\rm T}$ in GeV/c).

Electrons and photons are identified in the electromagnetic calorimeter (ECAL)
situated between the TPC and the superconducting coil. The ECAL is a sampling
calorimeter, made by lead plates and proportional wire chambers segmented in
$0.9^{\circ}\times 0.9^{\circ}$ projective towers. They are read out in three
sections in-depth and have a total thickness of around 22 radiation lengths.
Isolated photons are reconstructed with a relative energy resolution of
$0.18/\sqrt{E}+0.009$ ($E$ in GeV).

The iron return yoke constructed with 23 layers of streamer tubes is also used
as the hadron calorimeter (HCAL) for the detection of charged and neutral
hadrons. The relative energy resolution for hadrons is $0.85/\sqrt{E}$ ($E$ in
GeV). Muons are identified by their pattern in HCAL and by the muon chambers,
made by two double-layers of streamer tubes outside the HCAL.

The information from trackers and calorimeters is combined in an energy flow
algorithm Buskulic et al. (1995). This particle flow algorithm outputs a set
of charged and neutral particles, called the energy flow objects. They are
used in the jet reconstruction in this analysis.

##  III Dataset and Event Selection

The archived $e^{+}e^{-}$ annihilation data at a center-of-mass energy of 91.2
GeV were collected with the ALEPH detector at LEP in 1994. To analyze these
data, an MIT Open Data format Tripathee et al. (2017) was created and
validated in the two-particle correlation function analysis Badea et al.
(2019). Hadronic events are selected by requiring the sphericity axis Heister
et al. (2004) to have a polar angle in the laboratory reference frame
($\theta_{\text{lab}}$) between $7\pi/36$ and $29\pi/36$ to ensure that the
event is well contained within the detector. At least five tracks having a
minimum total energy of 15 GeV are also required to suppress electromagnetic
interactions Barate et al. (1998). The residual contamination from processes
such as $e^{+}e^{-}\rightarrow\tau^{+}\tau^{-}$ is expected to be less than
0.26% for these event selections Barate et al. (1998). Approximately 1.36
million $e^{+}e^{-}$ collisions resulting in the decay of a $Z^{0}$ boson to
quarks are analyzed.

High-quality tracks from particles are selected using requirements identical
to those in previous ALEPH analyses Barate et al. (1998) and are also required
to have a transverse momentum with respect to the beam axis ($p_{\rm T}^{\rm
lab}$) above 0.2 GeV/c and $|\cos{\theta_{\text{lab}}}|<0.94$ in the lab
frame. Secondary charged particles from neutral particle decays are suppressed
by $V^{0}$ reconstruction in the energy flow algorithm Barate et al. (1998).
In addition, the total visible energy from energy flow objects is required to
be smaller than 200 GeV in order to remove laser calibration events during the
run.

Event thrust distributions Farhi (1977) published by the ALEPH Collaboration
using a similar dataset Heister et al. (2004) were successfully reproduced
within uncertainties, affirming that the archived data is analyzed properly.

##  IV Simulated Samples

Archived Pythia 6.1 Sjöstrand et al. (2001) Monte Carlo (MC) simulation
samples which were produced with the detector conditions during the 1994 run
are used to correct for detector effects and to compare with the data. A set
of new Pythia events are also generated with Pythia version 8.303 Sjöstrand et
al. (2015) at a center-of-mass energy of 91.2 GeV. The Monash 2013 tune Skands
et al. (2014) is used, with the weak boson exchange and weak single boson
processes turned on. In order to select hadronic $Z^{0}$ decays, pure
electroweak events are rejected by requiring at last one hadron in the final
state particles. Sherpa samples are produced with version 2.2.5 Gleisberg et
al. (2009), with electron positron events generating 2 to 5 outgoing partons,
which are then showered into jets. The coupling constant $\alpha_{s}(M_{Z})$
is set to 0.1188 in the event generation. A set of Herwig Bellm et al. (2016)
samples is generated using version 7.2.2. In this sample, 2 to 3 outgoing
partons are specified in the hard process, and leptonic decays of $Z^{0}$
boson are turned off to increase the fraction of hadronic events. The order of
$\alpha$ coupling is set to 2, whereas the colored $\alpha_{s}$ order is set
to 0 in the sample generation. In order to demonstrate the potential
modification of the jet spectra and jet substructure observables from jet
quenching Bjorken (1982); Aad et al. (2010); Chatrchyan et al. (2011), a
sample is generated with the Pyquen generator Lokhtin and Snigirev (2006)
(version 1.5.3). The strength of jet quenching is set to be equivalent to that
in a minimum bias sample of PbPb collisions at 5.02 TeV and wide-angle
radiation is turned off during the sample generation.

##  V Jet reconstruction and performance

Experimental jets are clustered with the anti-kT algorithm Cacciari et al.
(2008) with a distance parameter $R=0.4$, using the FastJet package Cacciari
et al. (2012) (version 3.3.2). This resolution parameter was chosen since it
is widely used in the jet analyses in proton-proton and heavy-ion collisions
Aad et al. (2013); Aaboud et al. (2019); Sirunyan et al. (2021); Khachatryan
et al. (2017); Dokshitzer et al. (1997) carried out at the hadron colliders.
Moreover, the chosen value also gives us an opportunity to examine the shower
from quarks in detail. The distance measure $d_{ij}$ between jets $i$ and $j$
and $d_{iB}$ between jet $i$ and the beam is defined using the opening angle
$\theta_{ij}$:

$\displaystyle d_{ij}$
$\displaystyle=\min(E_{i}^{-2},E_{j}^{-2})\dfrac{1-\cos\theta_{ij}}{1-\cos R}$
$\displaystyle d_{iB}$ $\displaystyle=E_{i}^{-2},$

where $E_{i}$ is the energy of the $i$-th jet and $\theta_{ij}$ is the opening
angle of the $i$-th and $j$-th jets. The termination of the clustering process
is defined by the distance between the jet and the beam direction, $d_{iB}$.
In the data analysis, energy flow objects which are reconstructed using the
tracker and calorimeter information, are used for jet reconstruction. The mass
of the objects are assumed to be pion for charged hadron candidates, and
massless for photon candidates. Generator-level jets are clustered by
considering all visible final state particles by the ALEPH detector (i.e.,
excluding neutrinos). Reconstructed-level jets are clustered with all energy
flow candidates reconstructed by the tracker and calorimeters. In order to
avoid overlapping with the beam pipe, only jets between
$0.2\pi<\theta_{\text{jet}}<0.8\pi$ are considered in this measurement, where
$\theta_{\text{jet}}$ is the angle between jet momentum and beam axis.

The jets are calibrated following a multi-stage strategy. In the first stage,
simulated jet responses are corrected by using generated jets, clustered with
visible final state particles with a distance parameter of $R=0.4$, as the
reference. The correction factor is derived as a function of jet 3-momentum
$p_{\text{jet}}$ in 18 bins of jet direction $\theta_{\text{jet}}$. The
closure in the simulated sample is seen to be better than 0.5%.

The jets in data are corrected with the first stage correction from simulation
before deriving residual corrections characterizing the difference between
data and simulation. In the second stage, there are two steps in the residual
correction: relative residual correction, which equalizes the jet response
between the two halves of the detector, and absolute residual correction,
which calibrates the jet to an absolute scale. The jet relative correction is
performed using the leading jet distribution in different bins of
$\theta_{\text{jet}}$ and comparing the two sides of the detector. This
procedure aims to unify the jet energy scale in data between the two sides of
the detector. The correction is typically of scale up to 0.5% in data, while
the result on the simulated sample with the same procedure is consistent with
0.

A cross-check is carried out using the momentum balance of dijets in bins of
third jet activity related to the leading dijet: $E_{3}/(E_{1}+E_{2})$. The
residual scale is derived by extrapolating the third jet activity to zero. Due
to the larger uncertainty induced by the extrapolation procedure, this
procedure is not adopted as nominal. The corrections between the cross-check
and the nominal one are consistent within the uncertainty.

The absolute correction of jets is performed using multi-jet events with an
invariant mass around the $Z^{0}$ peak. Both the simulation correction and the
residual relative correction are applied before the derivation of the absolute
correction. The overall scale is fitted using up to $N$ high energy jets (with
$N$ ranging from 3 to 9), and with different requirements on the energy of the
$(N+1)$-th leading jet from 1 to 5 GeV as a handle to control potential bias
from the method. If there are less than $N$ jets reconstructed in the event,
the $(N+1)$-th jet has 0 energy and the event is included. The correction is
derived in bins of $\theta_{\text{jet}}$ and is up to 0.3%.

The jet energy resolution in simulation is measured using an empirical
function, in bins of $\theta_{\text{jet}}$:

$\displaystyle\sigma_{\text{jet}}=\sqrt{a_{0}+\dfrac{a_{1}}{p_{\text{jet}}}+\dfrac{a_{2}}{p_{\text{jet}}^{2}}+\dfrac{a_{3}}{p_{\text{jet}}^{3}}},$

where $a_{0}$ to $a_{3}$ are parameters to be fitted. Depending on the
direction of the jet, $a_{0}$ ranges from 0.005–0.010 GeV2, $a_{1}$ from
0.58–0.72 GeV3, $a_{2}$ up to 0.62 GeV4, while $a_{3}$ is negligible.

The difference between jet energy resolution between data and simulation is
measured using the momentum balance of the leading dijet in an event. The
magnitude of the third-leading jet is used as a handle to understand potential
systematic bias from the method, and it is constrained to be less than 3 GeV
in the final result. The difference is characterized as a scale factor, in
bins of $\theta_{\text{jet}}$. The resolution in data is found to be up to 5%
worse than that in simulation, depending on the direction of the jet Chen et
al. (2021).

##  VI Analysis

###  VI.1 Jet Grooming

The soft drop algorithm Larkoski et al. (2014); Dasgupta et al. (2013) is
proposed as a way to define a set of theoretically well-controlled observables
that are sensitive to structures inside a jet and has been used in various QCD
measurements in proton-proton and heavy-ion collisions Sirunyan et al.
(2018a); Acharya et al. (2020); Adam et al. (2020); Aad et al. (2020);
Sirunyan et al. (2018b, c). It is characterized by two free parameters
$z_{\text{cut}}$ and $\beta$. The constituents of a jet are first reclustered
using the Cambridge-Aachen algorithm Dokshitzer et al. (1997); Wobisch and
Wengler (1998), and then the clustering history is traced, and the algorithm
terminates when the soft drop condition is satisfied:

$\displaystyle z\equiv\dfrac{\min(E_{1},E_{2})}{E_{1}+E_{2}}\geq
z_{\text{cut}}\left(\dfrac{\theta_{12}}{R}\right)^{\beta},$

where $E_{1}$ and $E_{2}$ are the energies of the two branches, and
$\theta_{12}$ is their opening angle. Note the definition is adapted to using
the energy and opening angle (a true Great-Circle Distance), as opposed to
using the transverse energy and $\Delta
R\equiv\sqrt{(\Delta\eta)^{2}+(\Delta\phi)^{2}}$ as its counterpart in a
hadronic collider environment. When the algorithm terminates, the $z$ is
defined as $z_{G}$ Larkoski et al. (2017), and the opening angle is defined as
$R_{G}$.

In this measurement, in addition to the jet energy spectrum, the $z_{G}$ and
$R_{G}$ spectra in bins of jet energy for $z_{\text{cut}}=0.1,\beta=0$ are
reported. For the $z_{G}$ and $R_{G}$ distributions, the distributions are
self-normalized in order to decouple effects coming from jet energy migration
and the rest.

###  VI.2 Unfolding

In order to mitigate smearing effects from the detector, an unfolding is
performed on the measured spectra. The unfolding is done with a Bayes’
iterative method based on ref. D’Agostini (1995), as implemented in the
RooUnfold package. The method applies an iterative method with a corrected
error calculation which decreases sensitivity to the prior knowledge used in
the Bayes’ formulation.

The ideal number of iterations for the unfolding is obtained with a two-step
procedure. The generator level distribution from the simulation is first
combined with the smearing matrix to produce the idealized detector level
distribution. Pseudodatasets are then generated using this distribution with
the number of statistics compatible with that of the data. The ideal number of
iterations can then be obtained by comparing the unfolded distribution with
the generator level distribution. In the second step, data is first unfolded
using the number of iterations found from the first step, and the whole
procedure is repeated by using the unfolded distribution in place of the
generator level distribution. The numbers range from 7 to 17 and the exact
numbers can be found in ref. Chen et al. (2021).

The response matrix for the jet energy measurement is shown in Figure 1. The
measurements are done in bins of 1 GeV up to 45 GeV, after which the bin sizes
are progressively larger in order to accumulate enough entries in the tail. In
Figure 1, the matrix is row-normalized: the color scale indicates the
probability density of a generated (“Gen”) jet of given energy being smeared
into different reconstructed (“Reco”) jet energy. A tight correlation between
the reconstructed jet and the matched truth-level jet is observed. No
significant contribution from outliers is present in the correlation matrix.
The final result for this measurement starts from 10 GeV based on studies of
jet energy scales. The bins below 10 GeV are included in the unfolding
procedure as underflow bins, but they are not part of the final result.

Figure 1: Response matrix for jet $E$. The matrix is row-normalized: the color
scale indicates the probability density of a given generator level jet being
smeared into a different reconstructed level jet.

Due to non-negligible migration of jet energy, the substructure observables
are unfolded together with the jet energy $E$: $(E,z_{G})$, $(E,R_{G})$,
$(E,M/E)$ and $(E,M_{G}/E)$. The jet energy is binned in 5 GeV bins until the
kinematic limit, after which the bin size is increased due to the available
number of entries. Part of the response matrices are shown in Figs. 2, 3, 4
and 5. The full matrices are too large to be properly included in the main
paper body. The matrices are row-normalized. For most of the matrix, within
each block of jet energy, the response is mostly diagonal, and there are
significant contributions in the off-diagonal blocks due to jet energy
migration. For $z_{G}$ and $R_{G}$, there is an off-diagonal contribution
within each block at a few percent level. These correspond to cases where the
detector smearing causes the soft drop algorithm to terminate at a different
point in the declustering history. At very low energy below 10 GeV (not shown
in the partial matrices), the response is observed to be less diagonal. This
indicates that there might be a significant contribution of jets of
combinatorial origin.

Figure 2: Partial response matrix for $z_{G}$, in bins of jet $E$. Each block
is the sub-matrix for the indicated jet $E$ range. The matrix is row-
normalized across all jet $E$ ranges: the color scale indicates the
probability density of a given generator level (“Gen”) jet being smeared into
a different reconstructed level (“Reco”) jet. Figure 3: Partial response
matrix for $R_{G}$, in bins of jet $E$. Each block is the sub-matrix for the
indicated jet $E$ range. The matrix is row-normalized across all jet $E$
ranges: the color scale indicates the probability density of a given generator
level (“Gen”) jet being smeared into a different reconstructed level (“Reco”)
jet. Figure 4: Partial response matrix for $M_{G}/E$, in bins of $E$. Each
block is the sub-matrix for the indicated $E$ range. The matrix is row-
normalized across all $E$ ranges: the color scale indicates the probability
density of a given generator level (“Gen”) jet being smeared into a different
reconstructed level (“Reco”) jet. Figure 5: Partial response matrix for $M/E$,
in bins of $E$. Each block is the sub-matrix for the indicated $E$ range. The
matrix is row-normalized across all $E$ ranges: the color scale indicates the
probability density of a given generator level (“Gen”) jet being smeared into
a different reconstructed level (“Reco”) jet.

The statistical uncertainty from the unfolding machinery is validated using
pseudoexperiments. A set of 250 pseudodatasets are generated by combining the
unfolded data distribution and the smearing matrix. Each of the pseudodatasets
is unfolded, and the pull distribution between the unfolded pseudodataset and
the input distribution is examined. For all cases, the width of the pull
distribution is consistent with 1.

###  VI.3 Leading dijet selection

There is no detector acceptance close to the beamline, and therefore in order
to measure only the event-wide leading and subleading jet, referred to as
“leading dijet” in the later text, an additional selection is designed. The
total energy visible either inside the acceptance $0.2\pi<\theta<0.8\pi$ or
within radius $R$ of any jet above 5 GeV is calculated, and events with total
visible energy above 83 GeV are selected to ensure a purity of 99% of the
events having leading dijet inside the acceptance. The efficiency of such a
selection is about 60%. A correction is derived using simulated events to
correct for effects on the spectra due to this extra selection.

##  VII Systematic uncertainties

The majority of systematic uncertainties are evaluated by unfolding the data
with alternative versions of smearing matrices that encompass the systematic
variations we would like to examine, with a few exceptions, as discussed
individually below. Different sources are added in quadrature. These sources
are discussed in the following paragraphs:

Systematics related to hadronic event selections

The requirement on the number of charged particles passing all track
selections is changed from 5 to 6 in order to check the systematics related to
the hadronic event selection. The full analysis chain is rerun with the new
selection, and the results are compared to the nominal selection. The
deviation is negligible ($<0.1\%$).

Systematics related to jet reconstruction

The residual jet energy scale (“JEC”) is varied by $\pm 0.5\%$ which is the
maximum magnitude of the residual corrections. This uncertainty source ranges
from a few percent up to 20% in jet energy spectra, and it is the dominating
source of uncertainty at large $E_{\text{jet}}$. For jet substructure results,
due to the self-normalization, jet energy scale uncertainty, which manifests
itself mostly as an overall normalization in each energy bin, is small.

The jet energy is smeared by $\pm 2.5\%$, the difference between data and MC
simulation extracted from dijet balance studies, to construct alternative
response matrices in order to evaluate the uncertainty associated with the
resolution (“JER”) difference between data and simulation. It is the
dominating source in the mid-energy region and ranges from a few percent up to
15%.

The corrections to detector effects are extracted from the archived Monte
Carlo. We study the uncertainty coming from this choice of simulation with the
following procedure. We first reweight the simulated sample to match
predictions on certain observables from other event generators, such as Herwig
and Sherpa, before building the smearing matrix and repeating the whole
analysis. We considered the following target observables in the reweighting
studies: (1) reweight the jet constituent multiplicity with different energy
thresholds on the constituents (0.5 GeV, 1.0 GeV, and 1.5 GeV); (2) reweight
MC to match the substructure observable $z_{G}$ and $R_{G}$. The unfolded
results with reweighted MC changed by up to a few percent, and the maximum
difference from the nominal result is quoted as the systematic uncertainty.

For the angular and energy resolution of the subjets, the associated
uncertainties are determined in the following way. The directions of the
subjet momenta are smeared randomly based on the observed simulated
resolution. Although this choice is conservative, the overall effect is small.
The energy of the subjets is also smeared based on twice the nominal energy
resolution scale factor (+5%) obtained for the jets. The total uncertainty
from these two sources amounts to 3% for $z_{G}$ and $R_{G}$, and 5% for jet
mass. The total uncertainties from substructure modeling and MC reweighting
are shown in the “Model” category in the summary plots.

Contributions from fake or combinatorial jets are estimated by matching
reconstructed jets in simulation to generator-level jets. The fraction of
unmatched jets indicates the amount of potential fake jets. It ranges from
0–5% depending on the jet energy for inclusive jets.

A check on the results has been carried out by comparing the measured spectrum
from each side of the detector ($0.2\pi<\theta_{\text{jet}}<0.5\pi$ and
$0.5\pi<\theta_{\text{jet}}<0.8\pi$). The spectra are consistent with each
other.

Systematics related to unfolding

We examine the sensitivity to the choice of prior which is part of the input
to the unfolding. The difference between a flat prior and a prior with the
generator level distribution from the simulation is taken as the uncertainty.
It is grouped (“Unfold”) with the uncertainty associated with the number of
iterations done in the unfolding process. The number of iterations is varied
from the determined ideal number of iterations to that optimized on a separate
sample (ie., through Monte Carlo studies), and the difference in the unfolded
distribution is taken as uncertainty.

In order to assess for potential bias in the unfolding method itself, an
alternative unfolding method, the singular value decomposition (SVD) method
Hocker and Kartvelishvili (1996), is used. The SVD unfolding method employs a
regularization parameter, and the ideal parameter is determined using a method
similar to that used to determine the ideal number of iterations in the
nominal unfolding method. The difference between the SVD and the nominal is
quoted as a systematic uncertainty.

A final source of systematic uncertainty related to the unfolding is the
observed nonclosure on simulated events. This source of uncertainty is
negligible for jet $E$ and sizable, but still subdominant, for $z_{G}$ and
$R_{G}$ when $E$ is small.

Systematics related to leading dijet selection

The uncertainty on the extra leading dijet selection is evaluated by changing
the thresholds on the total visible energy to those corresponding to a leading
dijet purity of 98% and 99.5%. The difference to nominal is taken as the
systematic uncertainty. The effect from potential imperfect modeling of the
resolution of the total visible energy is investigated by smearing the
quantity event by event and repeating the analysis. The difference is added in
quadrature to the scale uncertainty.

Uncertainty that arises from the correction factor is evaluated through a
reweighting procedure on the simulated sample. The events are weighted by the
ratio of jet spectra for each jet between the nominal unfolded data and
simulation, and the weighted simulated sample is used to derive the
alternative correction factor. The difference is quoted as uncertainty.

A summary of the systematic uncertainties are shown in Fig. 6 for jet $E$,
Figs. 7 for leading dijet, and Figs. 8, 9, 11, and 10 for the jet substructure
observables.

Figure 6: Summary of the size of relative systematic uncertainty from
different sources for jet $E$.

Figure 7: Summary of the size of relative systematic uncertainty from
different sources for each of the two jet $E$ in leading dijet (left) and
leading dijet total energy (right). Figure 8: Summary of the size of relative
systematic uncertainty from different sources for $z_{G}$, in bins of jet $E$
from 10–15 GeV (left panel) to above 40 GeV (right panel). Figure 9: Summary
of the size of relative systematic uncertainty from different sources for
$R_{G}$, in bins of jet $E$ from 10–15 GeV (left panel) to above 40 GeV (right
panel). Figure 10: Summary of the size of relative systematic uncertainty from
different sources for $M_{G}/E$, in bins of jet $E$ from 10–15 GeV (left
panel) to above 40 GeV (right panel). Figure 11: Summary of the size of
relative systematic uncertainty from different sources for $M/E$, in bins of
jet $E$ from 10–15 GeV (left panel) to above 40 GeV (right panel).

##  VIII Results and discussions

The fully corrected jet spectra and jet substructure observables are presented
and compared with predictions from Pythia6 (red), Pythia8 (blue), Herwig
(green), and Sherpa (purple) event generators. The results are also compared
to analytical calculations with perturbative QCD. Finally, predictions from
the PYQUEN (gray) event generator, which added jet quenching effect to the
simulated $e^{+}e^{-}$ events, are also overlapped to illustrate the possible
modifications due to the presence of a strongly interacting medium.

Figure 12: Measured inclusive jet $E$ spectrum. The spectrum is normalized by
the number of events used to perform the measurement. Any overflow
contribution is included in the final bin. (Left Panel) The data is compared
with predictions from Pythia6 (red), Pythia8 (blue), Herwig (green), and
Sherpa (purple) generators. The predictions are normalized to have the same
area as the data. (Right Panel) The data spectra are compared with pQCD
calculations at NLO and NLL’+R resummation Neill et al. (2021).

The first unfolded inclusive jet energy spectrum is shown in Figure 12. The
distribution is peaked at around half of the $Z^{0}$ boson mass. This
structure is mainly coming from $Z^{0}\rightarrow q\bar{q}$, with the parton
shower of one of the outgoing (anti-)quark is almost fully captured by the
anti-$k_{T}$ algorithm with a resolution parameter of 0.4.

The spectrum decreases rapidly as one moves to lower jet energy, reaching a
minimum at around 20 GeV. Below that, the spectrum increases as we go to even
lower energy due to the contributions from soft radiation and combinatorial
jets. These features in the jet energy spectrum are captured by most of the
event generators. Pythia6, which was tuned to describe ALEPH data, gives the
best description of the energy spectrum. At low jet energy, Herwig gives the
worse description of the results. It over-predicts the jet spectrum at low jet
energy. Pyquen generator, which includes a large jet quenching effect
equivalent to that in the inelastic PbPb collisions, predicts a large
reduction of the population at around 45 GeV and a significant increase in the
number of jets at low jet energy.

The same data could also be compared to a next-to-leading-order perturbative
QCD calculation and a calculation based on a next-to-leading logarithmic
(NLL’) threshold and jet radius resummation (NLL’+R) Neill et al. (2021). The
NLO calculation predicts a sharper peak at large jet energy. The NLL’+R
calculation gives a better description of the inclusive jet energy spectrum,
which shows the importance to include the jet radius resummation.

Figure 13: Measured leading dijet jet $E$ spectrum. The spectrum is normalized
by the number of events used to perform the measurement. Any overflow
contribution is included in the final bin. The data is also compared with
predictions from Pythia6 (red), Pythia8 (blue), Herwig (green), and Sherpa
(purple) generators. The predictions are normalized to have the same area as
the data. (Right Panel) The data spectra are compared with pQCD calculations
at NLL’+R resummation Neill et al. (2021) Figure 14: Measured leading dijet
sum $E$ spectrum. The spectrum is normalized by the number of events used to
perform the measurement. Any overflow contribution is included in the final
bin. The data is also compared with predictions from Pythia6 (red), Pythia8
(blue), Herwig (green), and Sherpa (purple) generators. The predictions are
normalized to have the same area as the data.

In order to focus on the dominant energy flow and suppress the contributions
from soft radiation and combinatorial jets, the leading dijet energy is shown
in Figure 13. Due to the leading dijet selection, the rise at low jet energy
is suppressed. The results are consistent with predictions from Pythia6 within
the quoted systematical uncertainties. Pythia8, Herwig, Sherpa tend to
overpredict the spectrum at low jet energy. The dijet sum energy (Figure 14),
which is equivalent to the collision energy minus the sum of all the lower
energy jets, is more sensitive to the modeling of subleading jets. The levels
of (dis)agreement between the simulation and data for leading dijet energy and
the leading dijet sum energy are similar.

Figure 15: Measured $z_{G}$ spectra in bins of jet $E$. The spectrum is self-
normalized for each $E$ interval. The fraction of jets completely groomed away
are included in the first bin of each panel. The data is also compared with
predictions from Pythia6 (red), Pythia8 (blue), Herwig (green), and Sherpa
(purple) generators. Figure 16: Measured $R_{G}$ spectra in bins of jet $E$.
The spectrum is self-normalized for each $E$ interval. The fraction of jets
completely groomed away are included in the first bin of each panel. The data
is also compared with predictions from Pythia6 (red), Pythia8 (blue), Herwig
(green), and Sherpa (purple) generators.

Now we change our focus to the characterization of the substructure inside the
anti-$k_{T}$ jets. The groomed momentum sharing $z_{G}$ spectra and the
opening angle of subjets $R_{G}$ are presented in bins of jet energy. The
distributions are self-normalized for each energy bin to factor out effects
from an imperfect jet energy spectrum in the simulated samples. As shown in
Figure 15, the $z_{G}$ spectrum is falling as a function of $z_{G}$ value,
reaching a minimum at $z_{G}\sim 0.5$, which is similar to the data from
proton-proton and heavy-ion collisions Sirunyan et al. (2018a). The jets which
failed the grooming procedure are populated in the first bin $z_{G}=0$ in the
plot. At high jet energy $E>40$ GeV, Herwig predicts a flatter $z_{G}$ spectra
than data and over-predicts the jets with $z_{G}$ close to 0.5 while Pythia6
under-predicts the population at $z_{G}\sim 0.5$. This is similar to the
results from pp collisions at a higher jet transverse momentum Sirunyan et al.
(2018a). At low jet energy $E<40$ GeV, the data are consistent with
predictions from HERWIG 7. In most of the jet energy intervals, Pythia6,
Pythia8 and Sherpa under-predict the number of jets at large $z_{G}$ by up to
roughly 10% while they tend to over-predict the $z_{G}$ spectra at low
$z_{G}$.

The $R_{G}$ is the distance between the groomed subjets. The agreement between
$R_{G}$ spectra and event generators, as shown in Figure 16, is worse than
that observed in the $z_{G}$ spectra. At high jet energy, event generators
predict a slightly narrower core of the $R_{G}$ spectra compared to data. At
low jet energy, most event generators predict on average a larger separation
between subjects (larger $R_{G}$) than the unfolded data. The PYQUEN generator
that includes jet quenching effect predicts an even larger fraction of jets
with large $R_{G}$.

Figure 17: Measured $M/E$ and $M_{G}/E$ spectra in bins of jet $E$. The
spectrum is self-normalized for each $E$ interval. The fraction of jets
completely groomed away are included in the first bin of each panel. The data
is also compared with predictions from Pythia6 (red), Pythia8 (blue), Herwig
(green), and Sherpa (purple) generators.

The jet mass and groomed jet mass results are shown in Figure 17. A similar
trend to the $R_{G}$ result as a function of jet energy is observed: $M/E$ and
$M_{G}/E$ are generally smaller at larger jet energy compared to lower energy
results. In addition to the jet transverse width, the mass is also sensitive
to the scale of the initial parton. The difference in agreement between
$R_{G}$ and mass observables can shed light on the difference in agreement
between jet transverse width and the initial scale of the parton. The Pyquen
generator has a more pronounced effect on the mass which manifests as a shift
in the measured spectra, which is especially significant at low jet energy.
The difference between $M/E$ and $M_{G}/E$ is also of interest here, and it is
sensitive to the amount of wide-angle soft radiation in the jet. We leave the
explicit measurement of $\Delta M/E\equiv(M-M_{G})/E$ to future studies.

##  IX Summary

In summary, the first measurements of anti-$k_{T}$ jet energy spectrum and
substructure in hadronic $Z^{0}$ decays are presented. The archived
$e^{+}e^{-}$ annihilation data at a center-of-mass energy of 91.2 GeV were
collected with the ALEPH detector at LEP in 1994. Inclusive jet and dijet
energy spectra are presented using jets clustered with a distance parameter
$R=0.4$. In addition, jet substructure observables are analyzed as a function
of jet energy which includes groomed and ungroomed jet mass to jet energy
ratio, groomed momentum sharing, and groomed jet radius. The results are
compared with the perturbative QCD calculations and predictions from the
Sherpa, Herwig v7.1.5, Pythia6, Pythia8, and Pyquen event generators in order
to test the Monte Carlo modeling of initial state radiation, final state
radiation, parton shower, and hadronization effects. While the inclusive jet
and dijet energy spectra agree better with the predictions from the archived
Pythia6 MC, which was tuned to describe the ALEPH data, none of the event
generators gives a satisfactory description of the data. While the inclusive
jet spectrum can not be described by parton-level perturbative QCD calculation
at next-to-leading order, calculations based on next-to-leading logarithmic
threshold and jet radius resummation could give a good description of the
leading and inclusive jet energy spectra. These results provide the cleanest
test of the perturbative QCD calculation of jets and serve as an important
reference to the studies of jets in proton-proton, heavy-ion, and future
electron-ion collisions.

###### Acknowledgements.

The authors would like to thank the ALEPH Collaboration for their support and
foresight in archiving their data. We thank Joao Pires for providing the NLO
calculations of the jet energy spectrum. We appreciate the help from Felix
Ringer for providing the calculations with NLL’+R resummation. We would like
to thank the useful comments and suggestions from Roberto Tenchini, Guenther
Dissertori, Andrew Lakoski, Liliana Apolinario, Ben Nachman, Camelia Mironov
and Jing Wang. This work has been supported by the Department of Energy,
Office of Science, under Grant No. DE-SC0011088 (to Y.L., Y.C., M.P. and
T.S.), Grant No. DE-SC0012567 (to J.T.), Grant No. DE-SC0005131 (to A. Baty),
DOE Grants No. DE-SC0018117 and DE-FG02-03ER41244, and the Research
Corporation for Scientific Advancement (to D.V.P. and C.M.), the Ministry of
Education and the Ministry of Science and Technology of Taiwan (to P.C.), the
Harvard Frederick Sheldon Traveling Fellowship (to A. Badea).

## References

  * Abdesselam et al. (2011) A. Abdesselam et al., Eur. Phys. J. C 71, 1661 (2011), eprint 1012.5412.
  * Salam (2010) G. P. Salam, Eur. Phys. J. C 67, 637 (2010), eprint 0906.1833.
  * Altheimer et al. (2012) A. Altheimer et al., J. Phys. G 39, 063001 (2012), eprint 1201.0008.
  * Larkoski et al. (2020) A. J. Larkoski, I. Moult, and B. Nachman, Phys. Rept. 841, 1 (2020), eprint 1709.04464.
  * Kogler et al. (2019) R. Kogler et al., Rev. Mod. Phys. 91, 045003 (2019), eprint 1803.06991.
  * Connors et al. (2018) M. Connors, C. Nattrass, R. Reed, and S. Salur, Rev. Mod. Phys. 90, 025005 (2018), eprint 1705.01974.
  * Cao and Wang (2021) S. Cao and X.-N. Wang, Rept. Prog. Phys. 84, 024301 (2021), eprint 2002.04028.
  * Cunqueiro and Sickles (2022) L. Cunqueiro and A. M. Sickles, Prog. Part. Nucl. Phys. 124, 103940 (2022), eprint 2110.14490.
  * Sjöstrand et al. (2001) T. Sjöstrand, P. Eden, C. Friberg, L. Lonnblad, G. Miu, S. Mrenna, and E. Norrbin, Comput. Phys. Commun. 135, 238 (2001), eprint hep-ph/0010017.
  * Gleisberg et al. (2009) T. Gleisberg, S. Hoeche, F. Krauss, M. Schonherr, S. Schumann, F. Siegert, and J. Winter, JHEP 02, 007 (2009), eprint 0811.4622.
  * Reichelt et al. (2017) D. Reichelt, P. Richardson, and A. Siodmok, Eur. Phys. J. C77, 876 (2017), eprint 1708.01491.
  * Neill et al. (2021) D. Neill, F. Ringer, and N. Sato, JHEP 07, 041 (2021), eprint 2103.16573.
  * Andrews et al. (2020) H. A. Andrews et al., J. Phys. G 47, 065102 (2020), eprint 1808.03689.
  * Abdul Khalek et al. (2021) R. Abdul Khalek et al. (2021), eprint 2103.05419.
  * Decamp et al. (1990) D. Decamp et al. (ALEPH), Nucl. Instrum. Meth. A294, 121 (1990), [Erratum: Nucl. Instrum. Meth.A303,393(1991)].
  * Buskulic et al. (1995) D. Buskulic et al. (ALEPH), Nucl. Instrum. Meth. A 360, 481 (1995).
  * Tripathee et al. (2017) A. Tripathee, W. Xue, A. Larkoski, S. Marzani, and J. Thaler, Phys. Rev. D96, 074003 (2017), eprint 1704.05842.
  * Badea et al. (2019) A. Badea, A. Baty, P. Chang, G. M. Innocenti, M. Maggi, C. Mcginn, M. Peters, T.-A. Sheng, J. Thaler, and Y.-J. Lee, Phys. Rev. Lett. 123, 212002 (2019), eprint 1906.00489.
  * Heister et al. (2004) A. Heister et al. (ALEPH), Eur. Phys. J. C35, 457 (2004).
  * Barate et al. (1998) R. Barate et al. (ALEPH), Phys. Rept. 294, 1 (1998).
  * Farhi (1977) E. Farhi, Phys. Rev. Lett. 39, 1587 (1977).
  * 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, Comput. Phys. Commun. 191, 159 (2015), eprint 1410.3012.
  * Skands et al. (2014) P. Skands, S. Carrazza, and J. Rojo, Eur. Phys. J. C 74, 3024 (2014), eprint 1404.5630.
  * Bellm et al. (2016) J. Bellm et al., Eur. Phys. J. C76, 196 (2016), eprint 1512.01178.
  * Bjorken (1982) J. D. Bjorken (1982), eprint FERMILAB-PUB-82-059-THY.
  * Aad et al. (2010) G. Aad et al. (ATLAS), Phys. Rev. Lett. 105, 252303 (2010), eprint 1011.6182.
  * Chatrchyan et al. (2011) S. Chatrchyan et al. (CMS), Phys. Rev. C 84, 024906 (2011), eprint 1102.1957.
  * Lokhtin and Snigirev (2006) I. P. Lokhtin and A. M. Snigirev, Eur. Phys. J. C 45, 211 (2006), eprint hep-ph/0506189.
  * Cacciari et al. (2008) M. Cacciari, G. P. Salam, and G. Soyez, JHEP 04, 063 (2008), eprint 0802.1189.
  * Cacciari et al. (2012) M. Cacciari, G. P. Salam, and G. Soyez, Eur. Phys. J. C 72, 1896 (2012), eprint 1111.6097.
  * Aad et al. (2013) G. Aad et al. (ATLAS), Phys. Lett. B 719, 220 (2013), eprint 1208.1967.
  * Aaboud et al. (2019) M. Aaboud et al. (ATLAS), Phys. Lett. B 790, 108 (2019), eprint 1805.05635.
  * Sirunyan et al. (2021) A. M. Sirunyan et al. (CMS), JHEP 05, 284 (2021), eprint 2102.13080.
  * Khachatryan et al. (2017) V. Khachatryan et al. (CMS), Phys. Rev. C 96, 015202 (2017), eprint 1609.05383.
  * Dokshitzer et al. (1997) Y. L. Dokshitzer, G. D. Leder, S. Moretti, and B. R. Webber, JHEP 08, 001 (1997), eprint hep-ph/9707323.
  * Chen et al. (2021) Y. Chen, Y.-J. Lee, M. Maggi, P. Chang, Y.-T. Chien, C. McGinn, and D. Perepelitsa (2021), eprint 2108.04877.
  * Larkoski et al. (2014) A. J. Larkoski, S. Marzani, G. Soyez, and J. Thaler, JHEP 05, 146 (2014), eprint 1402.2657.
  * Dasgupta et al. (2013) M. Dasgupta, A. Fregoso, S. Marzani, and G. P. Salam, JHEP 09, 029 (2013), eprint 1307.0007.
  * Sirunyan et al. (2018a) A. M. Sirunyan et al. (CMS), Phys. Rev. Lett. 120, 142302 (2018a), eprint 1708.09429.
  * Acharya et al. (2020) S. Acharya et al. (ALICE), Phys. Lett. B 802, 135227 (2020), eprint 1905.02512.
  * Adam et al. (2020) J. Adam et al. (STAR), Phys. Lett. B 811, 135846 (2020), eprint 2003.02114.
  * Aad et al. (2020) G. Aad et al. (ATLAS), Phys. Rev. D 101, 052007 (2020), eprint 1912.09837.
  * Sirunyan et al. (2018b) A. M. Sirunyan et al. (CMS), Phys. Rev. D 98, 092014 (2018b), eprint 1808.07340.
  * Sirunyan et al. (2018c) A. M. Sirunyan et al. (CMS), JHEP 10, 161 (2018c), eprint 1805.05145.
  * Wobisch and Wengler (1998) M. Wobisch and T. Wengler, in _Workshop on Monte Carlo Generators for HERA Physics (Plenary Starting Meeting)_ (1998), pp. 270–279, eprint hep-ph/9907280.
  * Larkoski et al. (2017) A. Larkoski, S. Marzani, J. Thaler, A. Tripathee, and W. Xue, Phys. Rev. Lett. 119, 132003 (2017), eprint 1704.05066.
  * D’Agostini (1995) G. D’Agostini, Nucl. Instrum. Meth. A 362, 487 (1995).
  * Hocker and Kartvelishvili (1996) A. Hocker and V. Kartvelishvili, Nucl. Instrum. Meth. A 372, 469 (1996), eprint hep-ph/9509307.