CERN-TH-2022-217

# Deep Learning Searches for Vector-Like Leptons at the LHC
and Electron/Muon Colliders

António P. Morais<EMAIL_ADDRESS>Departamento de Física da Universidade de
Aveiro and Centre for Research and Development in Mathematics and Applications
(CIDMA), 3810-183 Aveiro, Portugal
Theoretical Physics Department, CERN, 1211 Geneva 23, Switzerland
António Onofre<EMAIL_ADDRESS>Departamento de Física da Universidade
do Minho, 4710-057 Braga, Portugal
Felipe F. Freitas<EMAIL_ADDRESS>Departamento de Física da Universidade
de Aveiro and Centre for Research and Development in Mathematics and
Applications (CIDMA), 3810-183 Aveiro, Portugal
João Gonçalves<EMAIL_ADDRESS>Departamento de Física da Universidade de
Aveiro and Centre for Research and Development in Mathematics and Applications
(CIDMA), 3810-183 Aveiro, Portugal
Roman Pasechnik<EMAIL_ADDRESS>Department of Astronomy and
Theoretical Physics, Lund University, Sölvegatan 14A, SE 223-62 Lund, Sweden
Rui Santos<EMAIL_ADDRESS>Centro de Física Teórica e Computacional,
Faculdade de Ciências, Universidade de Lisboa, 1749-016 Lisboa, Portugal
ISEL - Instituto Superior de Engenharia de Lisboa, Instituto Politécnico de
Lisboa, 1959-007 Lisboa, Portugal

###### Abstract

The discovery potential of both singlet and doublet vector-like leptons (VLLs)
at the Large Hadron Collider (LHC) as well as at the not-so-far future muon
and electron machines is explored. The focus is on a single production channel
for LHC direct searches while double production signatures are proposed for
the leptonic colliders. A Deep Learning algorithm to determine the discovery
(or exclusion) statistical significance at the LHC is employed. While doublet
VLLs can be probed up to masses of 1 TeV, their singlet counterparts have very
low cross sections and can hardly be tested beyond a few hundreds of GeV at
the LHC. This motivates a physics-case analysis in the context of leptonic
colliders where one obtains larger cross sections in VLL double production
channels, allowing to probe higher mass regimes otherwise inaccessible even to
the LHC high-luminosity upgrade.

Beyond Standard Model, $\mathrm{(g-2)}_{\mu}$ anomaly, Vector-like leptons,
Deep Learning, Lepton Colliders

##  I Introduction

The Standard Model (SM) of particle physics has, so far, served as the guide
towards unraveling the nature of all subatomic phenomena and its success is
now indisputable. However, it is clear that it lacks some of the necessary
ingredients to fully describe nature, such as an explanation for neutrino
masses, as firstly observed in neutrino oscillations by experiments like the
Super-Kamiokande Fukuda  _et al._ (1998) or a particle candidate to explain
experimental evidences for the existence of Dark Matter (DM) Bertone  _et al._
(2005). Recently, other discrepancies have also come into the forefront, such
as the measurement of the anomalous magnetic moment of the muon at Fermilab
Abi  _et al._ (2021), indicating a tension of $4.2\sigma$ with the SM
predictions, as well as the recent results coming from the LHCb experiment,
where hints of lepton flavour universality violation have been reported at a
$3.1\sigma$ significance Aaij  _et al._ (2021). Naturally, such deviations do
not meet the $5\sigma$ requirement to indicate the existence of new physics
(NP) and more data still needs to be collected if the experimental deviations
are confirmed. Together with the theoretical calculations in the framework of
the SM, an explanation of anomalies requires the addition of NP.

Finding well motivated NP scenarios that simultaneously address all the
aforementioned questions is of utmost importance and, in particular, Grand
Unified theories (GUTs) provide interesting avenues to follow. These classes
of models typically predict new particles present at the TeV scale (or not too
far from it), that may be probed at current and future colliders. Therefore,
the study of simplified frameworks that can later be mapped to the low-energy
limit of a certain UV complete framework can help in constraining their viable
parameter space. In particular, some of the authors have proposed a GUT
framework, which unifies all matter, Higgs and gauge sectors within the
$\mathrm{E}_{8}$ group, either at a conventional high-scale scenario Morais
_et al._ (2020); de Anda  _et al._ (2020), or a low-scale one Aranda  _et al._
(2021a, b). It turns out that, as a common prediction of any of such
proposals, the presence of new doublet Vector-Like Lepton (VLL) particles at
the TeV scale or below is expected. A phenomenological analysis of doublet
type VLLs, in the context of GUT model, was extensively discussed by some of
the authors in a previous work Freitas  _et al._ (2021). Naturally, such
description was bounded by the constraints and details of the underlying UV
theory and as such, the main VLL features can only be fully captured in more
simple phenomenological models, without the vast particle content and
constraints of GUT frameworks. Additionally, singlet-type VLLs, which are not
predicted in the context of the previously considered GUT, can also be
constructed and represent a significant par of this work. Studying the
differences between both types of VLL is of relevance not only in the context
of the LHC but also at future proposed leptonic colliders. Furthermore, if the
anomaly in the measurement of the magnetic moment of the muon is confirmed,
which can happen in about one year from the date this article is being
written, a possible explanation can arise in the form of TeV scale VLLs. In
this context, we refer to Yin and Yamaguchi (2020); Crivellin  _et al._
(2020); Bißmann  _et al._ (2021); Hiller  _et al._ (2020a, b); Guedes and
Santiago (2021) for further recent VLLs phenomenological studies.

With the above arguments in mind we study the phenomenological viability of
simple extensions of the SM enlarged with a single VLL and a right-handed
neutrino. For completeness, we consider two scenarios, the
$\mathrm{SU(2)_{L}}$ doublet and the singlet VLLs. Due to their vector-like
nature, they escape constraints from the fourth generation searches Eberhardt
_et al._ (2012), while non-vanishing couplings to the SM can induce
contributions not only to the muon $(g-2)_{\mu}$ anomaly, but also opening
tightly constrained channels such as $\mu\rightarrow e\gamma$ Zyla  _et al._
(2020).

This article is organized as follows. In Sec. II, we introduce the simplified
models, presenting the relevant theoretical details that will serve as the
basis behind the subsequent numerical simulations. In Sec. III we introduce
the numerical methodology for the studies we conduct in this work. In
particular, in Sec. III.1, we discuss the impact of the main flavour
constraints associated with the addition of a VLL, while in Sec. III.2 we
delve into collider analysis, both in the context of the LHC and future
electron (linear) or muon colliders. Finally, in Sec. IV we present our main
results and conclude in Sec. V.

##  II The models

In this article, we study collider phenomenology of both $\mathrm{SU(2)_{L}}$
doublet and singlet VLL extensions of the SM accompanied with new right-handed
neutrinos. We present the NP and SM fields’ quantum numbers in Tabs. 1 and 2
respectively.

Field | $\mathbf{SU(3)_{\text{C}}}$ | $\mathbf{SU(2)_{\text{L}}}$ | $\mathbf{U(1)_{\text{Y}}}$
---|---|---|---
$E_{\mathrm{L,R}}$ | 1 | 2 | $-1/2$
$\mathcal{E}_{\mathrm{L,R}}$ | 1 | 1 | $1$
$\nu_{\mathrm{R}}$ | 1 | 1 | $0$
Table 1:  Quantum numbers of the new exotic fermions. We adopt the notation of $E_{\mathrm{L,R}}$ being an $\mathrm{SU(2)}$ doublet, while $\mathcal{E}_{\mathrm{L,R}}$ is a singlet, and $\nu_{\mathrm{R}}$ is the new right-handed neutrino. Field | $\mathbf{SU(3)_{\text{C}}}$ | $\mathbf{SU(2)_{\text{L}}}$ | $\mathbf{U(1)_{\text{Y}}}$ | # of generations
---|---|---|---|---
$Q_{\mathrm{L}}$ | 3 | 2 | $1/6$ | 3
$L$ | 1 | 2 | $-1/2$ | 3
$d_{\mathrm{R}}$ | 3 | 1 | $-1/3$ | 3
$u_{\mathrm{R}}$ | 3 | 1 | $2/3$ | 3
$e_{\mathrm{R}}$ | 1 | 1 | $-1$ | 3
$\phi$ | $\bm{1}$ | $\bm{2}$ | $1/2$ | 1
Table 2:  Quantum numbers for the SM fields. Here, $Q_{\mathrm{L}}$,
$d_{\mathrm{R}}$, $u_{\mathrm{R}}$ are the quark fields, $L$ and
$e_{\mathrm{R}}$ are the lepton fields, and $\phi$ is the Higgs doublet.

The $\mathrm{SU(2)_{L}}$ doublets are defined as

$Q_{L}^{i}=\begin{bmatrix}u_{\mathrm{L}}\\\
d_{\mathrm{L}}\end{bmatrix}^{i},\quad L^{i}=\begin{bmatrix}\nu_{\mathrm{L}}\\\
e_{\mathrm{L}}\end{bmatrix}^{i},\quad
E_{\mathrm{L,R}}=\begin{bmatrix}\nu^{\prime}_{\mathrm{L,R}}\\\
e^{\prime}_{\mathrm{L,R}}\end{bmatrix},\quad\phi=\begin{bmatrix}\phi^{+}\\\
\phi^{0}\end{bmatrix},$ (1)

where $i$ is a generation index, $i=1,2,3$. The most general and
renormalizable Lagrangian density of the doublet VLL model, consistent with
the symmetries in Tab. 1, reads as

$\displaystyle\mathcal{L}_{d}=$
$\displaystyle\quantity(\Theta_{i}\bar{E}_{\mathrm{L}}e_{\mathrm{R}}^{i}\phi+\Upsilon\bar{E}_{\mathrm{L}}\nu_{\mathrm{R}}\tilde{\phi}+\Sigma_{i}\bar{L}^{i}\nu_{\mathrm{R}}\tilde{\phi}+\Omega\bar{E}_{\mathrm{R}}\nu_{R}\tilde{\phi}+\Pi_{ij}\bar{L}^{i}e_{\mathrm{R}}^{j}\phi+\mathrm{H.c.})+$
(2) $\displaystyle+M_{\rm E}\bar{E}_{\mathrm{L}}E_{\mathrm{R}}+M_{\rm
LE}^{i}\bar{L}_{i}E_{\mathrm{R}}+\frac{1}{2}M_{\nu_{\mathrm{R}}}\bar{\nu}_{R}\nu_{R}\,,$

where the Yukawa couplings are denoted by $\Theta$, $\Upsilon$, $\Sigma$,
$\Omega$ and $\Pi$, whereas bi-linear mass terms are indicated by $M_{E}$,
$M_{LE}$ and $M_{\nu_{\mathrm{R}}}$. For the singlet VLL case, the Lagrangian
is written as

$\displaystyle\mathcal{L}_{s}=$
$\displaystyle\quantity(\theta_{i}\bar{L}^{i}\mathcal{E}_{\mathrm{R}}\phi+\sigma^{i}\bar{L}_{i}\nu_{\mathrm{R}}\tilde{\phi}+\pi_{ij}\bar{L}^{i}e_{\mathrm{R}}^{j}\phi+\mathrm{H.c.})+M_{\rm
E}\bar{\mathcal{E}}_{\mathrm{L}}\mathcal{E}_{\mathrm{R}}+\frac{1}{2}M_{\nu_{\mathrm{R}}}\bar{\nu}_{R}\nu_{R}\,,$
(3)

where $\pi$, $\sigma$ and $\theta$ are the corresponding Yukawa couplings. The
mass matrices for both the singlet and the doublet scenarios, expressed in the
$\\{e_{L}^{1},e_{L}^{2},e_{L}^{2},e_{L}^{\prime}\\}\otimes\\{e_{R}^{1},e_{R}^{2},e_{R}^{2},e_{R}^{\prime}\\}$
basis, take the following form

$M_{L}^{\mathrm{doublet}}=\begin{bmatrix}\quantity(\dfrac{v\Pi}{\sqrt{2}})_{3\times
3}&\quantity(M_{\rm LE})_{3\times 1}\\\\[15.00002pt]
\quantity(\dfrac{v\Theta}{\sqrt{2}})_{1\times 3}&\quantity(M_{\rm E})_{1\times
1}\end{bmatrix},\quad
M_{L}^{\mathrm{singlet}}=\begin{bmatrix}\quantity(\dfrac{v\pi}{\sqrt{2}})_{3\times
3}&\quantity(\dfrac{v\theta}{\sqrt{2}})_{3\times 1}\\\\[15.00002pt]
\quantity(0)_{1\times 3}&\quantity(M_{\rm E})_{1\times 1}\end{bmatrix}\,,$ (4)

with $v=246\leavevmode\nobreak\ \mathrm{GeV}$ the electroweak symmetry
breaking (EWSB) Higgs doublet VEV $\left(0,\leavevmode\nobreak\
v\right)^{\top}$, and whose eigenstates describe both the chiral SM-like
charged leptons as well as the new VLLs.

For the neutrino sector the mass matrices can be expressed as

$\displaystyle M_{\nu}^{\mathrm{doublet}}=\begin{bmatrix}\quantity(0)_{3\times
3}&\quantity(\dfrac{\sigma v}{\sqrt{2}})_{3\times 1}&\quantity(0)_{3\times
1}&\quantity(M_{\rm LE})_{3\times 1}\\\ \quantity(\dfrac{\sigma
v}{\sqrt{2}})_{1\times 3}&M_{\nu_{\mathrm{R}}}&\quantity(\dfrac{\Upsilon
v}{\sqrt{2}})_{1\times 1}&\quantity(\dfrac{\Omega v}{\sqrt{2}})_{1\times 1}\\\
\quantity(0)_{1\times 3}&\quantity(\dfrac{\Upsilon v}{\sqrt{2}})_{1\times
1}&(0)_{1\times 1}&(M_{\rm E})_{1\times 1}&\\\ \quantity(M_{\rm LE})_{1\times
3}&\quantity(\dfrac{\Omega v}{\sqrt{2}})_{1\times 1}&(M_{\rm E})_{1\times
1}&\quantity(0)_{1\times 1}\end{bmatrix}\,,\qquad
M_{\nu}^{\mathrm{singlet}}=\begin{bmatrix}(0)_{3\times
3}&\quantity(\dfrac{\sigma v}{\sqrt{2}})_{3\times 1}\\\
\quantity(\dfrac{\sigma v}{\sqrt{2}})^{\dagger}_{1\times
3}&\quantity(M_{\nu_{\mathrm{R}}})_{1\times 1}\end{bmatrix}\,.$ (5)

where the
$\\{\nu_{L}^{1},\nu_{L}^{2},\nu_{L}^{3},\nu_{R},\nu_{L}^{{}^{\prime}},\nu_{R}^{{}^{\prime}}\\}\otimes\\{\nu_{L}^{1},\nu_{L}^{2},\nu_{L}^{3},\nu_{R},\nu_{L}^{{}^{\prime}},\nu_{R}^{{}^{\prime}}\\}$
basis was used for the doublet scenario whereas that of the singlet model was
chosen as
$\\{\nu_{L}^{1},\nu_{L}^{2},\nu_{L}^{3},\nu_{R}\\}\otimes\\{\nu_{L}^{1},\nu_{L}^{2},\nu_{L}^{3},\nu_{R}\\}$.

In this article, one considers that the lightest beyond-the-SM (BSM) neutrino
is in the KeV range and sterile, acting as missing energy in the detector.
Therefore, while the neutrino mass matrices above have the standard type-I
seesaw form, such a light sterile neutrino implies that $M_{\nu_{\mathrm{R}}}$
is equally small, which means that $\Sigma$, $\Upsilon$ and $\Omega$ (for the
doublet scenario) as well as $\sigma$ (for the singlet scenario), need to be
very small in order for both models to be consistent with the active neutrinos
mass scale. While a detailed discussion for the neutrino mass generation
mechanism is beyond the scope of this article, it is not a difficult task to
minimally extend both the doublet and singlet scenarios, e.g. by introducing
the dimension-five Weinberg operator, without affecting our analysis and
conclusions Strumia and Vissani (2006). However, the fact that $\Sigma$,
$\Upsilon$, $\Omega$ and $\sigma$ need to be very small is calling for the
introduction of an approximate global lepton number symmetry, which, for
convenience, we define as
$\mathrm{U}(1)_{e}\times\mathrm{U}(1)_{\mu}\times\mathrm{U}(1)_{\tau}$, such
that only leptons transform non-trivially according to the following quantum
numbers:

Field | $\mathrm{U(1)_{e}}$ | $\mathrm{U(1)_{\mu}}$ | $\mathrm{U(1)_{\tau}}$
---|---|---|---
$E_{\mathrm{L,R}}$ | $0$ | $1$ | $0$
$\mathcal{E}_{\mathrm{L,R}}$ | $0$ | $1$ | $0$
$\nu_{\mathrm{R}}$ | $0$ | $0$ | $0$
$L^{1}$ | $1$ | $0$ | $0$
$L^{2}$ | $0$ | $1$ | $0$
$L^{3}$ | $0$ | $0$ | $1$
$e_{\mathrm{R}}^{1}$ | $1$ | $0$ | $0$
$e_{\mathrm{R}}^{2}$ | $0$ | $1$ | $0$
$e_{\mathrm{R}}^{3}$ | $0$ | $0$ | $1$
Table 3:  Quantum numbers under the approximate lepton number symmetry.

Such an approximate symmetry constrains the model parameters as follows:

  * •

$\Sigma$, $\Upsilon$, $\Omega$ (for the doublet case) and $\sigma$ (for the
singlet case) must be tiny in consistency with the smallness of the active
neutrino masses;

  * •

the charged lepton-sector Yukawa couplings $\Pi$ (for the doublet model) and
$\pi$ (for the singlet model), are approximately diagonal and off-diagonal
elements can be ignored in our analysis;

  * •

The VLLs solely couple to muons such that only $\Theta_{2}$, $M_{{\rm LE},2}$
(for the doublet scenario) and $\theta_{2}$ (singlet scenario) are sizeable
and relevant for our numerical studies.

In what follows, we are essentially interested in lepton couplings to both
gauge and Higgs bosons. In particular, the focus is on the $W\nu
E(\mathcal{E})$, $Z^{0}E(\mathcal{E})\ell$ and $HE(\mathcal{E})\ell$ vertices.
All numerical computations are performed in the mass basis such that one needs
to rotate the gauge eigenstates via the bi-unitary transformations

$\displaystyle\mathrm{Leptons:}\hphantom{..}M_{L}^{\mathrm{diag}}=\quantity(U_{L}^{e})^{\dagger}M_{L}\quantity(U_{R}^{e}),$
(6)
$\displaystyle\mathrm{Neutrinos:}\hphantom{..}M_{\nu}^{\mathrm{diag}}=\quantity(U_{\nu})^{\dagger}M_{\nu}\quantity(U_{\nu}).$

The list of all Feynman rules that are relevant to the numerical analysis are
shown in Appendix A.

Before presenting the methodology behind the numerical calculations, it is
instructive to revisit the main experimental constraints for VLLs as well as
the most relevant couplings for this discussion. First, let us recall that
there are not many direct searches for VLLs, at least, in comparison with
their quark counterparts. In fact, LEP searches have lead to a lower bound on
the VLLs mass of $M_{\mathrm{VLL}}>100.6$ GeV Achard  _et al._ (2001), while
the most restrictive constraints are coming from CMS where direct searches for
doublet VLLs with sizeable couplings to taus are excluded at 95% Confidence
Level (CL) in the 120-790 GeV mass range Sirunyan  _et al._ (2019).
Additionally, there are searches for charged lepton resonances which constrain
the mass of a hypothetical VLL to be above a few hundred GeV, depending on the
underlying assumptions (114–176 GeV for VLL and 100–468 GeV for a heavy lepton
in the type-III inverse seesaw model Aad  _et al._ (2015), as well as 80-210
GeV in Ref. Chatrchyan  _et al._ (2012)). However, for the model under
consideration, the $[\mathrm{U}(1)]^{3}$ lepton number symmetry dictates that
the relevant couplings are those to muons such that the most restrictive
search in Sirunyan  _et al._ (2019) does not apply.

The left-mixing matrix that rotates from the flavour to the mass basis can be
parametrized by a single mixing angle, $\alpha$, as follows

$U_{L}^{e}=\begin{bmatrix}1&0&0&0\\\ 0&\cos{\alpha}&0&-\sin{\alpha}\\\
0&0&1&0\\\ 0&-\sin{\alpha}&0&-\cos{\alpha}\end{bmatrix}.$ (7)

With the
$\mathrm{U}(1)_{e}\times\mathrm{U}(1)_{\mu}\times\mathrm{U}(1)_{\tau}$
approximate symmetry the mass matrices for both the doublet and singlet cases
take a simple form

$M_{L}^{\mathrm{doublet}}=\begin{bmatrix}\dfrac{v\Pi_{11}}{\sqrt{2}}&0&0&0\\\
0&\dfrac{v\Pi_{22}}{\sqrt{2}}&0&M_{LE,2}\\\
0&0&\dfrac{v\Pi_{33}}{\sqrt{2}}&0\\\
0&\dfrac{v\Theta_{2}}{\sqrt{2}}&0&M_{E}\end{bmatrix}\,,\quad\quad
M_{L}^{\mathrm{singlet}}=\begin{bmatrix}\dfrac{v\pi_{11}}{\sqrt{2}}&0&0&0\\\
0&\dfrac{v\pi_{22}}{\sqrt{2}}&0&\dfrac{v\theta_{2}}{\sqrt{2}}\\\
0&0&\dfrac{v\pi_{33}}{\sqrt{2}}&0\\\ 0&0&0&M_{E}\end{bmatrix}\,.$ (8)

Besides lepton mixing terms one also considers a neutrino mixing with the
following structure,

$\displaystyle
U_{\nu}^{\mathrm{doublet}}=\begin{bmatrix}(U_{\mathrm{PMNS}})_{11}&(U_{\mathrm{PMNS}})_{12}&(U_{\mathrm{PMNS}})_{13}&0&0&0\\\
(U_{\mathrm{PMNS}})_{21}&(U_{\mathrm{PMNS}})_{22}&(U_{\mathrm{PMNS}})_{23}&0&0&0\\\
(U_{\mathrm{PMNS}})_{31}&(U_{\mathrm{PMNS}})_{32}&(U_{\mathrm{PMNS}})_{33}&0&0&0\\\
0&0&0&1&0&0\\\ 0&0&0&0&(U_{\nu})_{55}&(U_{\nu})_{56}\\\
0&0&0&0&(U_{\nu})_{65}&(U_{\nu})_{66}\\\ \end{bmatrix},$ (9) $\displaystyle
U_{\nu}^{\mathrm{singlet}}=\begin{bmatrix}(U_{\mathrm{PMNS}})_{11}&(U_{\mathrm{PMNS}})_{12}&(U_{\mathrm{PMNS}})_{13}&0\\\
(U_{\mathrm{PMNS}})_{21}&(U_{\mathrm{PMNS}})_{22}&(U_{\mathrm{PMNS}})_{23}&0\\\
(U_{\mathrm{PMNS}})_{31}&(U_{\mathrm{PMNS}})_{32}&(U_{\mathrm{PMNS}})_{33}&0\\\
0&0&0&1\end{bmatrix}$

where the $3\times 3$ SM block is fixed by the Pontecorvo–Maki–Nakagawa–Sakata
(PMNS) matrix Esteban  _et al._ (2020). Technically, and considering only such
a block, the PMNS matrix must in rigour be defined as
$U_{\mathrm{PMNS}}=U_{\nu}\cdot U_{L}^{\dagger}$ such that the actual mixing
is instead given as $(U_{\mathrm{PMNS}})_{i1}=(U_{\nu})_{i1}$,
$(U_{\mathrm{PMNS}})_{i3}=(U_{\nu})_{i3}$ and
$(U_{\mathrm{PMNS}})_{i2}=(U_{\nu})_{i2}\cos\alpha$ with $i=1,2,3$. However,
the benchmark scenarios considered in Sec. IV are such that $\cos\alpha\approx
1$ for both the doublet and the singlet models and therefore the
parametrization in Eq. 9 is consistent with a realistic lepton mixing.

The $2\times 2$ block in the doublet case represents the mixing between the
$\nu_{5}$ and $\nu_{6}$ mass eigenstates. These elements do not contribute to
the interaction vertices used in our analysis and therefore their size and
numerical values are not relevant for the discussion. For consistency, we
assume a generic $\mathcal{O}(1)$ mixing. It also follows from the approximate
lepton number symmetry that the mixing between the right-handed neutrino and
the remaining neutral leptons is negligible and can be ignored. This means
that $(U_{\nu})_{44}=1$, $(U_{\nu})_{4j}=(U_{\nu})_{j4}=0$ and the couplings
between the SM and BSM neutrinos can be neglected.

Let us note that if one requires consistency with the recently observed muon
$(g-2)_{\mu}$ anomaly Abi  _et al._ (2021), the preferred values of the cosine
of the mixing angle $\alpha$ must not be far from unity. On the other hand, a
small $\alpha$ (or equivalently, small Yukawa parameters $\theta_{2}$ and
$\Theta_{2}$) is preferable in order to keep tree-level modifications to the
Feynman rules for the $H\mu^{+}\mu-$ and $\mathrm{Z^{0}}\mu^{+}\mu^{-}$
vertices under control (a list of the rules are shown in appendix A). Indeed,
a large mixing angle would make the values of the branching fractions
$\mathrm{BR}(Z^{0}\rightarrow\mu^{+}\mu^{-})$ and
$\mathrm{BR}(H\rightarrow\mu^{+}\mu^{-})$ unacceptably large. Furthermore, and
noting that the VLL only couples to muons, one can see from del Aguila  _et
al._ (2008) that EW global fits constrain the mixing angle to be less than
0.034, which we impose in our numerical analysis. With such a small mixing,
lepton flavour violation observables can be kept under control and the VLL
single production cross-section is, for the singlet model, dominated by the
$W\nu_{\ell}\mathcal{E}$ vertex and, for the doublet model, by the
$W\nu_{L}^{\prime}E$ vertex, as we discuss in Sec. III.2. It follows from the
mixing structure in Eq. 9 (see also the Feynman rules in Eqs. 20, 21, 24 and
25) that interactions involving right-handed neutrinos are absent in both
scenarios. This is consistent with the approximate lepton number symmetry,
where one can assume that the smallness of the $\sigma$ and $\Sigma$ Yukawa
couplings, results in a negligible $W\nu_{R}\mathcal{E}$ and $W\nu_{R}E$
interaction strength. Therefore, in our numerical analysis, one can ignore
such interactions by safely setting it to zero.

Another important element to take into consideration is the VLL decay width.
In particular, it can be expressed as

$\displaystyle\Gamma(E(\mathcal{E})\rightarrow
W\nu_{\ell})=\frac{g^{2}U_{\mathrm{MIX}}^{2}m_{E,\mathcal{E}}^{2}}{64\pi
M_{W}^{2}}\quantity(1-\frac{M_{W}^{2}}{m_{E,\mathcal{E}}^{2}})^{2}\quantity(1+2\frac{M_{W}^{2}}{m_{E,\mathcal{E}}^{2}}).$
(10)

Here, $U_{\mathrm{MIX}}^{2}$ is a non-trivial combination of the neutrino and
charged lepton mixing matrices whose specific structure can be determined from
the Feynman rules in Eqs. 20, 21, 24 and 25. One notes that such a decay width
grows as a power law in $m_{E}$, becoming rather significant for large masses
as shown in Fig. 1. For example, taking the case where $U_{\mathrm{MIX}}=1$,
one can observe that for a mass of $2\leavevmode\nobreak\ \mathrm{TeV}$, the
decay width already surpasses the mass itself with $\Gamma=2705.4$ GeV. Such
an effect can be mitigated for smaller mixing angles as demonstrated in Fig.
1. Indeed, a decay width greater than its mass leads to conceptual issues with
respect to the interpretation of such large-width fields as particles. In
fact, within the context of Quantum Field Theory, a large value for the width
implies that the particle is highly non-local and the narrow-width
approximation is no longer valid calling for a more robust framework to
describe it (see e.g. Kuksa and Volchanskiy (2013); Kuksa (2006, 2009)). Such
a scenario is not relevant for the charged VLL as the mixing element
$U_{\mathrm{MIX}}$ must be small due to EW constraints involving Higgs and Z
boson decays. On the other hand, the VLL neutrino partners that reside in the
same doublet representation can also decay via a $W$ boson and a charged
lepton. In the limit of $m_{\nu_{L}^{\prime}}\gg m_{\ell}$, the expression for
the decay width is identical to that of the VLL (with the appropriate
changes),

$\displaystyle\Gamma(\nu_{L}^{\prime}\rightarrow
W\ell)=\frac{g^{2}U_{\mathrm{MIX}}^{2}m_{\nu_{L}^{\prime}}^{2}}{64\pi
M_{W}^{2}}\quantity(1-\frac{M_{W}^{2}}{m_{\nu_{L}^{\prime}}^{2}})^{2}\quantity(1+2\frac{M_{W}^{2}}{m_{\nu_{L}^{\prime}}^{2}}),$
(11)

where $U_{\mathrm{MIX}}$ can become sizable. Indeed, taking $\ell=\mu$, then,
in accordance with the Feynman rules, the mixing elements coming from the
left-chiral couplings are given as
$\cos{\alpha}\left[U_{\nu}^{\mathrm{doublet}}\right]_{55}$. Both
$\cos(\alpha)$ and $\left[U_{\nu}^{\mathrm{doublet}}\right]_{55}$ are of order
$\mathcal{O}(1)$ and therefore the $\nu_{L}^{\prime}\rightarrow W\ell$ decay
width becomes relevant for large masses. However, notice that the neutrino
mass is not independent of its charged VLL counterpart and thus cannot be made
arbitrarily small for the width not becoming too large. At leading order, both
the VLL and the heavy neutrino masses are controlled by the $M_{E}$ parameter
in Eq. 2 making them almost degenerate with a small deviation expected to be
induced at one-loop level Thomas and Wells (1998). Typically, such a mass
splitting is of the order of a few hundred MeV Thomas and Wells (1998), which
must not play a sizeable role. With this in mind, and without loss of
generality for this article’s analysis, we consider both the neutrino and its
corresponding charged VLL mass degenerate, ensuring that
$\Gamma(\nu_{L}^{\prime}\rightarrow W\ell)$ is never unacceptably large.
Although the widths must be kept under control, a sizeable value has
interesting experimental implications. For instance, it results in wider
distributions of kinematic variables such as the mass and the transverse
momentum ($p_{T}$), where large tails can extend into phase space regions not
populated by the SM background events.

Figure 1: Decay width as a function of the VLL’s mass, both displayed in GeV.
Each coloured line is representative of a different mixing element, whose
numerical value is displayed in the box to the right of the plot.

##  III VLLs: collider phenomenology and flavour constraints

###  III.1 Flavour constraints

It follows from interactions between the VLL and the muon that non-zero mixing
elements can trigger Lepton Flavour Violation (LFV) interactions, in
particular, those that relate to muon and tau decays. Therefore it is
important to confront all points considered in our numerical analysis with the
following LFV branching ratios (BR): BR($\mu\rightarrow e\gamma$),
BR($\tau\rightarrow e\gamma$), BR($\tau\rightarrow\mu\gamma$),
BR($\mu^{-}\rightarrow e^{-}e^{+}e^{-}$), BR($\tau^{-}\rightarrow
e^{-}e^{+}e^{-}$), BR($\tau^{-}\rightarrow\mu^{-}\mu^{+}\mu^{-}$),
BR($\tau^{-}\rightarrow e^{-}\mu^{+}\mu^{-}$),
BR($\tau^{-}\rightarrow\mu^{-}e^{+}e^{-}$),
BR($\tau^{-}\rightarrow\mu^{-}e^{+}\mu^{-}$), BR($\tau^{-}\rightarrow
e^{-}\mu^{+}e^{-}$), BR($\tau^{+}\rightarrow\pi^{0}e^{+}$) and
BR($\tau^{+}\rightarrow\pi^{0}\mu^{+}$). For such a purpose we use the latest
version of SPheno Porod and Staub (2012) to generate the Wilson coefficient
cards, which are then passed to flavio Straub (2018) in order to compute the
corresponding LFV observables. It is important to mention that SPheno also
computes the BRs i.e. flavio merely works as an extra layer of added scrutiny,
even if both programs’ numerical outputs (for the BRs) are well within each
other’s error and are therefore compatible. We use the 90% CL experimental
limits on the considered LFV observables as reported in the most recent issue
of the PDG review Zyla  _et al._ (2020).

Last but not least, let us comment that VLLs have long been proposed as a
potential explanation for the anomalous magnetic moment of the muon. However,
it follows from the smallness of the $\alpha$ mixing angle that neither the
singlet nor doublet models discussed in this article can successfully
accommodate such an anomaly as couplings of muons to VLLs and to Higgs or Z
bosons are too small. This was numerically verified and consistent with the
results obtained e.g. in Freitas  _et al._ (2014). Therefore, it is not
sufficient to extended the SM with one generation of VLLs to explain the muon
$g-2$ anomaly and additional new physics with less constrained couplings is
necessary. A possibility is to extend the scalar sector as proposed by various
authors Dermisek  _et al._ (2021); Kannike  _et al._ (2012); Jegerlehner and
Nyffeler (2009); Dermisek and Raval (2013); Crivellin  _et al._ (2018);
Bonilla  _et al._ (2021).

###  III.2 Collider analysis: LHC and muon colliders

The analysis techniques and Deep Learning algorithms used in the current study
for both the doublet and the singlet models were first implemented in Freitas
_et al._ (2021) and Bonilla  _et al._ (2021), respectively. Starting with the
doublet case, the channel we are probing is characterised by a final state
with 4 light jets, originating from the decays of two $W$ bosons, a charged
lepton, which we take to be the muon, and transverse missing energy (MET),
whose origin results from the undetected neutrinos in the final state. At
Leading Order (LO), the production diagram, at the LHC, can be seen in Fig. 2,
where the heavy $\nu_{L}^{\prime}$ neutrino belongs to a VLL doublet as shown
in Eq. 1. For this signal topology the main irreducible backgrounds are the
double top-quark production, $t\bar{t}$, $W+$jets (we consider up to 4 jets)
and Diboson plus jets. Notice that for the $t\bar{t}$ channel we consider that
one of the $W$’s decays in its hadronic channel with the other one decaying
into a muon and a neutrino. All backgrounds are generated at LO precision. For
the singlet scenario, the final state topology is characterized by one
isolated charged lepton from the $W$ boson decay, and transverse missing
energy associated with the undetected neutrinos. In particular, we focus on
the channel where the charged lepton is a muon. For this signal topology, we
consider the irreducible background together with all physics processes that
lead to the production of, at least, one muon and up to two jets in the final
state, $pp\rightarrow\mu^{-}\bar{\nu}_{\mu}$ (+jets)111In more practical
terms, this background is simulated in MadGraph via the commands generate p p
> mu+ vm && add process p p > mu+ vm j && add process p p > mu+ vm j j ..

For the case of the VLL singlet only one BSM coupling contributes to the decay
chain, in particular, the one containing a muon-specific neutrino, whereas for
the doublet scenario the couplings involving both the charged and neutral
components of the VLL doublet are present. Taking into account the Feynman
rules shown in appendix A, the strength of the relevant left and right chiral
charged current projections read as

$\displaystyle
g_{L}^{E\nu_{\mu}W}=\dfrac{g}{\sqrt{2}}\sin{\alpha}\left[U_{\nu}^{\mathrm{doublet}}\right]_{22},\qquad\qquad$
$\displaystyle
g_{L}^{\mathcal{E}\nu_{\mu}W}=\dfrac{g}{\sqrt{2}}\sin{\alpha}\left[U_{\nu}^{\mathrm{singlet}}\right]_{22},\qquad\qquad$
(12) $\displaystyle
g_{R}^{\nu_{L}^{\prime}\mu^{-}W^{+}}=\dfrac{g}{\sqrt{2}}\left[U_{\nu}^{\mathrm{doublet}}\right]_{55}\left[U_{R}^{e}\right]_{24},\qquad\qquad$
$\displaystyle
g_{L}^{\nu_{L}^{\prime}\mu^{-}W^{+}}=\dfrac{g}{\sqrt{2}}\sin{\alpha}\left[U_{\nu}^{\mathrm{doublet}}\right]_{55},$
$\displaystyle
g_{L}^{E\nu_{L}^{\prime}W}=\dfrac{g}{\sqrt{2}}\cos{\alpha}\left[U_{\nu}^{\mathrm{doublet}}\right]_{55},\qquad\qquad$
$\displaystyle
g_{R}^{E\nu_{L}^{\prime}W}=-\dfrac{g}{\sqrt{2}}\left[U^{e}_{R}\right]_{24}\left[U_{\nu}^{\mathrm{doublet}}\right]_{55}\,.$

This in turn means that the cross section for the doublet case will be much
larger than the one for the singlet case, because the strength of the dominant
contribution for the former, in particular the coupling between the VLL and
its neutrino counterpart, is proportional to $\cos{\alpha}$ while for the
latter all couplings are controlled by $\sin{\alpha}$. Notice that the single
production channel in Fig. 3 is valid for both the doublet and singlet
scenarios where, in the case of equal angles
$g_{L}^{E\nu_{\mu}W}=g_{L}^{\mathcal{E}\nu_{\mu}W}$. Therefore, it is
sufficient to study the singlet model and extract identical conclusions for
the doublet case.

Figure 2: Leading-order Feynman diagram for single production of the doublet
VLL. Here, $q$ and $\bar{q}$ correspond to quarks originating from the
colliding protons, $E$ represents the VLL and $\nu_{L}^{\prime}$ denotes the
VLL doublet partner. Both $W^{\prime}s$ in the hadronic channel with $j$
indicating first and second generation chiral quarks. Figure 3: Leading-order
Feynman diagram for VLL single production. $q$ and $\bar{q}$ corresponds to
quarks originating from the initial protons and $\nu_{\ell}$ are the SM
neutrinos.

The analysis follows a well-defined guideline. Both models are first
implemented in SARAH Staub (2014) to generate all interaction vertices as well
as all relevant files that interface with Monte-Carlo generators, namely,
MadGraph5 Alwall  _et al._ (2014) for quark-level matrix-element calculations
and Pythia8 Sjöstrand  _et al._ (2015) for hadronization and showering leading
to final state particles. In MadGraph5 we simulate proton-proton collisions at
a centre-of-mass energy $\sqrt{s}=14$ TeV, for a total of 250k events for the
signal and backgrounds. We also employ the default LO parton-distribution
function NNPDF2.3 Ball  _et al._ (2013) which fixes the evolution of the
strong coupling constant, $\alpha_{s}$. Fast-simulation of the ATLAS detector
is conducted, with Delphes de Favereau  _et al._ (2014). All kinematic
distributions are extracted with the help of ROOT Brun and Rademakers (1997).
At this stage, we also impose selection criteria, to maximise the signal
significance. In particular, all events must satisfy

$\displaystyle p_{T}(\mu^{-})>25\text{ GeV},$ $\displaystyle\rm{MET}>15\text{
GeV}\,,$ $\displaystyle\absolutevalue{\eta(\mu^{-})}\leq 2.5\,\rm{and}$
$\displaystyle\absolutevalue{\eta(j)}\leq 2.5\,.$

where $\eta$ represents the muon pseudo-rapidity222The pseudorapidity of a
particle is defined as $\eta=-\ln(\tan(\theta/2))$, where $\theta$ is its
polar angle.. Additionally, for the doublet topology, which involves jets in
the final state, we require them to be tagged as originating from light chiral
quarks, i.e., they cannot be tagged as b-jets. Note that the signal production
represented in Fig. 3, is characterized by having no jets in the final states.
While it is true that the lack of no jet candidates at any given event in a
hadron collider is perhaps too unrealistic, the absence of such selection can
also cause complications in the Deep Learning algorithms used later in the
analysis, since signal/background classes are incredibly unbalanced and may
lead to over-fitting problems.

Initiate NN models Top 5 selection: Act as F+M Select best model New NN models
“Daughters” “Mutated” NN models “Daughters” 200 epochsAfter 5 gen.F+M
traits$\mathcal{P}(M)=20\%$traindaughters

Figure 4: Flowchart representative of all iterative steps involved in the
genetic algorithm that we employ in this work.

The kinematics of final states for both signal and background events are
translated into tabular datasets and used as inputs for neural models, whose
job is to separate the signal from the background. The neural network is
constructed using Keras Chollet  _et al._ (2015). For classification of the
singlet model, we choose a set of 5 observables computed in the laboratory
frame. Low-level features include the muon observables
$\cos(\theta_{\mu^{-}})$, $\eta(\mu^{-})$, $\phi(\mu^{-})$ and
$p_{T}(\mu^{-})$ as well as the MET.

For the doublet model, a richer final state is present and, as such, one can
compute a much more complete list of observables. Do note that, for this
physics scenario, the final states contain 4 light-jets, which implies that
one cannot distinguish between the jets originating from $W^{+}$ and those
originating from $W^{-}$. Therefore, for reconstructed observables, we
consider all possible combinations of jets that can be used. These are $C$ =
$(j_{1},j_{2})$, $(j_{1},j_{3})$, $(j_{1},j_{4})$, $(j_{2},j_{3})$,
$(j_{2},j_{4})$ and $(j_{3},j_{4})$, where we define $j_{1}$ as the leading
jet, with the highest $p_{T}$, and $j_{4}$ as the subleading jet, with the
lowest $p_{T}$. The full list of observables that are used in the training are
shown in Tab. 4.

| Dimension-full | Dimensionless
---|---|---
Doublet VLL |  $p_{T}(\mu^{-})$, $p_{T}(j_{1})$, $p_{T}(j_{2})$, $p_{T}(j_{3})$ MET, $p_{T}(j_{4})$, $E(\mu^{-})$, $E(j_{1})$, $E(j_{2})$, $E(j_{3})$, $E(j_{4})$, $M(j_{1},j_{2})$, $M(j_{1},j_{3})$, $M(j_{a},j_{b})$, $p_{T}(j_{a},j_{b})$, $M(j_{a},j_{b},\mu^{-})$, $p_{T}(j_{a},j_{b},\mu^{-})$, $M(j_{a},j_{b},\mathrm{MET})$, $p_{T}(j_{a},j_{b},\mathrm{MET})$ |  $\cos(\theta_{j_{a}j_{b}})$, $\cos(\theta_{j_{a}\mu^{-}})$, $\eta(j_{1})$, $\eta(j_{2})$, $\eta(j_{3})$, $\eta(j_{4})$, $\eta(\mu^{-})$, $\phi(j_{1})$, $\phi(j_{2})$, $\phi(j_{3})$, $\eta(j_{4})$, $\eta(j_{a},j_{b})$, $\phi(j_{a},j_{b})$, $\eta(j_{a},j_{b},\mu^{-})$, $\eta(j_{a},j_{b},\mu^{-})$ |  $\Delta R(j_{a},j_{b})$, $\Delta R(j_{a},\mu^{-})$, $\Delta\phi(j_{a},j_{a})$, $\Delta\phi(j_{a},\mu^{-})$, $\Delta\theta(j_{a},j_{b})$, $\Delta\theta(j_{a},\mu^{-})$.

Table 4: Angular and kinematic distributions for the analysis of the doublet
production topology. All observables are computed in the laboratory frame of
reference. To simplify notation, we define $(j_{a},j_{b})$ as corresponding to
all jet combinations, as described in the text. Variables that involve
combinations of final states, e.g. $p_{T}(j_{a},j_{a},\mu^{-})$, are
reconstructed from the states indicated inside the parenthesis.

To optimize the neural network architecture, we employ a genetic algorithm
following the same steps as described in Freitas  _et al._ (2021); Bonilla
_et al._ (2021) and schematically represented in the diagram of Fig. 4. The
algorithm begins by first generating an arbitrary number of neural networks,
whose architecture is determined by randomly pooling a list of predefined
hyperparameters (number of layers, activation functions, regularisers, etc.).
We then train each individual network over the data for a given number of
epochs. From the trained networks, we pick the top five best performing ones.
From these, we create Father-Mother pairs, where we combine 50% of the
Father’s traits (that is, its hyperparameters) and 50% of the Mother traits to
construct new neural architectures, which we dub as Daughters. We also
consider the possibility of mutation, that is, after the Daughter networks
have been built, we consider that the hyperparameters may change to another,
with the probability of $\mathcal{P}(M)=20\%$. We then train the Daughter
networks for a given number of epochs and repeat the procedure for a set of
generations. Finally, we select the best performing network based on some
metric of choice. In our analysis, the algorithm is designed to maximise the
Asimov statistical significance, based on a earlier work of Elwood and Adam in
Elwood and Krücker (2018). The Asimov metric is defined as Cowan  _et al._
(2011)

$\mathcal{Z_{A}}=\Bigg{[}2\Bigg{(}(s+b)\ln\Bigg(\frac{(s+b)(b+\sigma_{b}^{2})}{b^{2}+(s+b)\sigma_{b}^{2}}\Bigg{missing})-\frac{b^{2}}{\sigma_{b}^{2}}\ln\Bigg(1+\frac{\sigma_{b}^{2}s}{b(b+\sigma_{b})}\Bigg{missing})\Bigg{)}\Bigg{]}^{1/2},$
(13)

where $s$ is the number of signal events, $b$ is the number of background
events and $\sigma_{b}$ is the uncertainty of the background. As part of the
optimization procedure, we consider the same list of hyperparameters as in our
previous work Freitas  _et al._ (2021):

  * •

number of hidden layers: 1 to 5

  * •

number of neurons in each layer: 256, 512, 1024 or 2048

  * •

kernel initializer: ’normal’,’he normal’,’he uniform’

  * •

L2 regularization penalty: 1e-3, 1e-5, 1e-7

  * •

activation function: ’relu’, ’elu’, ’tanh’, ’sigmoid’

  * •

optimizer: ’adam’, ’sgd’, ’adamax’, ’nadam’

Our evolutionary algorithm is initialized by building a set of ten NNs. The
parameters are chosen randomly from the previous lists. Each network is
trained up 200 epochs and if any improvement is not observed by at least 5
epochs, then the training stops. Note that a subset of the neural network
properties are not subject to optimisation but instead remain fixed during the
runs. Namely, we consider the following

  * •

The input data of the NNs are standardised, that is, input vectors have a mean
of zero and a standard deviation of 1. All observables were extracted from the
ROOT files and outputed into dataframes. The data is reshuffled and then
divided into training samples (80% of the total data) and validation samples
(the remainder 20%). We also consider cross-validation with a five-fold scheme
during training. The statistical significance is computed based on the trained
NN predictions of the validation data.

  * •

We employ a cyclic learning rate during the training phase with 0.01 initial
value and maximal value of 0.1.

  * •

The output of the NN is a vector of probabilities such that we define a signal
for an output with probability greater than 0.5, otherwise it is considered as
a background.

  * •

Batch size of 32768.

  * •

The best NN is selected based on the Asimov metric. The loss function for this
case is defined as the inverse of Eq. (13), such that when the loss function
is minimised, the Asimov significance is maximised. We also consider a fixed
value of the background uncertainty to $\sigma_{b}=10^{-1}$.

An additional consideration is the fact that our datasets are unbalanced,
which is in part a result of the selection criteria imposed on the final
states, reducing the allowed phase space for both signal and backgrounds.
Unbalanced datasets can lead to overfitted networks and, as such, must be
properly dealt with. In this work, we have utilised the Synthetic Minority
Oversampling Technique (SMOTE) Chawla  _et al._ (2011), which oversamples the
minority classes of our training data. Do note that this algorithm is only
employed to the training dataset and we do not perform any re-sampling of the
validation samples.

With the LHC being a proton collider, production of coloured particles by the
strong interaction is heavily favoured, compared to identical processes
involving electroweak bosons. As VLLs are colour singlets, they can only be
produced via this last interaction at LO. It is then worth mentioning that, in
addition to the obvious motivation of studying VLL production at the LHC
already during RUN3 and, later on, in its high luminosity phase (HL-LHC), it
is of utmost importance to understand the sensitivity with which these new
particles can be probed, at future colliders.

Figure 5: LO Feynman diagrams for VLL double production at a lepton/anti-
lepton collider via exchange of either a photon or a $\mathrm{Z^{0}}$ boson.

In particular, there has been an active discussion within the community about
$e^{+}e^{-}$ colliders like the Compact Linear Collider (CLIC) Charles  _et
al._ (2018); Boland  _et al._ (2016) or the International Linear Collider
(ILC) Bae (2013), and more recently, on the possibility of building a
$\mu^{-}\mu^{+}$ collider Shiltsev (2021); Muo (2021). Besides offering
cleaner environments, when compared to hadronic machines, production via
electroweak processes is favoured, hence, VLLs should have a higher chance of
being observed, in case they exist, in this type of machines. For this
purpose, we perform numerical and analytical computations with FeynCalc
Shtabovenko  _et al._ (2020) for the pair-production of VLLs whose tree-level
diagrams can be seen in Fig. 5. As previously stated, singlets have lower
couplings compared to doublets, hence we discuss the prospects for VLL
discovery at lepton colliders for the singlet scenario and not the doublet, as
this can be seen as the worst-case scenario. Do note, however, that the
conclusions one takes from the singlet model can be easily generalized to the
doublet model. That is, the only difference between the collider studies for
the two scenarios are the values of the couplings.

In conclusion, the single VLL production is favoured at the LHC for the
doublet case, due mainly to the strength of the couplings. Double production
of VLLs has the obvious disadvantage of the need to produce two heavy states
and a more elaborate final state to detect. Additionally, as it was shown in
previous works Freitas  _et al._ (2021); Bonilla  _et al._ (2021), double
production topologies are sub-leading when compared with single production
ones and, as such, are not considered in this work. If one wants to probe the
singlet model, and assuming that it cannot be done at the LHC, the future
lepton colliders, where single production is precluded at LO, may give us
enough VLLs in the double production channel to test the singlet VLL scenario.

##  IV Results

In this section, we present a numerical analysis focusing on the LHC and
future lepton colliders. In Eqs. 14 and 15 we show two distinct viable
scenarios (one for the doublet and one for the singlet model), where all
couplings are fixed and only the VLL mass varies. The colour coding, black for
the doublet model and green for the singlet case, merely serves to distinguish
each of the benchmark scenarios and identify them throughout the text whenever
necessary. We first consider a benchmark point where we set the VLL mass to
$M_{\rm VLL}=700$ GeV. The remaining BSM neutrinos present in the doublet
model share, at LO, the same mass with their $\mathrm{SU}(2)_{\mathrm{L}}$
doublet VLL counterpart and, as noted in LABEL:eq:EnuW, these neutrinos
efficiently decay into muons and $W^{\pm}$ bosons, with an interaction
strength of the order of the weak gauge coupling, where, for the considered
benchmarks,
$\left[U^{\mathrm{doublet}}_{\nu}\right]_{55}\approx\left[U_{R}^{e}\right]_{24}\approx
1$. Indeed, the right chiral component of the coupling is the dominating
factor, as the left-coupling is suppressed by a factor of $\sin(\alpha)$,
which is of $\mathcal{O}(10^{-2})$ in our numerical analysis. The selected
benchmark point, represented in black, for the doublet model, compatible with
all flavour observables discussed above reads as

$\displaystyle U^{e}_{L}=\begin{bmatrix}1&0&0&0\\\
0&0.9997545&0&-0.02215850\\\ 0&0&1&0\\\
0&-0.02215850&0&0.9997545\end{bmatrix},\quad
U^{e}_{R}=\begin{bmatrix}1&0&0&0\\\ 0&0.999842&0&0.0177691\\\ 0&0&-1&0\\\
0&0.0177691&0&0.999842\end{bmatrix},$ (14) $\displaystyle
U_{\nu}=\begin{bmatrix}U_{\mathrm{PMNS}}&0\\\
0&U_{\nu}^{\mathrm{BSM}}\end{bmatrix},\quad
U_{\nu}^{\mathrm{BSM}}=\begin{bmatrix}1&0&0\\\ 0&-0.960772&-0.277338\\\
0&0.277338&-0.960772\end{bmatrix},$

whereas for the singlet scenario, in green, one has

$\displaystyle{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}U^{e}_{L}=\begin{bmatrix}1&0&0&0\\\
0&0.9995324&0&-0.03057725\\\ 0&0&1&0\\\
0&-0.03057725&0&0.9995324\end{bmatrix}},\quad{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}U^{e}_{R}=\begin{bmatrix}1&0&0&0\\\
0&1&0&-1.57282\times 10^{-5}\\\ 0&0&-1&0\\\ 0&-1.57282\times
10^{-5}&0&-1\end{bmatrix}},$ (15) $\displaystyle
U_{\nu}=\begin{bmatrix}U_{\mathrm{PMNS}}&0\\\ 0&1\end{bmatrix}.$

The neutrino mixing is identical in the two considered cases and therefore is
always represented in black.

Figure 6: Production cross section for both doublet and singlet topologies, in
femtobarn (fb), as a function of the VLL’s mass, in GeV. The $y$-axis is shown
in logarithmic scale. The black and green curves correspond to the same colour
mixing scenarios in Eqs. 14 and 15.

It follows from LABEL:eq:EnuW and the numerical values shown above that, upon
applying the selection criteria described in Section III.2, the signal and
background cross sections are

$\displaystyle\sigma(\text{doublet}):0.19\hphantom{.}\mathrm{fb},$
$\displaystyle\sigma(\text{singlet}):5.94\times
10^{-6}\hphantom{.}\mathrm{fb},$
$\displaystyle\sigma(pp\rightarrow\mu^{-}\bar{\nu}_{\mu}\text{+ up to 2 jets):
}4.79\times 10^{6}\hphantom{.}\mathrm{fb},$
$\displaystyle\sigma(t\bar{t}):13616.16\hphantom{.}\mathrm{fb},$
$\displaystyle\sigma(W+\mathrm{jets}):102448.54\hphantom{.}\mathrm{fb},$
$\displaystyle\sigma(\mathrm{Diboson}+{\mathrm{jets}}):296.94\hphantom{.}\mathrm{fb},$

It is evident that the doublet model cross section is far larger than that of
the singlet. This is indeed expected as the dominant contribution in the
singlet model contains a $\sin{\alpha}$ suppression factor as one can see in
LABEL:eq:EnuW. The dependency of the production cross section with the VLL
mass is shown in Fig. 6. One can now estimate the total expected number of
events, given the cross section above, after event selection. Assuming the
target luminosity of the HL-LHC, $\mathcal{L}=3000$ $\mathrm{fb^{-1}}$ and
that the expected number of events is given by $N=\sigma\mathcal{L}$, we have
for the doublet case $N=570.0$ events, while for the singlet we obtain
$N=0.01782$ events. For the singlet model, this implies that the production
cross section is not large enough to generate one event at the LHC. For this
reason it is meaningless to present the statistical significance for the
singlet model for such a heavy VLL. However, we see from Fig. 6 that, for
lighter singlet VLLs, in particular for their masses of $100$ and
$200\leavevmode\nobreak\ \mathrm{GeV}$, one has $N\approx 30$ and $N\approx
3$, respectively. Thus, it becomes possible to produce them at the HL-LHC
motivating a further full analysis.

The kinematic features used in the Deep Learning analysis can be seen in Figs.
10 and 11 of Appendix B for the doublet and singlet models, respectively. For
the singlet model, the main variables allowing for a good discrimination
between the signal and the background are the transverse momentum of the muon
i.e., $p_{T}(\mu^{-})$ as well as the MET distribution. These distributions
are characterized by long tails at higher energies where the SM background is
no longer present. On the other hand, the angular distributions for the cosine
of the polar angle, as well as the azimuthal angle, offer the least
discriminating power. This follows from the fact that both signal and
background have a similar shape. The pseudo-rapidity distributions can also be
used in the discrimination since signal events tend to peak around $\eta\sim
0$ whereas the SM backgrounds spread out over the entire
$\absolutevalue{\eta}\leq 2.5$ region. In particular, for the $\eta(\mu^{-})$
distributions, the SM backgrounds spread uniformly in the entire range. For
the doublet model, we note that the kinematic distributions offer the best
discriminating power, with both the $p_{T}$ and mass distributions peaking at
higher energies than those of the main irreducible backgrounds. While
$\Delta\theta$ distributions closely follow the SM background prediction, the
$\Delta\phi$ and $\Delta R$ distributions are distinct from those of the SM,
with $\Delta\phi$ possessing a double peak structure near zero, whereas
$\Delta R$ have its maximum at zero.

With this information at hand, one can construct multi-dimensional
distributions to be used as inputs for a neural network that solves a
classification task. This is done via an evolution algorithm to optimise the
various hyper-parameters of the neural model, whose metric to be maximised is
the Asimov statistical significance. For completeness, we present our results
for different statistical models, some more conservative than others. In this
article we use the same measures as in Freitas  _et al._ (2021); Bonilla  _et
al._ (2021), which include:

  1. 1.

The Asimov significance, $\mathcal{Z}_{A}$, with 1% systematic uncertainty 333
The dominant backgrounds of the analysis ($t\bar{t}$ and W+jets) have
currently their cross sections measured with $\mathcal{O}$(2%) precision with
an important contribution from the integrated luminosity measurement which
strongly dominates the total uncertainty. When preparing for the European
Strategy for Particle Physics reports that both ATLAS and CMS contributed to,
a range of 1-1.5% was considered as a guideline for the luminosity uncertainty
of the HL-LHC CER (2018). Furthermore, in a recent update ATL (2022) this
recommendation was further improved to $0.83\%$, based on the full run-II
accumulated dataset. As such, a global of 1% in our backgrounds can be seen as
a realistic target for both run-III and the HL phase of the LHC., as defined
in Eq. (13);

  2. 2.

A less conservative version of the Asimov significance, which we dub as
$\mathcal{Z}(<1\%)$. In this measure, we assume that backgrounds are known
with an error of $10^{-3}$. Of all measures, this is the most lenient one and
typically offers the most significant results;

  3. 3.

The more traditional metric, $s/\sqrt{s+b}$.

With this in mind, we compute the significances, for the doublet scenario, in
a wide range of masses, from 100 to 1000 GeV, in steps of 100 GeV. In
particular, we plot the various metrics as a function of the neural network
score in Fig. 7, taking a VLL mass equal to 700 GeV, for illustration
purposes. For this mass point and, for an integrated luminosity of
$\mathcal{L}=3000$ $\mathrm{fb^{-1}}$, we obtain

$\displaystyle s/\sqrt{s+b}=0.61\sigma,$
$\displaystyle\mathcal{Z}(<1\%)=11.62\sigma,$
$\displaystyle\mathcal{Z}_{A}=0.063\sigma,$

(a)  (b)  (c)

Figure 7: Statistical significance as a function of the classifier score given
by a neural network for different metrics, assuming a doublet-type VLL with
$M_{\mathrm{VLL}}=700$ GeV and a collider luminosity $\mathcal{L}=3000$
$\mathrm{fb}^{-1}$. The computation of the statistical significance is made
using the best neural network that the evolution algorithm found. From left to
right we plot the significance $s/\sqrt{s+b}$ in (a), the adapted Asimov
significance where we assume a background uncertainty of $10^{-3}$ in (b), and
the Asimov significance with systematics of 1% in (c).

where we note that we can exclude (or claim an hypothetical discovery) for a
VLL with such a mass, since the $\mathcal{Z}(<1\%)$ metric offers a
statistical significance larger than $5\sigma$. However, one must keep in mind
that this metric is the least conservative of the those considered in this
work. On the other hand, the Asimov metric is the stringiest and most
conservative one resulting in a tiny significance for this particular point.
It is also important to study the role of different values of the luminosity.
As such, in Fig. 8, we show the evolution of the significance as a function of
the collider’s luminosity for a VLL with a mass of $700\leavevmode\nobreak\
\mathrm{GeV}$. In particular, we highlight with dashed black vertical lines
the target luminosities at Run III (300 $\mathrm{fb^{-1}}$) and at the HL-LHC
(3000 $\mathrm{fb^{-1}}$). Focusing on Run III, we obtain

$\displaystyle s/\sqrt{s+b}=0.19\sigma,$
$\displaystyle\mathcal{Z}(<1\%)=2.95\sigma,$
$\displaystyle\mathcal{Z}_{A}=0.059\sigma,$

such that, for a $700\leavevmode\nobreak\ \mathrm{GeV}$ VLL, one does not
expect any significant excess and therefore can not extract conclusions.

Figure 8: Statistical significance as a function of the collider luminosity,
for a fixed doublet VLL mass of $700\leavevmode\nobreak\ \mathrm{GeV}$ (mixing
in Eq. 14). Target luminosities of the HL-LHC (3000 $\mathrm{fb}^{-1}$) and
Run III (300 $\mathrm{fb}^{-1}$) are marked with vertical dashed black lines.
Both axes are shown in logarithmic scale.

For an inclusive picture we perform a mass scan in the range
$M_{\mathrm{VLL}}\in[100,1000]$ GeV, whose numerical values for the couplings
are fixed to those shown in Eq. 14. On the other hand, it follows from the
discussion in the previous paragraphs that only two example points are
selected for the singlet model, both featuring the green mixing matrices in
Eq. 15. We also vary the mass of the heavy left-handed neutrinos, present in
the doublet model, such that they share the same mass with their charged
$\mathrm{SU}(2)_{\mathrm{L}}$ doublet counterpart. In Tab. 5 our results for
the scan are summarized, where the VLL masses and the calculated significance
for both Run III and the HL-LHC upgrade are shown.

Mass of VLL (GeV) |  $300$ $\mathrm{fb^{-1}}$ |  $3000$ $\mathrm{fb^{-1}}$
---|---|---
$s/\sqrt{s+b}$ | $\mathcal{Z}(<1\%)$ | $\mathcal{Z}_{A}$ | $s/\sqrt{s+b}$ | $\mathcal{Z}(<1\%)$ | $\mathcal{Z}_{A}$
$100$ | $75.32$ | $81.76$ | $25.01$ | $238.17$ | $258.55$ | $26.13$
$200$ | $28.97$ | $29.92$ | $9.37$ | $91.62$ | $94.61$ | $9.81$
$300$ | $8.02$ | $8.09$ | $2.56$ | $25.35$ | $25.58$ | $2.68$
$400$ | $2.76$ | $2.77$ | $0.88$ | $8.73$ | $8.76$ | $0.92$
$500$ | $1.05$ | $2.70$ | $0.33$ | $3.32$ | $11.19$ | $0.35$
$600$ | $0.43$ | $2.82$ | $0.14$ | $1.37$ | $11.43$ | $0.14$
$700$ | $0.19$ | $2.95$ | $0.060$ | $0.61$ | $11.62$ | $0.063$
$800$ | $0.086$ | $1.58$ | $0.027$ | $0.27$ | $7.94$ | $0.029$
$900$ | $0.042$ | $2.98$ | $0.013$ | $0.13$ | $11.57$ | $0.014$
$1000$ | $0.021$ | $1.84$ | $0.0066$ | $0.066$ | $11.14$ | $0.0069$
$\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}{100}$ | $0.94$ | $1.33$ | $0.00094$ | $2.98$ | $4.22$ | $0.0023$
$\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}{200}$ | $0.25$ | $0.36$ | $0.00078$ | $0.80$ | $1.13$ | $0.0023$

Table 5: Signal significance for VLL single-production at the LHC calculated
with an evolution algorithm that maximises the Asimov metric. The last two
rows represent the singlet model whereas the remaining ten are benchmark
points of the doublet model. The colour coding is the same as in Eqs. 14 and
15.

Mass of VLL (GeV) |  $3000$ $\mathrm{fb^{-1}}$
---|---
$s/\sqrt{s+b}$ | $\mathcal{Z}(<1\%)$ | $\mathcal{Z}_{A}$
$100$ | $238.17$ | $258.55$ | $26.13$
$200$ | $91.62$ | $94.61$ | $9.81$
$300$ | $25.35$ | $25.58$ | $2.68$
$400$ | $8.73$ | $10.19$ | $0.92$
$500$ | $3.32$ | $11.14$ | $0.35$
$600$ | $1.37$ | $11.45$ | $0.14$
${\color[rgb]{1,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{1,0,0}700}$ | ${\color[rgb]{1,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{1,0,0}0.61}$ | ${\color[rgb]{1,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{1,0,0}11.62}$ | ${\color[rgb]{1,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{1,0,0}0.063}$
$800$ | $0.38$ | $11.58$ | $0.029$
$900$ | $0.21$ | $11.54$ | $0.014$
$1000$ | $0.26$ | $11.14$ | $0.0069$

Table 6: Signal significance for VLL single-production at the LHC calculated
with the neural network trained in the 700 GeV data, with the network details
shown in Tab. 9

We notice that, both for $\mathcal{L}=3000$ $\mathrm{fb^{-1}}$ or
$\mathcal{L}=300$ $\mathrm{fb^{-1}}$ as well as any of the mass values, we are
able to obtain significances greater than $5\sigma$ if the simplified Asimov
metric and $s/\sqrt{s+b}$ are considered. These results may indeed suggest
that, significances above 5$\sigma$ can be achievable even for larger VLL mass
values at the LHC. This is particularly relevant in the doublet model while
for singlet VLLs the significance quickly drops if we go beyond
$200\leavevmode\nobreak\ \mathrm{GeV}$. If only the Asimov significance,
$\mathcal{Z}_{A}$, which is the most conservative one, is considered we can
still probe doublet VLLs up to about $200\leavevmode\nobreak\ \mathrm{GeV}$
already at the LHC Run III. For the singlet scenario, considering the HL-LHC
program, we can not exclude/discover VLLs, with the hightest significance
being 4.22$\sigma$ for the $\mathcal{Z}(<1\%)$ metric. While bellow the
discovery threshold, it can still be regarded as a potential anomaly, whose
significance can be increased with the combination of additional channels.

Although LEP constraints have already excluded VLLs up to
$100.6\leavevmode\nobreak\ \mathrm{GeV}$, our results for
$100\leavevmode\nobreak\ \mathrm{GeV}$ are merely indicative of how large can
the significance become for small doublet masses, if our analysis technique is
employed.

As one can note, there is an overall trend of the significance dropping as the
mass of the VLL decreases. However, we note that for masses above 600 GeV, we
have consistently obtained signifcances close to $10\sigma$ for the
$\mathcal{Z}(<1\%)$ measure. Here, there are two main factors at play. First,
considering larger masses, kinematic distributions tend to peak at higher
energies, which is particularly relevant for mass distributions, enhancing the
neural networks discriminating power. Additionally, we employ the evolutionary
algorithm to every single point, meaning that each network is optimized to a
specific phase space region444The neural networks found to be optimized for
each point are shown in appendix C.. In a realistic search scenario we would
not know the mass of the VLL and therefore it would not be reasonable to set a
stronger preference in one of the various networks optimized towards different
VLL masses. However, the networks that we use are engineered to be generic
enough and can in principle be applied to distinct masses up to a certain
discrimination power. To illustrate this, we apply the 700 GeV network to all
scanned masses and present the results in Tab. 6. As one can see, the
numerical values of the significance changed for the more lenient metrics,
$s\sqrt{s+b}$ and $\mathcal{Z}(<1\%)$, whereas for the most conservative one
the significance remained the same. However, an increase was also experienced,
in particular, for a VLL mass of 800 GeV, where the lenient Asimov metric grew
from $7.84\sigma$ to $11.58\sigma$. These results indicate that, despite the
networks’ training in distinct phase space regions, they are versatile enough
to offer a good discrimination power, which also indicates absence of over-
fitting. In essence, the hypothetical observation of a statistical excess for
a given mass could motivate employing an optimized network to potentially
enhance such a signal/excess.

Figure 9: VLL pair production cross section, in fb, as a function of centre-
of-mass energy for a lepton collider, in GeV. Each individual box is
representative of a different mass of the singlet VLL, which is also given in
GeV.

The low cross-sections of the singlet model call for a different approach on
how to probe them. It is in this context that the near future electron and
muon colliders can offer new opportunities. As it was mentioned in Sec. III.2,
production of particles via electroweak processes is favoured in these
colliders. As such, it is instructive to understand how can a lepton collider
enhance the cross-section for the case of pair-produced VLLs via an
$s$-channel process. Note that the analysis that follows is independent of the
chosen collider, as the $s$-channel cross-section is independent of the mass
of the initial colliding particles. Therefore, all results shown here are
valid both for the electron and muon machines. At LO, the main diagrams
involved are shown in Fig. 5. Note that, since we are assuming non-zero
couplings between the muon/electron and the VLL, there are also $t$-channel
contributions of the form

$\raisebox{-52.00005pt}{\includegraphics{TChannel-
Photon.pdf}}\quad+\quad\raisebox{-52.00005pt}{\includegraphics{TChannel-
ZBoson.pdf}}.$ (16)

However, such contributions are sub-leading when compared to production via
the $s$-channel process, since these depend on mixing structures of the VLL
with the muon, whereas in the $s$-channel the interaction vertices feature two
leptons of the same flavour coupling directly to vectors via gauge
interactions. Furthermore, in the particular case of our models, the
approximate family symmetry that only allows couplings between muons and VLLs
results in vanishing $t$-channel contributions at electron colliders, but not
at muon machines. Note that this channel can be seen as a direct probe to the
flavour structure in the leptonic sector if both muon and electron colliders
become operational. It is then safe to neglect the $t$-channel contribution in
the remainder of this analysis, where we use the same mixing structure as
defined in Eq. 15. We then consider the $\mathcal{E}\to W\nu_{\ell}$ channel,
which, depending on how the $W$ bosons decay, can lead to the $2\nu_{\ell}+4j$
or $2\nu_{\ell}+2j+\ell\nu_{\ell}$ or
$2\nu_{\ell}+\ell\nu_{\ell}+\ell\nu_{\ell}$ final states. For the topologies
that involve jets as final states, the main backgrounds include diboson
production ($WW$, $W\mathrm{Z^{0}}$ and $\mathrm{Z^{0}Z^{0}}$) with associated
jets as well as tau pair production. For the purely leptonic final state
$2\nu_{\ell}+\ell\nu_{\ell}+\ell\nu_{\ell}$, diboson production is relevant.
Indeed, a lepton collider is a much cleaner environment when compared to a
hadronic machines in such a way that backgrounds involving jets can be safely
discarded. VLLs can also decay into $\mathrm{Z^{0}}$ bosons as
$E(\mathcal{E})\rightarrow\mu\mathrm{Z^{0}}$. Such a channel would be an
important test for a hypothetical VLL discovery provided that the final states
can contain at least 6 charged leptons. Besides being a very clean process, 6
lepton topologies are also expected to be small in the context of a SM
background that typically results from triboson production processes. A
detailed phenomenological analysis at lepton colliders using analogous deep
learning methods to those discussed in the context of the LHC is beyond the
scope of this article and is left for future work.

For a comprehensive understanding on how the cross section depends on a lepton
collider center-of-mass energy, $E_{\mathrm{CM}}$, we use FeynCalc to obtain
the following expression

$\displaystyle\sigma=\dfrac{\sqrt{\dfrac{E_{\mathrm{CM}}^{2}}{4}-M_{\mathcal{E}}^{2}}}{256\pi
E_{\mathrm{CM}}^{7}(M_{\mathrm{Z^{0}}}^{2}-E_{\mathrm{CM}}^{2})^{2}}($
$\displaystyle\alpha_{1}E_{\mathrm{CM}}^{8}+\alpha_{2}E_{\mathrm{CM}}^{6}M_{\mathcal{E}}^{2}+\alpha_{3}E_{\mathrm{CM}}^{6}M_{\mathrm{Z^{0}}}^{2}+\alpha_{4}E_{\mathrm{CM}}^{4}M_{\mathrm{Z^{0}}}^{2}M_{\mathcal{E}}^{2}+\alpha_{5}E_{\mathrm{CM}}^{4}M_{\mathrm{Z^{0}}}^{4}+$
(17)
$\displaystyle+\alpha_{6}E_{\mathrm{CM}}^{2}M_{\mathrm{Z^{0}}}^{4}M_{\mathcal{E}}^{2}),$

where $\alpha_{i}$ for $i=1,\dots,6$ are dimensionless constants, proportional
to the product of various couplings, i.e.
$\alpha_{j}=\alpha_{j}(U_{\mathrm{L}}^{e},U_{\mathrm{R}}^{e},g,g^{\prime},\theta_{W})$,
$M_{\mathrm{Z^{0}}}$ is the mass of the $\mathrm{Z^{0}}$ boson and
$M_{\mathcal{E}}$ the mass of the singlet VLL. With the numerical values in
Eq. 15 these constants read as

$\displaystyle\alpha_{1}=0.189325,\quad\alpha_{2}=-0.213435,\quad\alpha_{3}=-0.328608,\quad\alpha_{4}=0.403093,\quad\alpha_{5}=0.152126,\quad\alpha_{6}=-0.209655.$
(18)

In Fig. 9 we plot the corresponding cross-section as a function of the centre-
of-mass energy, for VLL masses in the range between $500\leavevmode\nobreak\
\mathrm{GeV}$ and $6\leavevmode\nobreak\ \mathrm{TeV}$. Do note that we are
studying the singlet scenario, where the decay width always remain bellow the
mass of the VLL. It is interesting to note that the cross sections are rather
large, above 13.23 $\mathrm{fb}$ for a VLL mass of $500\leavevmode\nobreak\
\mathrm{GeV}$ and $E_{\mathrm{CM}}=1.5$ TeV, dropping to
$0.095\leavevmode\nobreak\ \mathrm{fb}$ for a VLL mass of $6\mathrm{TeV}$ and
$E_{\mathrm{CM}}=14$ TeV. This increase of the cross section allows to probe
higher mass ranges than the ones at the reach of the LHC, in a much cleaner
environment. Immediately noticeable is the fact that the cross section hits a
maximum shortly after $E_{\mathrm{CM}}\sim 2M_{\mathrm{VLL}}$, with a
subsequent drop as $E_{\mathrm{CM}}$ increases. This implies that when the
collider beam energy is two times the mass of the VLL, the double production
of VLL is enhanced and the discovery potential is maximized. For lower masses,
the drop is more pronounced when compared to the high-mass regime, essentially
because we are taking $E_{\mathrm{CM}}$ in the range of $3-14$ TeV for low
masses, while for high masses it occurs for much higher centre-of-mass
energies (beyond $E_{\mathrm{CM}}=14$ TeV), and therefore not as relevant for
the proposed lepton colliders. In particular, for $M_{\rm VLL}=500$ GeV and at
$E_{\mathrm{CM}}=14$ TeV we have $\sigma=0.23$ $\mathrm{fb}$, while for the
same mass, at $E_{\mathrm{CM}}=3$ TeV, we have $\sigma=4.63$ $\mathrm{fb}$.

| 200 GeV | 700 GeV | 1200 GeV | 1700 GeV | 2200 GeV | 2700 GeV | 3200 GeV | 3700 GeV
---|---|---|---|---|---|---|---|---
$E_{\mathrm{CM}}=1.5$ TeV | 19.18 $\mathrm{fb}$ | 5.51 $\mathrm{fb}$ | $-$ | $-$ | $-$ | $-$ | $-$ | $-$
$E_{\mathrm{CM}}=3$ TeV | 5.00 $\mathrm{fb}$ | 4.21 $\mathrm{fb}$ | 2.49 $\mathrm{fb}$ | $-$ | $-$ | $-$ | $-$ | $-$
$E_{\mathrm{CM}}=10$ TeV | 0.46 $\mathrm{fb}$ | 0.45 $\mathrm{fb}$ | 0.44 $\mathrm{fb}$ | 0.41 $\mathrm{fb}$ | 0.39 $\mathrm{fb}$ | 0.35 $\mathrm{fb}$ | 0.31 $\mathrm{fb}$ | 0.26 $\mathrm{fb}$
$E_{\mathrm{CM}}=14$ TeV | 0.23 $\mathrm{fb}$ | 0.23 $\mathrm{fb}$ | 0.23 $\mathrm{fb}$ | 0.22 $\mathrm{fb}$ | 0.21 $\mathrm{fb}$ | 0.21 $\mathrm{fb}$ | 0.19 $\mathrm{fb}$ | 0.18 $\mathrm{fb}$
$E_{\mathrm{CM}}=14$ TeV ($e^{+}e^{-}\gamma\gamma$) | 17.95 $\mathrm{fb}$ | 0.54 $\mathrm{fb}$ | 0.095 $\mathrm{fb}$ | 0.027 $\mathrm{fb}$ | 0.0095 $\mathrm{fb}$ | 0.0037 $\mathrm{fb}$ | 0.001585 $\mathrm{fb}$ | 0.00069 $\mathrm{fb}$
$E_{\mathrm{CM}}=14$ TeV ($\mu^{+}\mu^{-}\gamma\gamma$) | 9.10 $\mathrm{fb}$ | 0.26 $\mathrm{fb}$ | 0.045 $\mathrm{fb}$ | 0.012 $\mathrm{fb}$ | 0.0043 $\mathrm{fb}$ | 0.0017 $\mathrm{fb}$ | 0.00070 $\mathrm{fb}$ | 0.00030 $\mathrm{fb}$

Table 7: Singlet VLL double production cross-section at lepton colliders. As
idicated in the top row of the table, all considered mass values correspond to
the green benchmark point in Eq. 15. We show the cross-section for a centre-
of-mass energy $E_{\mathrm{CM}}=1.5$ and 3 TeV, which are the target energies
for the CLIC collider, whereas $E_{\mathrm{CM}}=3$, 10 and 14 TeV correspond
to the current proposal for the future $\mu^{+}\mu^{-}$ collider. Points
marked with “$-$” indicate that there is not enough energy to pair-produce the
particles at that mass.

Fixing the centre-of-mass energy and looking at various mass points, as shown
in Tab. 7, we notice that for the proposed high-energy colliders, the
variation of the mass does not cause significant deviations in the cross-
section. For example, taking $E_{\mathrm{CM}}=14$ TeV, we observe that the
cross-section remains nearly constant for the displayed masses, ranging from
$0.23\leavevmode\nobreak\ \mathrm{fb}$ for a 200 GeV VLL to
$0.18\leavevmode\nobreak\ \mathrm{fb}$ for a $3.7\leavevmode\nobreak\
\mathrm{TeV}$ one.

We end this section with a comment about the cross-sections at the high energy
end of the colliders. As the energy grows the $s$-channel cross sections
decrease. However, there is an alternative process that grows with
$\ln^{2}(s/m_{f}^{2})$, where $f$ stands for the incoming fermion. This is
$e^{+}e^{-}\to e^{+}e^{-}E({\mathcal{E}})\bar{E}(\bar{\mathcal{E}})$ for an
electron-positron collider and
$\mu^{+}\mu^{-}\to\mu^{+}\mu^{-}E({\mathcal{E}})\bar{E}(\bar{\mathcal{E}})$
for a muon collider. The photon fusion processes have cross-sections that grow
with the centre-of-mass energy Budnev  _et al._ (1975); Frixione  _et al._
(1993) and although they are not competitive for the lower energies they
become dominant at high energies. We present in Tab. 7 the values of these
cross-sections for an energy of 14 TeV which shows that they can play an
important role for very high energy lepton colliders.

##  V Conclusions

In this article we have studied two simple SM extensions featuring, each of
them, a new vector-like lepton and a sterile neutrino. We have confronted the
cases of a doublet and a singlet VLL and discussed their collider
phenomenology, both at the LHC and future leptonic machines. For the former we
have employed Deep Learning techniques to compute the statistical significance
of a hypothetical discovery. In the selection of benchmark scenarios we have
required that the coupling of VLL to muons is consistent with flavour
constraints and the branching fractions of Higgs and Z bosons to muons. We
have also shown that, decay width of both VLL doublet components becomes
increasingly larger with growing exotic lepton masses. In our analysis we have
only considered scenarios where the width is smaller than the mass.

In the context of the LHC studies, we have performed Monte-Carlo simulations
to generate data for signal and background topologies. The signal is
characterized by the presence of a single isolated lepton and a substantial
amount of MET for the singlet case. For the doublet case, the signal involves
4 light jets, a charged lepton and a neutrino as final state particles. To
separate the signal from the background we constructed neural networks which
follow from an implementation of an evolution algorithm that maximises the
Asimov significance. We have shown that, for doublet VLLs, we can exclude
masses up to 1 TeV with more than five standard deviations. In particular, for
masses of $1\leavevmode\nobreak\ \mathrm{TeV}$, we obtain a significance of
$\mathcal{Z}(<1\%)=11.14\sigma$ for the high luminosity phase of the LHC, with
an integrated luminosity of $\mathcal{L}=3000$ $\mathrm{fb^{-1}}$. We have
also verified that one can already test the doublet VLL scenario at the Run
III of the LHC, which will deliver $\mathcal{L}=300$ $\mathrm{fb^{-1}}$ of
data. For such luminosity, one can test VLL masses up to about
$300\leavevmode\nobreak\ \mathrm{GeV}$ obtaining $s/\sqrt{s+b}=8.02\sigma$,
$\mathcal{Z}(<1\%)=8.09\sigma$ and $\mathcal{Z}_{A}=2.56\sigma$. For the
singlet scenario, production cross-sections are substantially smaller and the
expected number of events is usually zero, with an exception for VLL masses of
$100\leavevmode\nobreak\ \mathrm{GeV}$ and $200\leavevmode\nobreak\
\mathrm{GeV}$. While for the former one finds $s/\sqrt{s+b}=2.98\sigma$,
$\mathcal{Z}(<1\%)=4.22\sigma$ and $\mathcal{Z}_{A}=0.0023\sigma$ at
$\mathcal{L}=3000$ $\mathrm{fb^{-1}}$, the run III estimation gives
$s/\sqrt{s+b}=0.94\sigma$, $\mathcal{Z}(<1\%)=1.33\sigma$ and
$\mathcal{Z}_{A}=0.00094\sigma$ at $\mathcal{L}=300$ $\mathrm{fb^{-1}}$,
experiencing a great drop for a mass of $200\leavevmode\nobreak\
\mathrm{GeV}$. In particular, the latter scenario can not be excluded, with a
statistical significance bellow the discovery threshold.

Owning to the low production cross-sections of the singlet scenario, a
supplementary analysis was made within the context of lepton colliders. We
have performed numerical computations for the expected VLL pair-production
cross-section in the $s$-channel. We find that, in general, larger cross-
sections are obtained when compared to the LHC analysis, allowing for a much
wider range of masses to be probed with relevance for singlet VLLs. In
particular, we note that for luminosities of the order of $\mathrm{ab}^{-1}$
the, even for singlet VLLs with mass $3.7\leavevmode\nobreak\ \mathrm{TeV}$
one can expect $260$ events at $E_{\mathrm{CM}}=10$ TeV and $180$ events if
$E_{\mathrm{CM}}=14$ TeV. Furthermore, even for a VLL as heavy as
$6\leavevmode\nobreak\ \mathrm{TeV}$, a $14\leavevmode\nobreak\ \mathrm{TeV}$
lepton machine delivering a luminosity of $1\leavevmode\nobreak\
\mathrm{ab}^{-1}$ is expected to produce almost 100 events. With this in mind
we conclude that the study of VLL particles at future electron or muon
colliders is a rather relevant physics case scenario to be explored, allowing
to significantly extend the current reach of the LHC.

Acknowledgements. The authors would like to thank Celso Nishi for valuable
comments made to the initial version of this manuscript. J.G., F.F.F., and
A.P.M. are supported by the Center for Research and Development in Mathematics
and Applications (CIDMA) through the Portuguese Foundation for Science and
Technology (FCT - Fundação para a Ciência e a Tecnologia), references
UIDB/04106/2020 and UIDP/04106/2020. A.P.M., F.F.F., J.G. and R.S. are
supported by the project PTDC/FIS-PAR/31000/2017. A.P.M., F.F.F., J.G. are
also supported by the projects CERN/FIS-PAR/0021/2021. Additionally, A.P.M.
and J.G. are supported by CERN/FIS-PAR/0019/2021. J.G. is also directly funded
by FCT through the doctoral program grant with the reference 2021.04527.BD.
A.P.M. is also supported by national funds (OE), through FCT, I.P., in the
scope of the framework contract foreseen in the numbers 4, 5 and 6 of the
article 23, of the Decree-Law 57/2016, of August 29, changed by Law 57/2017,
of July 19. R.P. is supported in part by the Swedish Research Council grant,
contract number 2016-05996, as well as by the European Research Council (ERC)
under the European Union’s Horizon 2020 research and innovation programme
(grant agreement No 668679). R.S. is supported by CFTC-UL under FCT contracts
UIDB/00618/2020, UIDP/00618/2020, and by the projects CERN/FISPAR /0002/2017,
CERN/FIS-PAR/0014/2019 and by the HARMONIA project of the National Science
Centre, Poland, under contract UMO- 2015/18/M/ST2/00518. A.O. is supported by
the FCT project CERN/FIS-PAR/0029/2019.

##  Appendix A Feynman Rules

In this appendix, the list of Feynman rules relevant for the numerical
analysis is presented. All rules are displayed in the mass basis. To simplify
the notation, we define the chirality projection operators as
$\mathrm{P_{L}}=(1-\gamma_{5})/2$ and $\mathrm{P_{R}}=(1+\gamma_{5})/2$.
Notation wise, we define the following: $\theta_{W}$ is the Weinberg mixing
angle, $g$ the $\mathrm{SU}(2)_{\mathrm{L}}$ gauge coupling and $g^{\prime}$
the $\mathrm{U}(1)_{\mathrm{Y}}$ gauge coupling. The Latin indices
$i,j=(1,2,3,4)$ denote $(e,\mu,\tau,E(\mathcal{E}))$ for leptons, while for
neutrinos
$i,j=(1,2,3,4,5,6)=(\nu_{1},\nu_{2},\nu_{3},\nu_{4},\nu_{5},\nu_{6})$, with
$\nu_{5}$ and $\nu_{6}$ arising from the doublet representation of the VLL and
$\nu_{4}$ is the sterile neutrino. For the singlet model, there is no extra
left-handed neutrinos, hence the latin indices only run from 1 to 4.

For the doublet case, the Feynman rules read as

$\displaystyle\hskip 57.18988pt\raisebox{-33.5001pt}{\includegraphics{ZBoson-
doublet.pdf}}=i\quantity(g^{Z}_{L}\gamma_{\mu}\mathrm{P_{L}}+g^{Z}_{R}\gamma_{\mu}\mathrm{P_{R}}),$
(19) $\displaystyle
g_{L}^{Z}=\frac{1}{2}\quantity(g\cos{\theta_{W}}-g^{\prime}\sin{\theta_{W}})\quantity[\quantity(U_{L}^{e})^{*}_{j4}\quantity(U_{L}^{e})_{i4}+\sum_{a=1}^{3}\quantity(U_{L}^{e})^{*}_{ja}\quantity(U_{L}^{e})_{ia}],$
$\displaystyle
g_{R}^{Z}=-g^{\prime}\sin(\theta_{W})\sum_{a=1}^{3}\quantity(U_{R}^{e})^{*}_{ia}\quantity(U_{R}^{e})_{ja}+\frac{1}{2}\quantity(U_{R}^{e})^{*}_{i4}\quantity(U_{R}^{e})_{j4}\quantity(g\cos{\theta_{W}}-g^{\prime}\sin{\theta_{W}}),$

$\displaystyle\hskip 7.11317pt\raisebox{-33.5001pt}{\includegraphics{Wplus-
doublet.pdf}}=i\quantity(g^{W^{+}}_{L}\gamma_{\mu}\mathrm{P_{L}}+g^{W^{+}}_{R}\gamma_{\mu}\mathrm{P_{R}}),$
(20) $\displaystyle
g_{L}^{W^{+}}=-\frac{g}{\sqrt{2}}\quantity[\quantity(U_{L}^{e})^{*}_{j4}\quantity(U_{\nu})_{i5}+\sum_{a=1}^{3}\quantity(U_{L}^{e})^{*}_{ja}\quantity(U_{\nu})_{ia}],$
$\displaystyle
g_{R}^{W^{+}}=-\frac{g}{\sqrt{2}}\quantity(U_{\nu})^{*}_{i6}\quantity(U_{R}^{e})_{j4},$

$\displaystyle\hskip 7.11317pt\raisebox{-33.5001pt}{\includegraphics{Wminus-
doublet.pdf}}=i\quantity(g^{W^{-}}_{L}\gamma_{\mu}\mathrm{P_{L}}+g^{W^{-}}_{R}\gamma_{\mu}\mathrm{P_{R}}),$
(21) $\displaystyle
g_{L}^{W^{-}}=-\frac{g}{\sqrt{2}}\quantity[\quantity(U_{\nu})^{*}_{i5}\quantity(U_{L}^{e})_{j4}+\sum_{a=1}^{3}\quantity(U_{\nu})^{*}_{ia}\quantity(U_{L}^{e})_{ja}],$
$\displaystyle
g_{R}^{W^{-}}=-\frac{g}{\sqrt{2}}\quantity(U_{R}^{e})^{*}_{j4}\quantity(U_{\nu})_{i6},$

$\displaystyle\hskip 24.18501pt\raisebox{-33.5001pt}{\includegraphics{Higgs-
doublet.pdf}}=i\quantity(g^{H}_{L}\mathrm{P_{L}}+g^{H}_{R}\mathrm{P_{R}}),$
(22) $\displaystyle
g_{L}^{H}=\frac{1}{\sqrt{2}}\quantity[\sum_{a=1}^{3}\quantity(U_{L}^{e})^{*}_{j4}\quantity(U_{R}^{e})^{*}_{ia}\Theta_{a}+\sum_{a=1}^{3}\sum_{b=1}^{3}\quantity(U_{R}^{e})^{*}_{ib}\quantity(U_{L}^{e})^{*}_{ja}\Pi_{ab}],$
$\displaystyle
g_{R}^{H}=\frac{1}{\sqrt{2}}\quantity[\sum_{a=1}^{3}\quantity(U_{R}^{e})_{ja}\quantity(U_{L}^{e})_{i4}\Theta_{a}^{*}+\sum_{a=1}^{3}\sum_{b=1}^{3}\quantity(U_{L}^{e})_{ia}\quantity(U_{R}^{e})_{jb}\Pi_{ab}^{*}],$

Similarly, for the singlet case we have,

$\displaystyle\hskip 55.48277pt\raisebox{-33.5001pt}{\includegraphics{ZBoson-
doublet.pdf}}=i\quantity(g^{Z}_{L}\gamma_{\mu}\mathrm{P_{L}}+g^{Z}_{R}\gamma_{\mu}\mathrm{P_{R}}),$
(23) $\displaystyle
g_{L}^{Z}=-g^{\prime}\quantity(U^{e}_{L})^{*}_{j4}\quantity(U^{e}_{L})_{i4}\sin(\theta_{W})+\frac{1}{2}\quantity(g\cos{\theta_{W}}-g^{\prime}\sin{\theta_{W}})\sum_{a=1}^{3}\quantity(U_{L}^{e})^{*}_{ja}\quantity(U_{L}^{e})_{ia},$
$\displaystyle
g_{R}^{Z}=-g^{\prime}\sin(\theta_{W})\quantity[\quantity(U_{R}^{e})^{*}_{i4}\quantity(U_{R}^{e})_{j4}+\sum_{a=1}^{3}\quantity(U_{R}^{e})^{*}_{ia}\quantity(U_{R}^{e})_{ja}],$

$\displaystyle\raisebox{-33.5001pt}{\includegraphics{Wplus-
doublet.pdf}}=ig^{W^{+}}_{L}\gamma_{\mu}\mathrm{P_{L}},$ (24) $\displaystyle
g_{L}^{W^{+}}=-\frac{g}{\sqrt{2}}\sum_{a=1}^{3}\quantity(U_{L}^{e})^{*}_{ja}\quantity(U_{\nu})_{ia},$

$\displaystyle\raisebox{-33.5001pt}{\includegraphics{Wminus-
doublet.pdf}}=ig^{W^{-}}_{L}\gamma_{\mu}\mathrm{P_{L}},$ (25) $\displaystyle
g_{L}^{W^{-}}=-\frac{g}{\sqrt{2}}\sum_{a=1}^{3}\quantity(U_{\nu})^{*}_{ia}\quantity(U_{L}^{e})_{ja},$

$\displaystyle\hskip 22.76228pt\raisebox{-33.5001pt}{\includegraphics{Higgs-
doublet.pdf}}=i\quantity(g^{H}_{L}\mathrm{P_{L}}+g^{H}_{R}\mathrm{P_{R}}),$
(26) $\displaystyle
g_{L}^{H}=\frac{1}{\sqrt{2}}\quantity[\sum_{a=1}^{3}\quantity(U_{R}^{e})^{*}_{i4}\quantity(U_{L}^{e})^{*}_{ja}\theta_{a}+\sum_{a=1}^{3}\sum_{b=1}^{3}\quantity(U_{R}^{e})^{*}_{ib}\quantity(U_{L}^{e})^{*}_{ja}\Pi_{ab}],$
$\displaystyle
g_{R}^{H}=\frac{1}{\sqrt{2}}\quantity[\sum_{a=1}^{3}\quantity(U_{L}^{e})_{ia}\quantity(U_{R}^{e})_{j4}\theta_{a}^{*}+\sum_{a=1}^{3}\sum_{b=1}^{3}\quantity(U_{L}^{e})_{ia}\quantity(U_{R}^{e})_{jb}\Pi_{ab}^{*}].$

##  Appendix B Kinematic distributions for VLL single-production

Figure 10: A sample of the kinematic variables that the neural network uses
for classification. Data is simulated for a doublet VLL with mass of 700 GeV
and in the ATLAS detector. From left to right and top to bottom, the variables
are $p_{T}$ of the muon, $p_{T}$ of jet 1, 2, 3 and 4, mass distributions for
the combinations $(j_{1},j_{2},\mathrm{MET})$, $(j_{1},j_{3},\mathrm{MET})$,
$(j_{1},j_{4},\mathrm{MET})$, $(j_{2},j_{3},\mathrm{MET})$,
$(j_{2},j_{4},\mathrm{MET})$, $(j_{3},j_{4},\mathrm{MET})$, $\Delta\theta$
distributions for the combinations $(j_{1},j_{3})$, $(j_{1},j_{4})$,
$(j_{2},j_{3})$, $(j_{2},j_{4})$ and $(j_{3},j_{4})$, $\Delta\phi$
distributions for the combinations $(j_{1},j_{2})$, $(j_{1},j_{3})$,
$(j_{1},j_{4})$, $(j_{2},j_{3})$, $(j_{2},j_{4})$ and $\Delta R$ distributions
for the combinations $(j_{1},j_{2})$, $(j_{1},j_{3})$, $(j_{1},j_{4})$,
$(j_{2},j_{3})$ and $(j_{2},j_{4})$. In the $y$-axis, it is indicated that
events are normalized (NE). We consider 30 bins for all background and signal
histograms. Figure 11: Kinematic variables that the neural network uses for
classification. Data is simulated for a doublet VLL with mass of 700 GeV and
in the ATLAS detector. From left to right and top to bottom, the variables are
$\cos(\theta_{\mu^{-}})$, transverse momentum of the muon, pseudo-rapidity of
the muon, azimuthal angle of the muon and MET. In the $y$-axis, it is
indicated that events are normalized (NE). We consider 30 bins for all
background and signal histograms.

##  Appendix C Neural networks found via the genetic algorithm

Mass | Doublet
---|---
100 GeV |  Layers : 1 input + 1 hidden +1 output. Input and hidden layer with 1024 neurons, output layer with 4 neurons Regularizer : L1L2 (layer 1 and 2) and none (layer 3) Activation function : RELU (for layers 1 and 2) and Sigmoid (layer 3) Initializer : VarianceScaling, with uniform distribution (layer 1 and 2) in fan_in mode and GlorotUniform (layer 3) with null seed
200 GeV |  Layers : 1 input + 3 hidden + 1 output. Hidden and input layers with 1024 neurons each, output layer with 4 neurons Regularizer : L1L2 (layers 1 to 4) and none (for layer 5) Activation function : RELU (layers 1 to 4) and Sigmoid (layer 5) Initializer : VarianceScaling, with uniform distribution (layers 1-4) in fan_in mode and GlorotUniform (layer 5) with null seed
300 GeV |  Layers : 1 input + 4 hidden + 1 output. Input and hidden layers with 512 neurons, output layer with 4 neurons Regularizer : L1L2 (for layers 1-5) and none (for layer 6) Activation function : Sigmoid (for all layers) Initializer : VarianceScaling, with truncated normal distribution (for layers 1-5) in fan_in mode and GlorotUniform distribution (layer 6) with null seed
400 GeV |  Layers : 1 input + 2 hidden + 1 output. Hidden and input layers with 256 neurons each, output layer with 4 neurons Regularizer : L1L2 (for layers 1-3) and none (for layer 4) Activation function : RELU (layers 1 to 3) and Sigmoid (layer 4) Initializer : RandomNormal, with 0 mean and 0.05 standard deviation (for layers 1 to 3) and GlorotUniform distribution (layer 4) with null seed
500 GeV |  Layers : 1 input + 1 output. Input layer with 1024 neurons, output layer with 4 neurons Regularizer : L1L2 (layer 1) and none (layer 2) Activation function : Sigmoid (all layers) Initializer : VarianceScaling, with truncated normal distribution (layer 1) in fan_in mode and GlorotUniform distribution (layer 2) with null seed

Table 8: Architectures for the neural networks that the genetic algorithm
finds for the scanned masses (100 to 500 GeV) of the doublet VLL scenario.

Mass | Doublet
---|---
600 GeV |  Layers : 1 input + 4 hidden + 1 output. Hidden and input layers with 512 neurons each, output layer with 4 neurons Regularizer : L1L2 (for layers 1 to 5) and none (for layer 6) Activation function : Sigmoid (for all layers) Initializer : RandomNormal, with mean 0 and 0.05 standard deviation (for layers 1 to 5) and GlorotUniform distribution (layer 6) with null seed
700 GeV |  Layers : 1 input + 3 hidden + 1 output. Hidden and input layers with 512 neurons each, output layer with 4 neurons Regularizer : L1L2 (for layers 1 to 4) and none (for layer 5) Activation function : Sigmoid (for all layers) Initializer : VarianceScaling, with uniform distribution (for layers 1 to 4) in fan_in mode and GlorotUniform distribution (layer 5) with null seed
800 GeV |  Layers : 1 input + 3 hidden + 1 output. Hidden and input layers with 256 neurons each, output layer with 4 neurons Regularizer : L1L2 (for layers 1 and 4) and none (for layer 5) Activation function : ELU (for layers 1 to 4) and Sigmoid (layer 5) Initializer : VarianceScaling, with truncated normal distribution (for layers 1 to 4) in fan_in mode and GlorotUniform distribution (layer 5) with null seed
900 GeV |  Layers : 1 input + 1 hidden + 1 output. Hidden and input layers with 256 neurons each, output layer with 4 neurons Regularizer : L1L2 (for layers 1 and 2) and none (for layer 3) Activation function : Tanh (for layers 1 and 2) and Sigmoid (layer 3) Initializer : RandomNormal, with mean 0 and 0.05 standard deviation (for layers 1 and 2) and GlorotUniform distribution (layer 3) with null seed
1000 GeV |  Layers : 1 input + 2 hidden + 1 output. Hidden and input layers with 2048 neurons each, output layer with 4 neurons Regularizer : L1L2 (for layers 1 to 3) and none (for layer 4) Activation function : ELU (for layers 1 to 3) and Sigmoid (layer 4) Initializer : VarianceScaling, with uniform distribution (for layers 1 to 3) in fan_in mode and GlorotUniform distribution (layer 4) with null seed

Table 9: Architectures for the neural networks that the genetic algorithm
finds for the scanned masses (600 to 1000 GeV) of the doublet VLL scenario.

Mass | Singlet
---|---
100 GeV |  Layers : 1 input + 1 hidden + 1 output. Hidden and input layers with 512 neurons each, output layer with 3 neurons Regularizer : L1L2 (for layers 1 to 4) and none (for layer 5) Activation function : sigmoid (for all layers) Initializer : RandomNormal, with 0 mean and 0.05 standard deviation (for layers 1 and 2) and VarianceScaling, with uniform distribution in fan_in mode (last layer)
200 GeV |  Layers : 1 input + 1 hidden + 1 output. Hidden and input layers with 2048 neurons each, output layer with 3 neurons Regularizer : L1L2 (for layers 1 and 2) and none (for layer 3) Activation function : sigmoid (for all layers) Initializer : VarianceScaling, with uniform distribution (layer 1 and 2) in fan_in mode and fan_avg mode (last layer)

Table 10: Architectures for the neural networks that the genetic algorithm
finds for the scanned masses (100 and 200 GeV) of the singlet VLL scenario.

##  Appendix D ROC plots for Doublet VLLs

(a)

(b)

(c)

(d)

(e)

Figure 12: Receiver operator characteristic (ROC) plots for the neural
networks of the doublet VLLs from the masses of 100 GeV (top) to 500 GeV
(bottom).

(a)

(b)

(c)

(d)

(e)

Figure 13: Receiver operator characteristic (ROC) plots for the neural
networks of the doublet VLLs from the masses of 600 GeV (top) to 1000 GeV
(bottom).

## References

  * Fukuda  _et al._ (1998) Y. Fukuda _et al._ (Super-Kamiokande), “Evidence for oscillation of atmospheric neutrinos,” Phys. Rev. Lett. 81, 1562–1567 (1998), arXiv:hep-ex/9807003 .
  * Bertone  _et al._ (2005) Gianfranco Bertone, Dan Hooper, and Joseph Silk, “Particle dark matter: Evidence, candidates and constraints,” Phys. Rept. 405, 279–390 (2005), arXiv:hep-ph/0404175 .
  * Abi  _et al._ (2021) B. Abi _et al._ (Muon g-2), “Measurement of the Positive Muon Anomalous Magnetic Moment to 0.46 ppm,” Phys. Rev. Lett. 126, 141801 (2021), arXiv:2104.03281 [hep-ex] .
  * Aaij  _et al._ (2021) Roel Aaij _et al._ (LHCb), “Test of lepton universality in beauty-quark decays,” (2021), arXiv:2103.11769 [hep-ex] .
  * Morais  _et al._ (2020) António P. Morais, Roman Pasechnik, and Werner Porod, “Prospects for new physics from gauge left-right-colour-family grand unification hypothesis,” Eur. Phys. J. C 80, 1162 (2020), arXiv:2001.06383 [hep-ph] .
  * de Anda  _et al._ (2020) Francisco J. de Anda, Alfredo Aranda, António P. Morais, and Roman Pasechnik, “Gauge couplings evolution from the Standard Model, through Pati-Salam theory, into $E_{8}$ unification of families and forces,” (2020), arXiv:2011.13902 [hep-ph] .
  * Aranda  _et al._ (2021a) Alfredo Aranda, Francisco J. de Anda, António P. Morais, and Roman Pasechnik, “A Different Take on $E_{8}$ Unification,” (2021a), arXiv:2107.05421 [hep-ph] .
  * Aranda  _et al._ (2021b) Alfredo Aranda, Francisco J. de Anda, António P. Morais, and Roman Pasechnik, “Sculpting the Standard Model from low-scale Gauge-Higgs-Matter $E_{8}$ Grand Unification in ten dimensions,” (2021b), arXiv:2107.05495 [hep-ph] .
  * Freitas  _et al._ (2021) Felipe F. Freitas, João Gonçalves, António P. Morais, and Roman Pasechnik, “Phenomenology of vector-like leptons with Deep Learning at the Large Hadron Collider,” JHEP 01, 076 (2021), arXiv:2010.01307 [hep-ph] .
  * Yin and Yamaguchi (2020) Wen Yin and Masahiro Yamaguchi, “Muon $g-2$ at multi-TeV muon collider,” (2020), arXiv:2012.03928 [hep-ph] .
  * Crivellin  _et al._ (2020) Andreas Crivellin, Fiona Kirk, Claudio Andrea Manzari, and Marc Montull, “Global Electroweak Fit and Vector-Like Leptons in Light of the Cabibbo Angle Anomaly,” JHEP 12, 166 (2020), arXiv:2008.01113 [hep-ph] .
  * Bißmann  _et al._ (2021) Stefan Bißmann, Gudrun Hiller, Clara Hormigos-Feliu, and Daniel F. Litim, “Multi-lepton signatures of vector-like leptons with flavor,” Eur. Phys. J. C 81, 101 (2021), arXiv:2011.12964 [hep-ph] .
  * Hiller  _et al._ (2020a) Gudrun Hiller, Clara Hormigos-Feliu, Daniel F. Litim, and Tom Steudtner, “Anomalous magnetic moments from asymptotic safety,” Phys. Rev. D 102, 071901 (2020a), arXiv:1910.14062 [hep-ph] .
  * Hiller  _et al._ (2020b) Gudrun Hiller, Clara Hormigos-Feliu, Daniel F. Litim, and Tom Steudtner, “Model Building from Asymptotic Safety with Higgs and Flavor Portals,” Phys. Rev. D 102, 095023 (2020b), arXiv:2008.08606 [hep-ph] .
  * Guedes and Santiago (2021) Guilherme Guedes and José Santiago, “New leptons with exotic decays: collider limits and dark matter complementarity,” (2021), arXiv:2107.03429 [hep-ph] .
  * Eberhardt  _et al._ (2012) Otto Eberhardt, Geoffrey Herbert, Heiko Lacker, Alexander Lenz, Andreas Menzel, Ulrich Nierste, and Martin Wiebusch, “Impact of a Higgs boson at a mass of 126 GeV on the standard model with three and four fermion generations,” Phys. Rev. Lett. 109, 241802 (2012), arXiv:1209.1101 [hep-ph] .
  * Zyla  _et al._ (2020) P. A. Zyla _et al._ (Particle Data Group), “Review of Particle Physics,” PTEP 2020, 083C01 (2020).
  * Strumia and Vissani (2006) Alessandro Strumia and Francesco Vissani, “Neutrino masses and mixings and…” (2006), arXiv:hep-ph/0606054 .
  * Achard  _et al._ (2001) P. Achard _et al._ (L3), “Search for heavy neutral and charged leptons in $e^{+}e^{-}$ annihilation at LEP,” Phys. Lett. B 517, 75–85 (2001), arXiv:hep-ex/0107015 .
  * Sirunyan  _et al._ (2019) Albert M Sirunyan _et al._ (CMS), “Search for vector-like leptons in multilepton final states in proton-proton collisions at $\sqrt{s}$ = 13 TeV,” Phys. Rev. D 100, 052003 (2019), arXiv:1905.10853 [hep-ex] .
  * Aad  _et al._ (2015) Georges Aad _et al._ (ATLAS), “Search for heavy lepton resonances decaying to a $Z$ boson and a lepton in $pp$ collisions at $\sqrt{s}=8$ TeV with the ATLAS detector,” JHEP 09, 108 (2015), arXiv:1506.01291 [hep-ex] .
  * Chatrchyan  _et al._ (2012) Serguei Chatrchyan _et al._ (CMS), “Search for Heavy Lepton Partners of Neutrinos in Proton-Proton Collisions in the Context of the Type III Seesaw Mechanism,” Phys. Lett. B 718, 348–368 (2012), arXiv:1210.1797 [hep-ex] .
  * Esteban  _et al._ (2020) Ivan Esteban, M. C. Gonzalez-Garcia, Michele Maltoni, Thomas Schwetz, and Albert Zhou, “The fate of hints: updated global analysis of three-flavor neutrino oscillations,” JHEP 09, 178 (2020), arXiv:2007.14792 [hep-ph] .
  * del Aguila  _et al._ (2008) F. del Aguila, J. de Blas, and M. Perez-Victoria, “Effects of new leptons in Electroweak Precision Data,” Phys. Rev. D 78, 013010 (2008), arXiv:0803.4008 [hep-ph] .
  * Kuksa and Volchanskiy (2013) V. Kuksa and N. Volchanskiy, “Factorization in the model of unstable particles with smeared mass,” Central Eur. J. Phys. 11, 182–194 (2013), arXiv:1109.1541 [hep-ph] .
  * Kuksa (2006) V. I. Kuksa, “Convolution formula for decay rate,” Phys. Lett. B 633, 545–549 (2006), [Erratum: Phys.Lett.B 664, 315 (2008)], arXiv:hep-ph/0508164 .
  * Kuksa (2009) V. I. Kuksa, “Factorization method in the model of unstable particles with a smeared mass,” Phys. Atom. Nucl. 72, 1063–1073 (2009), arXiv:0902.4892 [hep-ph] .
  * Thomas and Wells (1998) Scott D. Thomas and James D. Wells, “Phenomenology of Massive Vectorlike Doublet Leptons,” Phys. Rev. Lett. 81, 34–37 (1998), arXiv:hep-ph/9804359 .
  * Porod and Staub (2012) W. Porod and F. Staub, “SPheno 3.1: Extensions including flavour, CP-phases and models beyond the MSSM,” Comput. Phys. Commun. 183, 2458–2469 (2012), arXiv:1104.1573 [hep-ph] .
  * Straub (2018) David M. Straub, “flavio: a Python package for flavour and precision phenomenology in the Standard Model and beyond,” (2018), arXiv:1810.08132 [hep-ph] .
  * Freitas  _et al._ (2014) Ayres Freitas, Joseph Lykken, Stefan Kell, and Susanne Westhoff, “Testing the Muon g-2 Anomaly at the LHC,” JHEP 05, 145 (2014), [Erratum: JHEP 09, 155 (2014)], arXiv:1402.7065 [hep-ph] .
  * Dermisek  _et al._ (2021) Radovan Dermisek, Keith Hermanek, and Navin McGinnis, “Muon $g-2$ in two Higgs doublet models with vectorlike leptons,” (2021), arXiv:2103.05645 [hep-ph] .
  * Kannike  _et al._ (2012) Kristjan Kannike, Martti Raidal, David M. Straub, and Alessandro Strumia, “Anthropic solution to the magnetic muon anomaly: the charged see-saw,” JHEP 02, 106 (2012), [Erratum: JHEP 10, 136 (2012)], arXiv:1111.2551 [hep-ph] .
  * Jegerlehner and Nyffeler (2009) Fred Jegerlehner and Andreas Nyffeler, “The Muon g-2,” Phys. Rept. 477, 1–110 (2009), arXiv:0902.3360 [hep-ph] .
  * Dermisek and Raval (2013) Radovan Dermisek and Aditi Raval, “Explanation of the Muon g-2 Anomaly with Vectorlike Leptons and its Implications for Higgs Decays,” Phys. Rev. D 88, 013017 (2013), arXiv:1305.3522 [hep-ph] .
  * Crivellin  _et al._ (2018) Andreas Crivellin, Martin Hoferichter, and Philipp Schmidt-Wellenburg, “Combined explanations of $(g-2)_{\mu,e}$ and implications for a large muon EDM,” Phys. Rev. D 98, 113002 (2018), arXiv:1807.11484 [hep-ph] .
  * Bonilla  _et al._ (2021) Cesar Bonilla, A. E. Cárcamo Hernández, João Gonçalves, Felipe F. Freitas, Antonio P. Morais, and R. Pasechnik, “Collider signatures of vector-like fermions from a flavor symmetric 2HDM,” (2021), arXiv:2107.14165 [hep-ph] .
  * Staub (2014) Florian Staub, “SARAH 4 : A tool for (not only SUSY) model builders,” Comput. Phys. Commun. 185, 1773–1790 (2014), arXiv:1309.7223 [hep-ph] .
  * 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 (2014), arXiv:1405.0301 [hep-ph] .
  * Sjöstrand  _et al._ (2015) Torbjörn Sjöstrand, Stefan Ask, Jesper R. Christiansen, Richard Corke, Nishita Desai, Philip Ilten, Stephen Mrenna, Stefan Prestel, Christine O. Rasmussen, and Peter Z. Skands, “An introduction to PYTHIA 8.2,” Comput. Phys. Commun. 191, 159–177 (2015), arXiv:1410.3012 [hep-ph] .
  * Ball  _et al._ (2013) Richard D. Ball, Valerio Bertone, Stefano Carrazza, Luigi Del Debbio, Stefano Forte, Alberto Guffanti, Nathan P. Hartland, and Juan Rojo (NNPDF), “Parton distributions with QED corrections,” Nucl. Phys. B 877, 290–320 (2013), arXiv:1308.0598 [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 (2014), arXiv:1307.6346 [hep-ex] .
  * Brun and Rademakers (1997) R. Brun and F. Rademakers, “ROOT: An object oriented data analysis framework,” Nucl. Instrum. Meth. A 389, 81–86 (1997).
  * Chollet  _et al._ (2015) François Chollet _et al._ , “Keras,” https://keras.io (2015).
  * Elwood and Krücker (2018) Adam Elwood and Dirk Krücker, “Direct optimisation of the discovery significance when training neural networks to search for new physics in particle colliders,” (2018), arXiv:1806.00322 [hep-ex] .
  * Cowan  _et al._ (2011) Glen Cowan, Kyle Cranmer, Eilam Gross, and Ofer Vitells, “Asymptotic formulae for likelihood-based tests of new physics,” Eur. Phys. J. C 71, 1554 (2011), [Erratum: Eur.Phys.J.C 73, 2501 (2013)], arXiv:1007.1727 [physics.data-an] .
  * Chawla  _et al._ (2011) N. V. Chawla, K. W. Bowyer, L. O. Hall, and W. P. Kegelmeyer, “SMOTE: Synthetic Minority Over-sampling Technique,” arXiv e-prints , arXiv:1106.1813 (2011), arXiv:1106.1813 [cs.AI] .
  * Charles  _et al._ (2018) T. K. Charles _et al._ (CLICdp, CLIC), “The Compact Linear Collider (CLIC) - 2018 Summary Report,”  2/2018 (2018), 10.23731/CYRM-2018-002, arXiv:1812.06018 [physics.acc-ph] .
  * Boland  _et al._ (2016) M J Boland _et al._ (CLIC, CLICdp), “Updated baseline for a staged Compact Linear Collider,”  (2016), 10.5170/CERN-2016-004, arXiv:1608.07537 [physics.acc-ph] .
  * Bae (2013) “The International Linear Collider Technical Design Report - Volume 2: Physics,” (2013), arXiv:1306.6352 [hep-ph] .
  * Shiltsev (2021) Vladimir Shiltsev, “Muon colliders –opening new horizons for particle physics,” https://indico.cern.ch/event/1022802/contributions/4293870/attachments/2218112/3755618/MuonCollider_UChi_03292021_v1.pdf (2021).
  * Muo (2021) “Muon collider forum kickoff meeting,” https://indico.fnal.gov/event/47038/ (2021).
  * Shtabovenko  _et al._ (2020) Vladyslav Shtabovenko, Rolf Mertig, and Frederik Orellana, “FeynCalc 9.3: New features and improvements,” Comput. Phys. Commun. 256, 107478 (2020), arXiv:2001.04407 [hep-ph] .
  * CER (2018) “Perspectives on the determination of systematic uncertainties at hl-lhc,” https://cds.cern.ch/record/2642427/files/ATL-PHYS-SLIDE-2018-893.pdf (2018).
  * ATL (2022) “Luminosity determination in $pp$ collisions at $\sqrt{s}=13$ TeV using the ATLAS detector at the LHC,” (2022), arXiv:2212.09379 [hep-ex] .
  * Budnev  _et al._ (1975) V. M. Budnev, I. F. Ginzburg, G. V. Meledin, and V. G. Serbo, “The Two photon particle production mechanism. Physical problems. Applications. Equivalent photon approximation,” Phys. Rept. 15, 181–281 (1975).
  * Frixione  _et al._ (1993) Stefano Frixione, Michelangelo L. Mangano, Paolo Nason, and Giovanni Ridolfi, “Improving the Weizsacker-Williams approximation in electron - proton collisions,” Phys. Lett. B 319, 339–345 (1993), arXiv:hep-ph/9310350 .