# Search for single vector-like $B$ quark production
in hadronic final states at the LHC

Bingfang Yang1111Email<EMAIL_ADDRESS>Zejun Li1222Email:
<EMAIL_ADDRESS>Xinglong Jia1333Email<EMAIL_ADDRESS>Stefano Moretti2,3444Email<EMAIL_ADDRESS><EMAIL_ADDRESS>and Liangliang Shang1,2555Email:
<EMAIL_ADDRESS><EMAIL_ADDRESS>1School of Physics,
Henan Normal University, Xinxiang 453007, PR China
2Department of Physics and Astronomy, Uppsala University, Box 516, SE-751 20
Uppsala, Sweden
3School of Physics and Astronomy, University of Southampton, Highfield,
Southampton SO17 1BJ, UK

###### Abstract

In this paper, we study the discovery potential of a Vector-Like $B$ quark
(VLB) via the process $pp\to B(\to bZ)j\to b(Z\to\nu_{l}\bar{\nu_{l}})j$ at
the Large Hadron Collider (LHC) with $\sqrt{s}=14$ TeV. In the framework of a
simplified model, we perform a scan over its parameter space and test its
viability following a Monte Carlo analysis developed to include all production
and decay dynamics. We use cut-and-count combined with Extreme Gradient
Boosting (XGBoost) methods to classify the signal and background events in
order to improve the efficiency of signal identification and background
rejection. We find that this approach can reduce background events
significantly while the signal retention rate is much higher than that of
traditional methods, thereby improving the VLB discovery potential. We then
calculate the exclusion and discovery capabilities for VLBs and find that the
advantages of the cut-and-count plus XGBoost method especially lie in the
high-mass region, i.e., $m_{B}>1500\text{ GeV}$. We finally obtain the
following LHC results in terms of the coupling and chiral structure of a
singlet heavy VLB interactions: (i) for $g^{\ast}$=0.2 and $R_{L}=0$ with 3000
fb-1, the $B$ quark mass can be be excluded (discovered) up to 3000 GeV (2500
GeV); (ii) for $g^{\ast}$=0.2 and $R_{L}=0.5$ with 3000 fb-1, the exclusion
(discovery) region can reach up to 4750 GeV (4250 GeV).

##  I Introduction

The gauge hierarchy problem Wilson:1970ag  caused by the quadratic divergence
of the Higgs mass is one of the unsolved problems in the Standard Model (SM)
of particle physics. For the past few decades, solving this problem has become
a guiding principle for establishing New Physics (NP) models. Amongst these,
Vector-Like Quarks (VLQs) are often introduced to alleviate this problem since
an extra SM-like chiral fermion has been disfavored by the Higgs data
collected at the Large Hadron Collider (LHC) ATLAS:2019nkf ; CMS:2018gwt .
VLQs are a type of hypothetical fermions with spin 1/2 and color-triplets,
which left- and right-handed components transform with the same properties
under the SM Electro-Weak (EW) symmetry group DeSimone:2012fs . There are
seven possible multiplets of VLQs, depending on how they couple to the SM
quarks Aguilar-Saavedra:2013qpa . These are EW singlets $[(T),(B)]$, doublets
$[(X,T),(T,B),(B,Y)$] and triplets $[(X,T,B),(T,B,Y)]$, where the fields
$X,Y,T,B$ have charges $+5/3e$, $-4/3e$, $+2/3e$ and $-1/3e$, respectively.
Moreover, VLQs can explain other anomalies observed in experiments
Branco:2022gja , such as the $W$ mass anomaly He:2022zjz ; Cao:2022mif . Many
phenomenological analyses on VLQs have also been carried out Vignaroli:2012sf
; Vignaroli:2012nf ; Yang:2019uea ; Liu:2015kmo ; Nutter:2012an ;
Deandrea:2021vje ; Buckley:2020wzk ; Gong:2019zws ; Wang:2020ips ;
Yang:2022wfa ; Han:2022rxy ; Han:2021kcr ; Shang:2022tkr ;
Cacciapaglia:2018qep ; Carvalho:2018jkq ; CarcamoHernandez:2021yev ;
Deandrea:2017rqp ; Cacciapaglia:2018lld ; Cacciapaglia:2017gzh ; King:2020mau
; Han:2024xyv ; Banerjee:2023upj ; Hernandez:2021tii ; Bhardwaj:2022wfz ;
Bardhan:2022sif ; Bhardwaj:2022nko . While most such analyses where concerned
with pair production of VLQs, as the limits on their masses grow stronger,
single production of VLQs has also become of phenomenological interest
Vignaroli:2012sf ; Vignaroli:2012nf ; Liu:2024hvp ; Canbay:2023vmj ;
Han:2023jzm ; Han:2022rxy ; Han:2022jcp ; Han:2022exz ; Yang:2022wfa ;
Han:2022iqh ; Han:2021lpg ; Tian:2021nmj ; Yang:2021btv ; Deandrea:2021vje ;
Cetinkaya:2020yjf ; Zhou:2020byj ; Cacciapaglia:2018qep ; Carvalho:2018jkq ;
Yang:2019uea ; Liu:2016jho ; Chen:2016yfv ; Aguilar-Saavedra:2013qpa ;
Barcelo:2011wu ; Shang:2024wwy ; Liu:2019jgp . In particular, in Ref.
Han:2022iqh , the authors have studied single production of Vector-Like $B$
quarks (VLBs) via the $bZ$ channel with $Z\to l^{+}l^{-}$ at the LHC. Herein,
we will complement that study as we will be dealing with the same production
process via the $bZ$ channel but with $Z\to\nu\bar{\nu}$, again at the LHC but
also in the presence of the up-to-date constraint $m_{B}>1500\text{ GeV}$.

In experiment, the ATLAS and CMS collaborations have implemented many analyses
in the search for VLQs at the LHC. For VLBs, which this work focuses on, the
CMS collaboration has performed a search via pair production in the fully
hadronic final state using Run 2 data with 138 fb-1 and set lower limits on
the VLB mass equal to 1570 GeV for BR$(B\to bH)$ = 100% and 1540 GeV for
BR$(B\to bZ)$ = 100% CMS:2020ttz , i.e., in terms of VLB Branching Ratios
(BRs). For the singlet case, BR$(B\to bH)=$BR$(B\to bZ)$= 25% and BR$(B\to
tW)=50\%$, the lower limit on the VLB mass is set to about 1060 GeV. For the
$(B,Y)$ doublet case BR$(B\to bH)=$BR$(B\to bZ)$= 50% and BR$(B\to tW)=0\%$,
the lower limit on the VLB mass is set to 1500 GeV. As intimated, although
VLBs can be pair produced through strong interactions, the single production
process through EW interactions becomes more important when VLBs are heavy,
because of phase space effects. However, we also notice that the single VLQ
production is model dependent, unlike the double VLQ one. The ATLAS
collaboration has performed a search for VLB single production decaying to a
$Z$ boson and a Higgs boson using Run 2 data with 139 fb-1 ATLAS:2023qqf . For
an isospin singlet VLB, the search has excluded values of the $c_{W}$
coupling666Note that, following Ref. ATLAS:2023qqf , $c_{W}$, $c_{Z}$, $c_{H}$
represents the coupling constants for interactions between the VLB and $W$,
$Z$ and $H$ bosons, respectively. For singlet $B$ and $m_{B}>1\text{ TeV}$,
$c_{W}\approx\kappa$, $c_{Z}\approx\kappa/(\sqrt{2}\cos\theta_{W})$, where
$\theta_{W}$ is the Weinberg angle, and $c_{H}=m_{B}\kappa/(\sqrt{2}v)$, where
$v$ is the Higgs vacuum expectation value, defined as $v\equiv 2m_{W}/g_{W}$,
with $g_{W}$ representing the $SU(2)_{L}$ weak coupling. Here, $\kappa$
corresponds to $\kappa_{B}$ in Ref. Buchkremer:2013bha , where $\zeta_{i}=1$
and $c_{W}$ and $c_{H}$ do not contain the factor $g_{W}/2$. Finally, in our
paper, $\kappa$ corresponds to $g^{*}$ when $R_{L}=0$. greater than 0.45 for
$m_{B}$ between 1.0 and 1.2 TeV. For $1.2\text{ TeV}<m_{B}<2.0\text{ TeV}$,
$c_{W}$ values in the range $0.5-0.65$ have been excluded. For a VLB as part
of a $(B,Y)$ doublet, the smallest excluded $c_{Z}$ coupling values range
between 0.3 and 0.5 for $1.2\text{ TeV}<m_{B}<2.0\text{ TeV}$.

The NP processes can be rare and only differ slightly from the SM yield, a
fact which represents a great challenge for physicists when analyzing data.
Furthermore, data collected by LHC experiments are complex and highly
dimensional whereas traditional data analysis techniques use a sequence of
boolean decisions followed by statistical analysis on the selected data, both
of which are based on the distribution of a single observed quantity motivated
by physics considerations, which is not easily extended to the required higher
dimensions Guest:2018yhq . As a result, physicists have done a lot of work to
explore the application of Machine Learning (ML), including Deep Learning
(DL), in particle physics and have developed a series of useful tools
Kasieczka:2019dbj ; Kagan:2020yrm ; Chiang:2022lsn . From these studies, we
can see that ML has a natural advantage in identifying signals, wherein the
Extreme Gradient Boosting (XGBoost) Chen_2016  algorithm, stemming from
Boosted Decision Trees (BDTs), has a unique advantage in finding NP effects
Alvestad:2021sje . For these reasons, we will adopt a XGBoost algorithm to
improve the discovery potential of VLBs at the LHC, targeting not only Run 3
datasets, i.e., with 300 fb-1 of luminosity, but also those of a future High-
Luminosity LHC (HL-LHC) Gianotti:2002xx , i.e., with 3000 fb-1 of luminosity.

We organize our paper as follows. In Sec. 2, we briefly describe the effective
Lagrangian of the singlet VLB. In Sec. 3, we present the event generation and
analysis of signal and backgrounds. In Sec. 4, we explore the exclusion and
discovery potential at the (HL-) LHC. Finally, we summarize our results in
Sec. 5.

##  II The effective Lagrangian

In a simplified model including an $SU(2)$ singlet VLB, the Lagrangian that
describes how the VLB couples to SM quarks and EW bosons can be expressed as
follows Buchkremer:2013bha 777In Ref. Buchkremer:2013bha , the effective
Lagrangian for a singlet VLB is not clearly specified, yet, it can be obtained
by incorporating the coefficients for the singlet VLB from Eq. (3.2) in Ref.
Buchkremer:2013bha . This Lagrangian is presented in the singlet VLB model of
FeynRules Alloul:2013bka . However, there is a missing factor of $\sqrt{2}$ in
the Lagrangian and the charge of $W$ is also inaccurately described on the
website feynrules . Nevertheless, the provided codes are accurate, and we
employed these in our work.:

$\begin{split}\mathcal{L}_{B}=&\frac{gg^{*}}{\sqrt{2}}\left\\{\sqrt{\frac{R_{L}}{1+R_{L}}}\frac{1}{\sqrt{2}}\left[\bar{B}_{L}W_{\mu}^{-}\gamma^{\mu}u_{L}\right]+\sqrt{\frac{1}{1+R_{L}}}\frac{1}{\sqrt{2}}\left[\bar{B}_{L}W_{\mu}^{-}\gamma^{\mu}t_{L}\right]\right.\\\
&\left.+\sqrt{\frac{R_{L}}{1+R_{L}}}\frac{1}{2\cos\theta_{W}}\left[\bar{B}_{L}Z_{\mu}\gamma^{\mu}d_{L}\right]+\sqrt{\frac{1}{1+R_{L}}}\frac{1}{2\cos\theta_{W}}\left[\bar{B}_{L}Z_{\mu}\gamma^{\mu}b_{L}\right]\right.\\\
&\left.-\sqrt{\frac{R_{L}}{1+R_{L}}}\frac{m_{B}}{2m_{W}}\left[\bar{B}_{R}Hd_{L}\right]-\sqrt{\frac{1}{1+R_{L}}}\frac{m_{B}}{2m_{W}}\left[\bar{B}_{R}Hb_{L}\right]\right\\}+h.c.\end{split}$
(1)

Here, $g$ is the $SU(2)_{L}$ gauge coupling constant, $g^{\ast}$ parameterizes
the single VLB production in association with SM quarks whereaa $m_{W}$ is the
mass of the $W$ boson. Further, $R_{L}$ is the rate of the VLB decaying into
the first and the third generation quarks. The decay widths of the VLB into SM
quarks can be estimated as:

$\displaystyle\Gamma(B\to Zd_{i})$
$\displaystyle=\frac{g^{2}g^{*2}c^{2}_{d_{i}}}{512m^{2}_{Z}m^{3}_{B}\cos^{2}\theta_{W}}\left(m^{4}_{B}+m^{4}_{d_{i}}+m^{2}_{B}m^{2}_{Z}-2m^{4}_{Z}+m^{2}_{d_{i}}(-2m^{2}_{B}+m^{2}_{Z})\right)$
(2)
$\displaystyle\times\sqrt{m^{4}_{d_{i}}+\left(m^{2}_{B}-m^{2}_{Z}\right)^{2}-2m^{2}_{d_{i}}(m^{2}_{B}+m^{2}_{Z})},$
$\displaystyle\Gamma(B\to W^{-}u_{i})$
$\displaystyle=\frac{g^{2}g^{*2}c^{2}_{u_{i}}}{256m^{2}_{W}m^{3}_{B}}\left(m^{4}_{B}+m^{4}_{u_{i}}+m^{2}_{u_{i}}m^{2}_{W}-2m^{4}_{W}+m^{2}_{B}(-2m^{2}_{u_{i}}+m^{2}_{W})\right)$
$\displaystyle\times\sqrt{m^{4}_{B}+\left(m^{2}_{u_{i}}-m^{2}_{W}\right)^{2}-2m^{2}_{B}(m^{2}_{u_{i}}+m^{2}_{W})},$
$\displaystyle\Gamma(B\to Hd_{i})$
$\displaystyle=\frac{g^{2}g^{*2}c^{2}_{d_{i}}}{512m_{W}^{2}m_{B}}\left(m^{2}_{B}+m^{2}_{d_{i}}-m^{2}_{H}\right)\sqrt{m^{4}_{d_{i}}+\left(m^{2}_{B}-m^{2}_{H}\right)^{2}-2m^{2}_{d_{i}}\left(m^{2}_{B}+m^{2}_{H}\right).}$

Here, $i=1,3$ denotes the generation of SM quarks,
$c^{2}_{d_{i},u_{i}}=R_{L}/(1+R_{L})$ for $i=1$,
$c^{2}_{d_{i},u_{i}}=1/(1+R_{L})$ for $i=3$. $m_{u_{i}}$ and $m_{d_{i}}$ are
masses of the up-type and down-type SM quarks, respectively. Furthermore,
$m_{H}$ is the mass of the SM-like Higgs boson. If $R_{L}=0$, the VLB is
assumed to couple only to the third generation of SM quarks Aguilar-
Saavedra:2013qpa  and the Lagrangian can be simplified as:

$\begin{split}\mathcal{L}_{B}=&\frac{gg^{*}}{\sqrt{2}}\left\\{\frac{1}{\sqrt{2}}\left[\bar{B}W_{\mu}^{-}\gamma^{\mu}t_{L}\right]+\frac{1}{2\cos\theta_{W}}\left[\bar{B}Z_{\mu}^{-}\gamma^{\mu}b_{L}\right]\right.\\\
&\left.-\frac{m_{B}}{2m_{W}}\left[\bar{B}_{R}Hb_{L}\right]-\frac{m_{B}}{2m_{W}}\left[\bar{B}_{L}Hb_{R}\right]\right\\}+h.c.\end{split}$
(3)

For $m_{B}>1\text{ TeV}$, the decay BR of the VLB tends to be such that
BR$(B\to tW)$:BR$(B\to bZ)$:BR$(B\to bH)$ = 2:1:1, as expected from the
Goldstone equivalence theorem He:1992nga ; He:1993yd .

Figure 1: The BRs of the VLB as a function of $R_{L}$ with $m_{B}$=1500 GeV
(left) and the width-to-mass ratios of the VLB as a function of $m_{B}$
(right).

In Fig. 1, we show the BRs of our VLB as a function of $R_{L}$ (left) for
$m_{B}$=1500 GeV and the width-to-mass ratios as a function of $m_{B}$
(right). The BRs are independent of the coupling coefficient $g^{\ast}$
whereas the width-to-mass ratios are approximately independent of $R_{L}$ when
$m_{B}>1\text{ TeV}$. For a fixed $g^{\ast}$, the width-to-mass ratios are
proportional to $m_{B}^{2}$. We can see that the VLB mainly decays to the
third generation quarks for $R_{L}<0.5$. However, when $R_{L}>1$, the VLB
decay to the first generation quarks gradually becomes the dominant decay
mode. In this work, we primarily focus on the VLB coupling to the third
generation quarks, hence, we consider two cases: $R_{L}=0$ and $R_{L}=0.5$.
Since the constraints from the EW Precision Observables (EWPOs) on the
coupling coefficient $g^{\ast}$ are of $\mathcal{O}(10^{-1})$ Cao:2022mif , we
consider here the range $g^{\ast}\leq 0.5$.

Figure 2: Feynman diagrams for single VLB (in red) production with subsequent
decay into $bZ$ at the LHC.

We explore the discovery potential of single VLB production through the
channel $pp\to B(\to bZ)j\to b(Z\to\nu_{l}\bar{\nu_{l}})j$, where $b$ denotes
the $b$-tagged jet and $j$ stands for light-flavor jets. We show the typical
Leading Order (LO) Feynman diagrams in Fig. 2. From this figure, we can see
that the signal events are characterized by two jets with at least one
$b$-tagged jet among these, along with significant missing (transverse)
energy. Therefore, the primary SM backgrounds are $pp\to Zjjj$, $pp\to bZj$,
$pp\to ZZ$ and $pp\to ZH$. Among these, $pp\to bZj$ is an irreducible
background while the others are considered reducible backgrounds. The latter
must be taken into account due to jets exhibiting lower tagging efficiencies
and higher mistagging rates compared to leptons at hadron colliders. We have
also evaluated further backgrounds including $pp\to W^{\pm}Z$ and determined
that their impact can be disregarded according to the selection criteria that
will be discussed later.

Figure 3: Cross sections of the signal process $pp\to B(\to bZ)j\to
b(Z\to\nu_{l}\bar{\nu_{l}})j$ as a function of $m_{B}$ at the 14 TeV LHC with
$g^{\ast}$ = 0.1, 0.2, 0.5 for $R_{L}$=0 (left) and $R_{L}$=0.5 (right). Here,
the conjugate processes have been included.

We calculate the LO cross sections of the signal process by MadGraph5_aMC@NLO
Alwall:2014hca  with the NNPDF23NLO Parton Distribution Functions (PDFs) NNPDF
and default renormalization and factorization scales. The cross sections at 14
TeV LHC are depicted as a function of $m_{B}$ in Fig. 3, where we take
$g^{\ast}$ = 0.1, 0.2 and 0.5 as well as $R_{L}=0$ and $0.5$ as examples. The
input SM parameters are taken as follows888Here, $m_{W}$ is calculated as
$\sqrt{\frac{M_{Z}^{2}}{2}+\sqrt{\frac{M_{Z}^{4}}{4}-\frac{\alpha_{\text{EW}}\pi
M_{Z}^{2}}{G_{F}\sqrt{2}}}}$ within MadGraph5_aMC@NLO.
ParticleDataGroup:2022pth :

$m_{b}=4.18\text{ GeV}$, $m_{Z}=91.1876\text{ GeV}$, $m_{H}=125.25\text{
GeV}$,

$G_{F}=1.166378\times 10^{-5}\text{ GeV}^{-2}$, $\sin^{2}\theta_{W}=0.231$,
$\alpha_{\text{EM}}(m_{Z})=1/128$.

From Fig. 3, we can see that the signal cross sections decrease rapidly with
increasing $m_{B}$. For a fixed $m_{B}$, the cross sections increase rapidly
with increasing $g^{\ast}$ due to $\sigma\propto g^{\ast 2}$. Furthermore, for
equivalent values of the coupling $g^{*}$ and the mass $m_{B}$, the cross
section for $R_{L}=0.5$ is larger than that for $R_{L}=0$. This is because the
former scenario involves $B$ mixing with the first generation quarks, which
significantly enlarges the single VLB production rate due to the increased sea
quark PDFs. However, we limit ourselves to studying cases where $R_{L}<0.5$
because there are stringent constraints from non-LHC flavor physics data on
VLQs mixing significantly with lighter quarks Belfatto:2023tbv ;
Crivellin:2022rhw ; Balaji:2021lpr ; Branco:2021vhs ; Vatsyayan:2020jan ;
Botella:2016ibj ; Ishiwata:2015cga ; Alok:2015iha ; Cacciapaglia:2015ixa ;
Botella:2012ju ; Cacciapaglia:2010vn . We list the LO cross sections and the
$K$ factors representing the QCD corrections for the background processes in
Tab. 1.

Processes | $Zjjj$ | $bZj$ | $ZZ$ | $ZH$
---|---|---|---|---
$\sigma_{\text{LO}}$(fb) | 22.5 | 5.8 | 0.1 | 0.1
$K$ factor | 1.2 Alwall:2014hca  | 1.2 Campbell:2005zv  | 1.8 Campbell:1999ah  | 1.3 Alwall:2014hca
Table 1: The LO cross sections and the $K$ factors representing the QCD
corrections for the background processes at the 14 TeV LHC. Here, the cross
section of $bZj$ includes the contribution from the conjugate process of
$pp\to\bar{b}Zj$.

##  III Event generation and analysis

We conduct a detailed detector simulation in this work. The signal model file
is obtained from FeynRules feynruls  and, again, parton-level events are
generated using MadGraph5_aMC$@$NLO with the NNPDF23NLO PDF set. Indeed, the
factorization and renormalization scales are still set to their default values
in MadGraph5_aMC$@$NLO.

These events are then input into Pythia 8.3 Bierlich:2022pfr  for parton
showering and hadronization. Subsequently, fast detector simulations are
performed using Delphes 3.4.2 deFavereau:2013fsa  with the built-in detector
configurations of the LHC Run 3 and HL-LHC website_hllhc . Jets are clustered
using FastJet Cacciari:2011ma  employing the anti-$kt$ algorithm
Cacciari:2005hq . Settings for the cone radius of the jet tagging algorithm
and isolated leptons are implemented in Delphes with default values. Then,
both signal and background events are analyzed using MadAnalysis 5
Conte:2012fm . Finally, we scan the VLB parameter space and connect these
programs with the help of the EasyScan_HEP package Shang:2023gfy . To quantify
the observability of the signal, we utilize the following formula to calculate
the expected discovery and exclusion significance, without considering
systematic uncertainties for simplicity Burns:2011xf :

$\mathcal{Z}_{\rm
disc}=\sqrt{2\left[\left(s+b\right)\ln\left(1+\frac{s}{b}\right)-s\right]},$
(4)

and

$\mathcal{Z}_{\rm
excl}=\sqrt{2\left[s-b\ln\left(1+\frac{s}{b}\right)\right]},$ (5)

where $s$ and $b$ are the number of signal ($s$) and total background ($b$)
events, respectively, after applying the selection criteria that will be
discussed later. The exclusion significance is denoted by $\mathcal{Z}_{\rm
excl}$ = 2 whereas the discovery significance is denoted by $\mathcal{Z}_{\rm
disc}$ = 5.

To accommodate contributions from higher-order QCD corrections, we adjust the
cross sections of the dominant backgrounds from the LO to the Next-to-LO (NLO)
or Next-to-NLO (NNLO) using $K$ factors, as detailed in Tab. 1. We assume that
kinematic distributions are minimally impacted by these higher-order effects.
Therefore, for simplicity, we scale the kinematic distributions, which we will
discuss later, using constant bin-independent $K$ factors. It is then worth
noting that, for the signal, we maintain the LO cross section, ensuring that
the exclusion and discovery potentials depicted in this analysis remain
conservative. The inclusion of higher-order QCD corrections could indeed
bolster these potentials.

For our analysis, we establish as reference values for the couplings
$g^{*}=0.2$ and $R_{L}=0$. Furthermore, along the VLB mass direction, our
Benchmark Points (BPs) are located at the $m_{B}$ values of 1500 GeV and 2000
GeV. However, we will later showcase the excluded $2\sigma$ and $5\sigma$
discovery thresholds in the $g^{*}-m_{B}$ plane, considering two different
values of $R_{L}$, namely 0 and 0.5 and more finely sifted $m_{B}$ values.

Figure 4: Normalized distributions for the signals ($m_{B}=1500$ GeV and 2000
GeV are selected as examples) and backgrounds at the 14 TeV LHC for $g^{\ast}$
= 0.2 and $R_{L}=0$.

Due to the large mass of the VLB, its decay products are highly boosted, this
resulting in the transverse momentum ($p_{T}$) peaks of the signals being
larger than those of the (SM) backgrounds. In Fig. 4, we present several
differential distributions for signals and backgrounds at the LHC, including
the transverse momentum spectrum of the leading $b$-tagged jet
($p_{T}^{b_{1}}$), the distribution of total transverse hadronic energy
$H_{T}$ (defined as $H_{T}\equiv\sum_{\text{hadronic
particles}}|\vec{p}_{T}|$) and the that of the missing transverse energy
$\not{E}_{T}$. Additionally, the signal process involves an outgoing light
quark, leading to a light-flavor forward jet within the detector. Therefore,
we also display kinematic distributions showcasing distinct spatial
characteristics, such as the pseudorapidity distribution of the leading light-
flavor jet $\eta_{j_{1}}$ and the distribution of spatial separation $\Delta
R$ between the leading $b$-tagged jet and the leading light-flavor jet. The
spatial separation is defined as $\Delta
R=\sqrt{\Delta\Phi^{2}+\Delta\eta^{2}}$, representing the separation in the
rapidity ($\eta$)–azimuth ($\phi$) plane.

Then, we employ the XGBoost algorithm to enhance the discovery potential of
the VLB at the LHC. We select four signal BPs with $m_{B}$ values of 1500 GeV,
2000 GeV, 3000 GeV and 4000 GeV while fixing $g^{*}=0.1$ and $R_{L}=0$
throughout, as kinematic distributions are minimally affected by the coupling
$g^{*}$ and the mixing parameter $R_{L}$. The background events are generated
for the following processes here:

  * •

$pp\to Z(\to\nu_{l}\bar{\nu_{l}})jjj$,

  * •

$pp\to bZ(\to\nu_{l}\bar{\nu_{l}})j$,

  * •

$pp\to Z(\to\nu_{l}\bar{\nu_{l}})Z(\to b\bar{b}/jj)$,

  * •

$pp\to Z(\to\nu_{l}\bar{\nu_{l}})H(\to b\bar{b})$.

Considering the characteristics of the signal distributions, for the purpose
of optimising the Monte Carlo (MC) generation efficiency, events in this
analysis must satisfy $p_{T}^{j,b}>200\text{ GeV}$ and $\not{E}_{T}>400\text{
GeV}$ at the parton level999We have indeed checked that these are not biasing
the detector level events, see Cut-2 and -4 below., implemented using
MadGraph5_aMC$@$NLO, and $N(j)\geq 1$ and $N(b)\geq 1$ at the detector level,
implemented using MadAnalysis 5.

Each event from the signal BPs and the four types of backgrounds corresponds
to a 10-dimensional event vector constructed by 9-dimensional input features
and a 1-dimensional target variable. The input features consist of the number
of $b$-tagged jets $N(b)$, the number of light-flavor jets $N(j)$, the
transverse momentum ($p_{T}^{b_{1}}$) and pseudorapidity ($\eta_{b_{1}}$) of
the leading $b$-tagged jet, the transverse momentum ($p_{T}^{j_{1}}$) and
pseudorapidity ($\eta_{j_{1}}$) of the leading light-flavor jet, the spatial
separation between the leading $b$-jet and light-jet $\Delta R(b_{1},j_{1})$,
the total transverse hadronic energy $H_{T}$ and the missing energy
$\not{E}_{T}$, obtained with the assistance of expert mode in MadAnalysis 5.
The target variable is set as 1 for signal events and 0 for background events.
We employ $8\times 10^{5}$ event vectors, with equal contributions from each
signal BP and background process, as training data. This dataset enables the
XGBoost model (hereafter referred to as the X-model) to learn the
relationships between features and how to map these onto target variables.
Additionally, we use $1.6\times 10^{4}$ event vectors, with equal
contributions from each signal BP and background process, as testing data to
evaluate the performance of the trained X-model on new signal and background
events different from those in the training data.

There are hyperparameters used to fine-tune the performance of the X-model.
The learning rate $\alpha$ controls the step size during the minimization of a
logistic loss function $\mathcal{F}$:

$\mathcal{F}=-\frac{1}{N}\sum_{i=1}^{N}\left[y_{i}\ln(p_{i})+(1-y_{i})\ln(1-p_{i})\right].$
(6)

Here, $N$ represents the number of events used to train the X-model, $y_{i}$
denotes the target variable for each event and $p_{i}$ indicates the predicted
probability of the signal class for the $i$-th event. The fractions of event
vectors ($r_{\text{s}}$) and features ($r_{\text{f}}$) in each event vector
can be randomly selected for each tree. This randomness is introduced into the
X-model to enhance its generality. Additionally, the parameter $\gamma$
prevents overfitting by penalizing excessively large and deep trees. If
$\mathcal{F}>\gamma$, a further partition on a leaf node of the tree is
required. We specify the learning rate to range between 0.01 and 0.5. The
ranges for the fractions of event vectors and features are established as
$r_{\text{s,f}}\in[0.5,1.0]$. Furthermore, the ranges for the number of trees
$n_{\text{t}}$ and the maximum depth of each tree $n_{\text{d}}$ are set as
[50, 201] and [3, 21], respectively. A summary of these hyperparameter ranges
and the optimal set of hyperparameters is provided in Tab. 2.

Hyperparameter | $\alpha$ | $\gamma$ | $n_{\text{t}}$ | $n_{\text{d}}$ | $r_{\text{s}}$ | $r_{\text{f}}$
---|---|---|---|---|---|---
Range | [0.01, 0.3] | [1, 5] | [20, 201] | [3, 21] | [0.5, 1] | [0.5, 1]
Optimal value | 0.2 | 1.3 | 183 | 15 | 0.96 | 0.99
Table 2: The hyperparameter ranges and optimal values considered in this work.

When assessing the performance of a set of hyperparameters, each $p_{i}$ and
$p_{i+1}$ serves as threshold probabilities denoted by $\mathcal{T}_{j}$
($j=1,\dots,N+1$). For every $\mathcal{T}_{j}$, if $p_{i}>\mathcal{T}_{j}$
($i=1,\dots,N$), the corresponding event is classified as a signal, otherwise,
it is classified as background. Consequently, each set of hyperparameters is
associated with a Receiver Operating Characteristic (ROC) curve, which
illustrates the True Positive Rate (TPR) against the False Positive Rate (FPR)
at various threshold settings. Here, $j$ ranges from 1 to $N+1$ and the Area
Under the ROC Curve (AUC) FAWCETT2006861  determines the effectiveness of the
hyperparameters. Bayesian optimization, as introduced by snoek2012practical ,
is employed to discover the optimal set of hyperparameters for this study. The
objective function chosen for this optimization is the AUC value. During each
iteration of Bayesian optimization, the training event vectors are randomly
split into a 4:1 ratio. The former portion is utilized to train the X-model
using XGBoost, while the latter is reserved for guiding the selection of the
next set of hyperparameters to evaluate via Bayesian optimization.
Subsequently, the optimal set of hyperparameters is identified based on
maximizing the AUC value. Once this optimal set is determined, it represents
the desired configuration for the X-model, which we refer to as the trained
X-model. The trained X-model is then applied to the testing data as discussed
earlier. Subsequently, the threshold probability closest to the point (1,0) on
the resulting ROC curve is used as the selection criterion of the XGBoost
algorithm to distinguish signal from background events used to achieve
discovery or exclusion. We denote this threshold probability as
$\mathcal{T}^{*}$. Ultimately, we utilize the trained X-model to generate the
ROC curve and probability distributions of X-model outputs for events in the
testing data. These plots, along with the selection criterion
$\mathcal{T}^{*}$, are depicted in Fig. 5.

Figure 5: Model performance on testing data. The left is the ROC curve while
AUC represents the area under the curve. The right is the score distributions
relevant to testing data. The green line stands for the signal while the blue
line stands for all the backgrounds. The vertical red line stands for the
threshold we choose.

##  IV LHC Exclusion and Discovery Potential

Cuts | Signals | Backgrounds
---|---|---
|  1500 GeV |  2000 GeV |  |  $Zjjj$ |  $bZj$ |  $ZZ$ |  $ZH$
Initial |  0.079 |  0.024 |  |  26.94 |  6.26 |  0.14 |  0.11
Cut-1 |  0.056 |  0.017 |  |  2.26 |  4.55 |  0.017 |  0.046
Cut-2 |  0.030 |  0.011 |  |  0.24 |  0.77 |  0.0009 |  0.013
Cut-3 |  0.015 |  0.006 |  |  0.07 |  0.13 |  0.0002 |  0.0046
Cut-4 |  0.010 |  0.005 |  |  0.02 |  0.04 |  $4.26\times 10^{-6}$ |  0.0001
Efficiency |  12.9% |  20.8% |  |  $7.4\times 10^{-4}$ |  $6.8\times 10^{-3}$ |  $3.0\times 10^{-5}$ |  $9.0\times 10^{-4}$
Table 3: Cut flow of the cross sections (in fb) for the signals with two
typical VLB quark masses and backgrounds at the 14 TeV LHC using the cut-and-
count method. Here we take the parameters $g^{\ast}=0.2$ and $R_{L}=0$.

Drawing from the kinematic distributions illustrated in Fig. 4, we implement
the subsequent kinematic selection criteria on the events to distinguish the
signal from the backgrounds:

  1. Cut-1:

There are at least two jets ($N(j,b)\geq 2$) and at least one $b$-tagged jet
among these ($N(b)\geq 1$).

  2. Cut-2:

The transverse momentum of the leading $b$-tagged jet is required to be
greater than $500\text{ GeV}$ ($p_{T}^{b_{1}}>500\text{ GeV}$).

  3. Cut-3:

The pseudorapidity of the leading light-flavor jet is required to be greater
than $1.8$ ($|\eta_{j_{1}}|>1.8$). Additionally, the spatial separation
between the leading $b$-tagged jet and the leading light-flavor jet is
required to be greater than $2.0$ to reduce background from the process $pp\to
bZj$, which contains the $b$ quark in the initial state ($\Delta
R(j_{1},b_{1})>2.0$).

  4. Cut-4:

The total transverse hadronic energy is required to be greater than $800\text{
GeV}$ ($H_{T}>800\text{ GeV}$) and the transverse missing energy is required
to be greater than $500\text{ GeV}$ ($\not{E}_{T}>500\text{ GeV}$).

Cuts | Signals | Backgrounds
---|---|---
|  1500 GeV |  2000 GeV |  |  $Zjjj$ |  $bZj$ |  $ZZ$ |  $ZH$
Initial |  0.079 |  0.024 |  |  26.94 |  6.26 |  0.14 |  0.11
Basic cuts |  0.056 |  0.017 |  |  2.26 |  4.55 |  0.017 |  0.046
$\mathcal{T}^{*}$ |  0.006 |  0.005 |  |  0.0075 |  0.019 |  0.000005 |  0.00002
Efficiency |  7.5% |  20.8% |  |  $2.8\times 10^{-4}$ |  $3.0\times 10^{-3}$ |  $3.6\times 10^{-5}$ |  $1.8\times 10^{-4}$
Table 4: Cut flow of the cross sections (in fb) for the signals with two
typical VLB quark masses and backgrounds at the 14 TeV LHC using the XGBoost
method. Here we take the parameters $g^{\ast}=0.2$ and $R_{L}=0$.

We present the cross sections of two representative signals, with $m_{B}$
values set at 1500 GeV and 2000 GeV, along with the corresponding backgrounds.
These are depicted after applying the selection criteria detailed in Tab. 3
and Tab. 4 for the cut-and-count method and XGBoost method, respectively. To
construct these tables, we have simulated $10^{5}$ events at the detector
level for each signal and background process.

From Tab. 3, we can see that, among the four types of backgrounds, the
dominant ones are the $Zjjj$ and $bZj$ prior to any cuts. The first two cuts,
pertaining to the number of final jets and transverse momentum of the leading
$b$-tagged jets, significantly cut events from $Zjjj$ and $ZZ$, achieving a
cumulative reduction rate of approximately 99%. Subsequently, the third cut,
addressing pseudorapidities and spatial separations, efficiently cuts $bZj$
events with a reduction rate of about 83%. Furthermore, imposing stringent
requirements on $H_{T}$ and $\not{E}_{T}$ successfully cuts overall around 70%
of background events while retaining 70% of signal ones. Altogether, at the
end of the cut flow, all backgrounds undergo notable suppression whereas the
signals maintain a relatively favorable efficiency. In the end, the primary
background contributions originate from the $Zjjj$ and $bZj$ SM processes,
which remain comparable in magnitude to the signal remnants, exhibiting cross
sections of $0.02\text{ fb}$ and $0.04\text{ fb}$, respectively. Finally,
compared to the cut-and-count method, from Tab. 4, we can see that the signal
significance is clearly improved further by the XBoost method.

To determine the exclusion and discovery reach at the 14 TeV LHC, using the
two methods adopted here (i.e., the traditional cut-and-count one and the ML
one), we perform calculations spanning the interval $m_{B}=1500\text{ GeV}$ to
3000 GeV. BPs are now selected at intervals of 250 GeV in the VLB mass, with
$R_{L}=0$ and $g^{*}=0.2$. For each signal BP, we simulate $10^{5}$ events,
while for each background, we simulate $10^{6}$ events. It has been verified
that different values of $R_{L}$ and $g^{*}$ yield the same efficiency when
$m_{B}$ remains constant.

In Fig. 6, we show the $2\sigma$ (exclusion) and $5\sigma$ (discovery) lines
on the $g^{\ast}-m_{B}$ plane for $R_{L}=0$ at the 14 TeV LHC. We can see that
VLBs can be excluded in the regions of $g^{\ast}\in\left[0.35,0.50\right]$ and
$m_{B}\in\left[1500\text{ GeV},2300\text{ GeV}\right]$ assuming an integrated
luminosity of 300 fb-1 (Run 3). With a value for the latter of 3000 fb-1 (HL-
LHC), instead, the exclusion regions can be extended to
$g^{\ast}\in\left[0.18,0.50\right]$ and $m_{B}\in\left[1500\text{
GeV},3000\text{ GeV}\right]$. Accordingly, the discovery regions are
$g^{\ast}\in\left[0.30,0.50\right]$ and $m_{B}\in\left[1500\text{
GeV},2500\text{ GeV}\right]$.

Fig. 7 is the same as Fig. 6, but for $R_{L}=0.5$. We can see that the
exclusion regions are $g^{\ast}\in\left[0.10,0.50\right]$ and
$m_{B}\in\left[1500\text{ GeV},4050\text{ GeV}\right]$ with 300 fb-1, which
can be extended to $g^{\ast}\in\left[0.05,0.50\right]$ and
$m_{B}\in\left[1500\text{ GeV},4750\text{ GeV}\right]$ with 3000 fb-1. The
discovery regions are $g^{\ast}\in\left[0.16,0.50\right]$ and
$m_{B}\in\left[1500\text{ GeV},3500\text{ GeV}\right]$ with 300 fb-1 plus
$g^{\ast}\in\left[0.08,0.50\right]$ and $m_{B}\in\left[1500\text{
GeV},4250\text{ GeV}\right]$ with 3000 fb-1.

Figure 6: Exclusion (left) and discovery (right) potential for the VLB signal
on the $g^{\ast}-m_{B}$ planes for $R_{L}=0$ at the 14 TeV LHC, considering
the two methods adopted in our analysis. Three typical luminosities, 300 fb-1,
1000 fb-1 and 3000 fb-1, are illustrated. The dashed lines denote the width-
to-mass ratios $\Gamma_{B}/m_{B}$.

Figure 7: Same as Fig. 6, but for $R_{L}=0.5$.

##  V Summary

In this work, we have utilized a simplified model approach to test the
exclusion and discovery potential of a VLB signal at both the current LHC (Run
3) and a future HL-LHC. We have employed both a standard cut-and-count method
and a ML based XGBoost algorithm to achieve this goal. Specifically, we have
focused on the process $pp\to B(\to bZ)j\to b(Z\to\nu_{l}\bar{\nu_{l}})j$ as
signal.

Initially, we have generated signal and background events using state-of-the-
art computing tools. Subsequently, we have applied the cut-and-count method
and the XGBoost algorithm to classify these events, effectively distinguishing
between signal and background. For the XGBoost approach, we have trained the
ML model using MC data and utilized Bayesian optimization to determine the
optimal combination of hyperparameters. We have then employed this trained
model to classify events and determine cut efficiencies via selection
thresholds. Finally, we have evaluated the exclusion and discovery
capabilities of the signal over both reducible and irreducible backgrounds
across different LHC luminosities (and fixed $\sqrt{s}=14$ TeV).

Our analysis has revealed significant sensitivities to this VLB process across
various LHC stages, spanning from current to future machine configurations.
These appear to improve significantly upon more traditional methods used in
phenomenological analyses.

LHC | Decay | Exclusion | Discovery | Reference
---|---|---|---|---
|  | $g^{\ast}$ |  $m_{B}$ (GeV) | $g^{\ast}$ |  $m_{B}$ (GeV) |
$R_{L}=0$ | $B\to tW$ | [0.2,0.5] | [1800,2400] | [0.3,0.5] | [1750,2100] | Han:2022jcp
$R_{L}=0$ | $B\to bZ$ | [0.2,0.5] | [1300,1800] | [0.3,0.5] | [1300,1500] | Han:2022iqh
$R_{L}=0$ | $B\to bZ$ | [0.2,0.5] | [1500,3000] | [0.3,0.5] | [1500,2500] | this work
$R_{L}$=0.5 | $B\to bZ$ | [0.2,0.5] | [3750,4750] | [0.3,0.5] | [3700,4250]
Table 5: Comparison of exclusion and discovery capabilities of different VLB
signals with $\sqrt{s}=14$ TeV and $\mathcal{L}=3000$ fb-1 from different
works.

Indeed, for comparison, we finish by listing previous results on searching for
a VLB at the HL-LHC in Tab. 5, also including $B\to tW$ decays, against ours.
Since the advantages of XGBoost method lie in the high-mass region, we focus
for this comparison over the region $m_{B}>$1500 GeV. (Recall that the VLB
coupling strength $g^{\ast}$ is also described as $\kappa_{B}$ in other
references and their relation is $g^{\ast}=\kappa_{B}$ Buchkremer:2013bha .)
By comparing these results, we can see that the XGBoost algorithm does improve
significantly the sensitivity of searches for a VLB with the same luminosity,
a result which makes the exploration of VLB signals ever more promising also
from the experimental side. However, it should be noted that, due to the
inherent characteristics of the ML algorithm, certain limits to its
applicability exist. For instance, when the mass of the VLB greatly exceeds
the upper bound of the training set (in this work, 4000 GeV), the predictive
capability of the model significantly diminishes. In forthcoming research, our
primary focus will be on enhancing the model to address this issue.

## Acknowledgements

We thank Luca Panizzi for helpful discussions. This work of BY, ZL, XJ and LS
is supported by the Natural Science Foundation of Henan Province under Grant
No. 232300421217, the National Research Project Cultivation Foundation of
Henan Normal University under Grant No. 2021PL10, the China Scholarship
Council under Grant No. 202208410277 and also powered by the High Performance
Computing Center of Henan Normal University. The work of SM is supported in
part through the NExT Institute, the Knut and Alice Wallenberg Foundation
under the Grant No. KAW 2017.0100 (SHIFT) and the STFC Consolidated Grant No.
ST/L000296/1.

## References

  * (1) K. G. Wilson, “The Renormalization Group and Strong Interactions,” Phys. Rev. D 3 (1971) 1818.
  * (2) ATLAS Collaboration, G. Aad et al., “Combined measurements of Higgs boson production and decay using up to $80$ fb-1 of proton-proton collision data at $\sqrt{s}=$ 13 TeV collected with the ATLAS experiment,” Phys. Rev. D 101 no. 1, (2020) 012002, arXiv:1909.02845 [hep-ex].
  * (3) CMS Collaboration, A. M. Sirunyan et al., “Measurement and interpretation of differential cross sections for Higgs boson production at $\sqrt{s}=$ 13 TeV,” Phys. Lett. B 792 (2019) 369–396, arXiv:1812.06504 [hep-ex].
  * (4) A. De Simone, O. Matsedonskyi, R. Rattazzi, and A. Wulzer, “A First Top Partner Hunter’s Guide,” JHEP 04 (2013) 004, arXiv:1211.5663 [hep-ph].
  * (5) J. A. Aguilar-Saavedra, R. Benbrik, S. Heinemeyer, and M. Pérez-Victoria, “Handbook of vectorlike quarks: Mixing and single production,” Phys. Rev. D 88 no. 9, (2013) 094010, arXiv:1306.0572 [hep-ph].
  * (6) G. C. Branco and M. N. Rebelo, “Vector-like Quarks,” PoS DISCRETE2020-2021 (2022) 004, arXiv:2208.07235 [hep-ph].
  * (7) S. P. He, “Leptoquark and vector-like quark extended model for simultaneous explanation of $W$ boson mass and muon g–2 anomalies*,” Chin. Phys. C 47 no. 4, (2023) 043102, arXiv:2205.02088 [hep-ph].
  * (8) J. J. Cao, M. L., L. L. Shang, S. Y. Wang, and B. F. Yang, “Interpreting the $W$-mass anomaly in vectorlike quark models,” Phys. Rev. D 106 no. 5, (2022) 055042, arXiv:2204.09477 [hep-ph].
  * (9) N. Vignaroli, “Discovering the composite Higgs through the decay of a heavy fermion,” JHEP 07 (2012) 158, arXiv:1204.0468 [hep-ph].
  * (10) N. Vignaroli, “Early discovery of top partners and test of the Higgs nature,” Phys. Rev. D 86 (2012) 075017, arXiv:1207.0830 [hep-ph].
  * (11) B. F. Yang, B. F. Hou, H. Y. Zhang, and N. Liu, “Single production of vectorlike top partner decaying to $Wb$ in the leptonic channel at $ep$ colliders in the littlest Higgs model with $T$-parity,” Phys. Rev. D 99 no. 9, (2019) 095002, arXiv:1904.07434 [hep-ph].
  * (12) N. Liu, L. Wu, B. F. Yang, and M. C. Zhang, “Single top partner production in the Higgs to diphoton channel in the Littlest Higgs Model with $T$-parity,” Phys. Lett. B 753 (2016) 664–669, arXiv:1508.07116 [hep-ph].
  * (13) J. Nutter, R. Schwienhorst, D. G. E. Walker, and J. Yu, “Single Top Production as a Probe of $B^{\prime}$ Quarks,” Phys. Rev. D 86 (2012) 094006, arXiv:1207.5179 [hep-ph].
  * (14) A. Deandrea, T. Flacke, B. Fuks, L. Panizzi, and H. Shao, “Single production of vector-like quarks: the effects of large width, interference and NLO corrections,” JHEP 08 (2021) 107, arXiv:2105.08745 [hep-ph]. [Erratum: JHEP 11, 028 (2022)].
  * (15) A. Buckley, J. M. Butterworth, L. Corpe, D. Huang, and P. Sun, “New sensitivity of current LHC measurements to vector-like quarks,” SciPost Phys. 9 no. 5, (2020) 069, arXiv:2006.07172 [hep-ph].
  * (16) X. Gong, C. X. Yue, and Y. C. Guo, “Search for vector-like bottom quark via $Zb$ production at the LHC,” Phys. Lett. B 793 (2019) 175–180, arXiv:1904.08206 [hep-ph].
  * (17) D. H. Wang, L. Wu, and M. C. Zhang, “Hunting for top partner with a new signature at the LHC,” Phys. Rev. D 103 no. 11, (2021) 115017, arXiv:2007.09722 [hep-ph].
  * (18) B. F. Yang, X. L. Sima, S. Y. Wang, and L. L. Shang, “Single vectorlike top quark production in the $tZ$ channel at high energy pp colliders,” Phys. Rev. D 105 no. 9, (2022) 096010.
  * (19) J. Z. Han, S. Xu, H. Q. Song, and Y. J. Wang, “Single production of vectorlike $B$ quarks decaying into $bh$ at the CLIC,” Nucl. Phys. B 985 (2022) 116030.
  * (20) L. Han and J. F. Shen, “Search for single production of vector-like $B$ quark decaying to $bZ$ at future linear colliders,” Eur. Phys. J. C 81 no. 5, (2021) 463.
  * (21) L. L. Shang, C. P. Chen, S. Y. Wang, and B. F. Yang, “Single production of vector-like $B$ quark decaying into $bZ$ at future $ep$ colliders,” Nucl. Phys. B 984 (2022) 115977.
  * (22) G. Cacciapaglia, A. Carvalho, A. Deandrea, T. Flacke, B. Fuks, D. Majumder, L. Panizzi, and H.-S. Shao, “Next-to-leading-order predictions for single vector-like quark production at the LHC,” Phys. Lett. B 793 (2019) 206–211, arXiv:1811.05055 [hep-ph].
  * (23) A. Carvalho, S. Moretti, D. O’Brien, L. Panizzi, and H. Prager, “Single production of vectorlike quarks with large width at the Large Hadron Collider,” Phys. Rev. D 98 no. 1, (2018) 015029, arXiv:1805.06402 [hep-ph].
  * (24) A. E. Cárcamo Hernández, S. F. King, and H. Lee, “$Z$ mediated flavor changing neutral currents with a fourth vectorlike family,” Phys. Rev. D 105 no. 1, (2022) 015021, arXiv:2110.07630 [hep-ph].
  * (25) A. Deandrea and A. M. Iyer, “Vectorlike quarks and heavy colored bosons at the LHC,” Phys. Rev. D 97 no. 5, (2018) 055002, arXiv:1710.01515 [hep-ph].
  * (26) G. Cacciapaglia, A. Deandrea, N. Gaur, D. Harada, Y. Okada, and L. Panizzi, “The LHC potential of Vector-like quark doublets,” JHEP 11 (2018) 055, arXiv:1806.01024 [hep-ph].
  * (27) G. Cacciapaglia, H. Cai, A. Carvalho, A. Deandrea, T. Flacke, B. Fuks, D. Majumder, and H.-S. Shao, “Probing vector-like quark models with Higgs-boson pair production,” JHEP 07 (2017) 005, arXiv:1703.10614 [hep-ph].
  * (28) S. J. D. King, S. F. King, S. Moretti, and S. J. Rowley, “Discovering the origin of Yukawa couplings at the LHC with a singlet Higgs and vector-like quarks,” JHEP 21 (2020) 144, arXiv:2102.06091 [hep-ph].
  * (29) J.-Z. Han, Y.-B. Liu, and S.-Y. Xu, “Pair production of the singlet vector-like B quark at the CLIC,” Eur. Phys. J. C 84 no. 1, (2024) 61, arXiv:2401.11423 [hep-ph].
  * (30) A. Banerjee, V. Ellajosyula, and L. Panizzi, “Heavy vector-like quarks decaying to exotic scalars: a case study with triplets,” JHEP 01 (2024) 187, arXiv:2311.17877 [hep-ph].
  * (31) A. E. C. Hernández, S. F. King, and H. Lee, “Fermion mass hierarchies from vectorlike families with an extended 2HDM and a possible explanation for the electron and muon anomalous magnetic moments,” Phys. Rev. D 103 no. 11, (2021) 115024, arXiv:2101.05819 [hep-ph].
  * (32) A. Bhardwaj, K. Bhide, T. Mandal, S. Mitra, and C. Neeraj, “Discovery prospects of a vectorlike top partner decaying to a singlet boson,” Phys. Rev. D 106 no. 7, (2022) 075024, arXiv:2204.09005 [hep-ph].
  * (33) J. Bardhan, T. Mandal, S. Mitra, and C. Neeraj, “Machine learning-enhanced search for a vectorlike singlet B quark decaying to a singlet scalar or pseudoscalar,” Phys. Rev. D 107 no. 11, (2023) 115001, arXiv:2212.02442 [hep-ph].
  * (34) A. Bhardwaj, T. Mandal, S. Mitra, and C. Neeraj, “Roadmap to explore vectorlike quarks decaying to a new scalar or pseudoscalar,” Phys. Rev. D 106 no. 9, (2022) 095014, arXiv:2203.13753 [hep-ph].
  * (35) Y.-B. Liu, B. Hu, and C.-Z. Li, “Single production of vectorlike quarks with charge 5/3 at the 14 TeV LHC,” arXiv:2402.01248 [hep-ph].
  * (36) A. C. Canbay and O. Cakir, “Investigating the single production of vectorlike quarks decaying into a top quark and W boson through hadronic channels at the HL-LHC,” Phys. Rev. D 108 no. 9, (2023) 095006, arXiv:2307.12883 [hep-ph].
  * (37) J.-Z. Han, S. Xu, W.-J. Mao, and H.-Q. Song, “Single production of vectorlike $X$ quarks at the HL-LHC,” Nucl. Phys. B 992 (2023) 116235.
  * (38) J. Z. Han, Y. B. Liu, L. Xing, and S. Xu, “Search for single production of vectorlike $B$-quarks at the LHC*,” Chin. Phys. C 46 no. 10, (2022) 103103, arXiv:2208.06845 [hep-ph].
  * (39) L. Han, L.-F. Du, and Y.-B. Liu, “Single production of vectorlike T quark at future high-energy linear e+e- collider,” Phys. Rev. D 105 no. 11, (2022) 115032, arXiv:2206.06562 [hep-ph].
  * (40) J. Z. Han, J. Yang, S. Xu, and H. K. Wang, “Search for single production of vector-like $B$ quark via bZ channel at the HL-LHC,” Nucl. Phys. B 975 (2022) 115672.
  * (41) J.-Z. Han, J. Yang, S. Xu, and H.-K. Wang, “Single production of vectorlike $B$ quarks at the CLIC,” Phys. Rev. D 105 no. 1, (2022) 015005, arXiv:2112.15044 [hep-ph].
  * (42) X.-Y. Tian, L.-F. Du, and Y.-B. Liu, “Search for single production of vectorlike top partners through $th$ channel at the HE-LHC and FCC-hh,” Eur. Phys. J. C 81 no. 7, (2021) 594.
  * (43) B. Yang, M. Wang, H. Bi, and L. Shang, “Single production of vectorlike $T$ quark decaying into $Wb$ at the LHC and the future $pp$ colliders,” Phys. Rev. D 103 no. 3, (2021) 036006.
  * (44) V. Cetinkaya, A. Ozansoy, V. Ari, O. M. Ozsimsek, and O. Cakir, “Single production of vectorlike Y quarks at the HL-LHC,” Nucl. Phys. B 973 (2021) 115580, arXiv:2012.15308 [hep-ph].
  * (45) H. Zhou and N. Liu, “Probing a new decay of a vectorlike top partner mediated by a heavy Majorana neutrino via single production,” Phys. Rev. D 101 no. 11, (2020) 115028, arXiv:2006.05296 [hep-ph].
  * (46) Y.-B. Liu, “Search for single production of the heavy vectorlike $T$ quark with $T\to th$ and $h\to\gamma\gamma$ at the high-luminosity LHC,” Phys. Rev. D 95 no. 3, (2017) 035013, arXiv:1612.05851 [hep-ph].
  * (47) C.-H. Chen and T. Nomura, “Single production of $X_{\pm 5/3}$ and $Y_{\mp 4/3}$ vectorlike quarks at the LHC,” Phys. Rev. D 94 no. 3, (2016) 035001, arXiv:1603.05837 [hep-ph].
  * (48) R. Barcelo, A. Carmona, M. Chala, M. Masip, and J. Santiago, “Single Vectorlike Quark Production at the LHC,” Nucl. Phys. B 857 (2012) 172–184, arXiv:1110.5914 [hep-ph].
  * (49) L. Shang, Y. Yan, S. Moretti, and B. Yang, “Single production of an exotic vector-like $Y$ quark at future high energy $pp$ colliders,” arXiv:2401.00770 [hep-ph].
  * (50) Y.-B. Liu and S. Moretti, “Search for single production of a top quark partner via the $T\to th$ and $h\to WW^{\ast}$ channels at the LHC,” Phys. Rev. D 100 no. 1, (2019) 015025, arXiv:1902.03022 [hep-ph].
  * (51) CMS Collaboration, A. M. Sirunyan et al., “A search for bottom-type, vector-like quark pair production in a fully hadronic final state in proton-proton collisions at $\sqrt{s}=$ 13 TeV,” Phys. Rev. D 102 (2020) 112004, arXiv:2008.09835 [hep-ex].
  * (52) ATLAS Collaboration, G. Aad et al., “Search for single vector-like $B$ quark production and decay via $B\rightarrow bH(b\bar{b})$ in $pp$ collisions at $\sqrt{s}=13$ TeV with the ATLAS detector,” arXiv:2308.02595 [hep-ex].
  * (53) M. Buchkremer, G. Cacciapaglia, A. Deandrea, and L. Panizzi, “Model Independent Framework for Searches of Top Partners,” Nucl. Phys. B 876 (2013) 376–417, arXiv:1305.4172 [hep-ph].
  * (54) D. Guest, K. Cranmer, and D. Whiteson, “Deep Learning and its Application to LHC Physics,” Ann. Rev. Nucl. Part. Sci. 68 (2018) 161–181, arXiv:1806.11484 [hep-ex].
  * (55) A. Butter et al., “The Machine Learning landscape of top taggers,” SciPost Phys. 7 (2019) 014, arXiv:1902.09914 [hep-ph].
  * (56) M. Kagan, “Image-Based Jet Analysis,” arXiv:2012.09719 [physics.data-an].
  * (57) C. W. Chiang, D. Shih, and S. F. Wei, “VBF vs. GGF Higgs with Full-Event Deep Learning: Towards a Decay-Agnostic Tagger,” Phys. Rev. D 107 no. 1, (2023) 016014, arXiv:2209.05518 [hep-ph].
  * (58) T. Chen and C. Guestrin, “XGBoost: A scalable tree boosting system,” in Proceedings of the 22nd ACM SIGKDD International Conference on Knowledge Discovery and Data Mining, KDD ’16, p. 785–794.  Association for Computing Machinery, New York, NY, USA, 2016.
  * (59) D. Alvestad et al., “Beyond cuts in small signal scenarios: Enhanced sneutrino detectability using machine learning,” Eur. Phys. J. C 83 no. 5, (2023) 379, arXiv:2108.03125 [hep-ph].
  * (60) F. Gianotti et al., “Physics potential and experimental challenges of the LHC luminosity upgrade,” Eur. Phys. J. C 39 (2005) 293–333, arXiv:hep-ph/0204087.
  * (61) A. Alloul, N. D. Christensen, C. Degrande, C. Duhr, and B. Fuks, “FeynRules 2.0 - A complete toolbox for tree-level phenomenology,” Comput. Phys. Commun. 185 (2014) 2250–2300, arXiv:1310.1921 [hep-ph].
  * (62) “Singlet b model (vlq),” https://feynrules.irmp.ucl.ac.be/wiki/VLQ_bsingletvl, 2014.  Accessed: February 2, 2024.
  * (63) H. J. He, Y. P. Kuang, and X. Y. Li, “On the precise formulation of equivalence theorem,” Phys. Rev. Lett. 69 (1992) 2619–2622.
  * (64) H. J. He, Y. P. Kuang, and X. Y. Li, “Further investigation on the precise formulation of the equivalence theorem,” Phys. Rev. D 49 (1994) 4842–4872.
  * (65) 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 (2014) 079, arXiv:1405.0301 [hep-ph].
  * (66) NNPDF Collaboration, R. D. Ball et al., “Parton distributions from high-precision collider data,” Eur. Phys. J. C 77 no. 10, (2017) 663, arXiv:1706.00428 [hep-ph].
  * (67) Particle Data Group Collaboration, R. L. Workman et al., “Review of Particle Physics,” PTEP 2022 (2022) 083C01.
  * (68) B. Belfatto and S. Trifinopoulos, “Cabibbo angle anomalies and oblique corrections: The remarkable role of the vectorlike quark doublet,” Phys. Rev. D 108 no. 3, (2023) 035022, arXiv:2302.14097 [hep-ph].
  * (69) A. Crivellin, M. Kirk, T. Kitahara, and F. Mescia, “Global fit of modified quark couplings to EW gauge bosons and vector-like quarks in light of the Cabibbo angle anomaly,” JHEP 03 (2023) 234, arXiv:2212.06862 [hep-ph].
  * (70) S. Balaji, “Asymmetry in flavour changing electromagnetic transitions of vector-like quarks,” JHEP 05 (2022) 015, arXiv:2110.05473 [hep-ph].
  * (71) G. C. Branco, J. T. Penedo, P. M. F. Pereira, M. N. Rebelo, and J. I. Silva-Marcos, “Addressing the CKM unitarity problem with a vector-like up quark,” JHEP 07 (2021) 099, arXiv:2103.13409 [hep-ph].
  * (72) D. Vatsyayan and A. Kundu, “Constraints on the quark mixing matrix with vector-like quarks,” Nucl. Phys. B 960 (2020) 115208, arXiv:2007.02327 [hep-ph].
  * (73) F. J. Botella, G. C. Branco, M. Nebot, M. N. Rebelo, and J. I. Silva-Marcos, “Vector-like Quarks at the Origin of Light Quark Masses and Mixing,” Eur. Phys. J. C 77 no. 6, (2017) 408, arXiv:1610.03018 [hep-ph].
  * (74) K. Ishiwata, Z. Ligeti, and M. B. Wise, “New Vector-Like Fermions and Flavor Physics,” JHEP 10 (2015) 027, arXiv:1506.03484 [hep-ph].
  * (75) A. K. Alok, S. Banerjee, D. Kumar, S. U. Sankar, and D. London, “New-physics signals of a model with a vector-singlet up-type quark,” Phys. Rev. D 92 (2015) 013002, arXiv:1504.00517 [hep-ph].
  * (76) G. Cacciapaglia, A. Deandrea, N. Gaur, D. Harada, Y. Okada, and L. Panizzi, “Interplay of vector-like top partner multiplets in a realistic mixing set-up,” JHEP 09 (2015) 012, arXiv:1502.00370 [hep-ph].
  * (77) F. J. Botella, G. C. Branco, and M. Nebot, “The Hunt for New Physics in the Flavour Sector with up vector-like quarks,” JHEP 12 (2012) 040, arXiv:1207.4440 [hep-ph].
  * (78) G. Cacciapaglia, A. Deandrea, D. Harada, and Y. Okada, “Bounds and Decays of New Heavy Vector-like Top Partners,” JHEP 11 (2010) 159, arXiv:1007.2933 [hep-ph].
  * (79) J. M. Campbell, R. K. Ellis, F. Maltoni, and S. Willenbrock, “Production of a $Z$ boson and two jets with one heavy-quark tag,” Phys. Rev. D 73 (2006) 054007, arXiv:hep-ph/0510362. [Erratum: Phys.Rev.D 77, 019903 (2008)].
  * (80) J. M. Campbell and R. K. Ellis, “An Update on vector boson pair production at hadron colliders,” Phys. Rev. D 60 (1999) 113006, arXiv:hep-ph/9905386.
  * (81) A. Alloul, N. D. Christensen, C. Degrande, C. Duhr, and B. Fuks, “FeynRules 2.0 - A complete toolbox for tree-level phenomenology,” Comput. Phys. Commun. 185 (2014) 2250–2300, arXiv:1310.1921 [hep-ph].
  * (82) C. Bierlich et al., “A comprehensive guide to the physics and usage of PYTHIA 8.3” arXiv:2203.11601 [hep-ph].
  * (83) DELPHES 3 Collaboration, J. de Favereau, C. Delaere, P. Demin, A. Giammanco, V. Lemaître, A. Mertens, and M. Selvaggi, “DELPHES 3, A modular framework for fast simulation of a generic collider experiment,” JHEP 02 (2014) 057, arXiv:1307.6346 [hep-ex].
  * (84) CERN Collaboration, M. Selvaggi, “Delphes cards for LHC Run-III, and HL-LHC.” December 6th, 2017.  https://github.com/delphes/delphes/blob/master/cards/delphes_card_HLLHC.tcl. Accessed on 2024-2-25.
  * (85) M. Cacciari, G. P. Salam, and G. Soyez, “FastJet User Manual,” Eur. Phys. J. C 72 (2012) 1896, arXiv:1111.6097 [hep-ph].
  * (86) M. Cacciari and G. P. Salam, “Dispelling the $N^{3}$ myth for the $k_{t}$ jet-finder,” Phys. Lett. B 641 (2006) 57–61, arXiv:hep-ph/0512210.
  * (87) E. Conte, B. Fuks, and G. Serret, “MadAnalysis 5, A User-Friendly Framework for Collider Phenomenology,” Comput. Phys. Commun. 184 (2013) 222–256, arXiv:1206.1599 [hep-ph].
  * (88) L. L. Shang and Y. Zhang, “EasyScan_HEP: a tool for connecting programs to scan the parameter space of physics models,” arXiv:2304.03636 [hep-ph].
  * (89) E. Burns and W. Fisher, “Testing the approximations described in ’Asymptotic formulae for likelihood-based tests of new physics’,” arXiv:1110.5002 [hep-ex].
  * (90) T. Fawcett, “An introduction to roc analysis,” Pattern Recognition Letters 27 no. 8, (2006) 861–874.
  * (91) J. Snoek, H. Larochelle, and R. P. Adams, “Practical Bayesian Optimization of Machine Learning Algorithms,” arXiv:1206.2944 [stat.ML].


§ INTRODUCTION

The study of microbial communities involves analyzing their biodiversity. The two axes of diversity typically examined are alpha diversity and beta diversity. Alpha diversity measures the biodiversity within a single community, while beta diversity measures the diversity between different communities. Measures of beta diversity allow the comparison of microbial communities and facilitate downstream analysis such as clustering and visualization.

The Bray-Curtis dissimilarity index <cit.> is commonly used to measure beta diversity by comparing microbial abundances. This index does not use external information about the compared operational taxonomic units (OTUs). The formula is simple and easy to calculate, and thus it is included in many studies.

UniFrac <cit.> is a distance metric over microbial abundance vectors that takes into account the phylogenetic relatedness between OTUs. It uses a user-provided phylogenetic tree and returns the Jaccard similarity (or dissimilarity) of the samples over the branches of the tree. This way, two different OTUs are in agreement over the part of the tree that they share, starting from the root and until they diverge to separate branches, and in disagreement over the parts of the tree past their divergence. The closer the OTUs are on the tree, the higher their similarity. Two flavors of UniFrac that are commonly used are weighted and unweighted. In weighted UniFrac, every species contributes to the total similarity relative to its abundance, while in unweighted UniFrac every species has the same weight.

The advantage of UniFrac is in the differential distances between the OTUs, according to their taxonomic relation. Similar environments in different geographies may harbor similar makeups of microbial families, orders, classes, etc. even if the individual species are different. These similarities are lost when using Bray-Curtis, but are kept when using UniFrac, allowing an informative comparison of environments.
Over time, improvements and variations to the UniFrac metric were developed. Generalized UniFrac <cit.>, Variance-Adjusted Weighted UniFrac <cit.>, Information UniFrac and Ratio UniFrac <cit.> were created to broaden the kinds of variance that can be measured. Performance improvements such as Fast UniFrac <cit.> were introduced to allow UniFrac to run on larger datasets. Recent advancements have introduced the Striped UniFrac <cit.> algorithm and GPU support <cit.>, which increased the speed of computation to make it applicable to datasets with over 100K samples.

These modern implementations come with specific challenges. R and Python-based implementations require having these languages installed, installing third-party libraries, and having some level of proficiency with these languages. Some of those third-party libraries, such as OpenMP and HDF5, require setting up the environment in which the program runs in a specific way. Bioconda’s implementation works with the BIOM and QZA formats, which are specific to that framework. Some of the implementations do not support Windows. Each such constraint adds friction on the way to getting UniFrac in the analysis pipeline. Scientists may give up on using UniFrac in their research if it requires more energy or expertise than what they have.

To demonstrate the growth potential for UniFrac, we searched Google Scholar for papers that mention Bray-Curtis without mentioning UniFrac, published in 2022 with source “nature”[Query: "bray curtis" -unifrac source:nature] (which included journals from the Nature Publishing Group) and got 562 results. Searching for papers that mention UniFrac, published in 2022 with source “nature”[Query: unifrac source:nature] yielded 282 results. Under the assumption that UniFrac is relevant to most of the studies mentioning Bray-Curtis, this suggests a potential to double the use of UniFrac.

Reviewing the latest 20 results from the latter search (papers that mention UniFrac), we found that UniFrac was used in studies with sample sizes ranging from 5 to 3211 samples. In 19 out of the 20 papers, sample sizes ranged from 5 to 247 samples, with an average of 69±66. Interestingly, 19 out of 20 of these papers used only the original variations of UniFrac - weighted and unweighted (supplementary table <ref>).

We infer that an inclusive implementation of UniFrac, one that is easily accessible, can substantially expand the usage of the UniFrac metric. From the surveyed studies we learn that supporting the weighted and unweighted variations, with reasonable performance on up to a few thousands of samples, can satisfy most currently published studies.

Here, we set out to create a more accessible implementation of the UniFrac metric. We created a toolkit that is optimized for ease of access and ease of use, and can be used on any platform without installation.

§ RESULTS

§.§ A standalone implementation

We evaluated our implementation against criteria that may affect a method’s usability (table <ref>). We refer to installation processes as “builtin” and “third party”, where a builtin installation is one that uses a language’s builtin package installation pipeline. Third-party installation mechanisms introduce more moving parts to the process, and are typically tested and maintained by smaller communities. This can lead to more fragile installation processes which take more debugging and therefore might exclude users.

Implementation     Supported operating systems     Required programming language     Required pre-existing software     Installation

FrackyFrac     All     None     None     None

Phyloseq <cit.>     All     R     R libraries     Third party, with BioConductor

QIIME2 <cit.>     Linux, macOS     Anaconda Python     >400 Python libraries     Builtin, with a conda configuration file

Scikit-bio     Linux, macOS     Python     Python libraries     Builtin, with pip

Bioconda-binaries (SSU) <cit.>     Linux, macOS     Anaconda Python     Python libraries     Builtin, with conda

Comparison of UniFrac implementations.

§.§ Performance

We ran FrackyFrac on the OTU tables of samples from the Earth Microbiome Project <cit.> (the “silva” reference was used). To compare run times on different input sizes, subsamples of the data were taken arbitrarily. We compared FrackyFrac with SSU which is the implementation used in QIIME2, and with Phyloseq. Scikit-bio was left out because QIIME2 is a faster Python library. FrackyFrac’s computation times were higher than SSU’s by a factor of about two on a single core, and by a factor of 1.4 on eight cores. Figure <ref> shows the performance of the compared implementations.

UniFrac computation times on Earth Microbiome Project data. (A) Computation times on 10, 20, 50, 100, 200, 500, 1000 and 2000 samples using a single core. (B) Computation times on 1000, 2000, 5000, 10000 and 23000 (all) samples using 8 cores.

§.§ Validation

FrackyFrac’s result distances from the previous section were compared with SSU’s and were found to be equal. The published repository contains a directory named ‘testdata’ with runnable tests and comparison scripts. Figure <ref> shows the separation of samples from several hosts from the Earth Microbiome Project, using FrackyFrac’s weighted and unweighted UniFrac.

PCoA of samples from different host species. Distances were calculated FrackyFrac’s unweighted (A) and weighted (B) UniFrac, on Earth Microbiome Project samples.

§.§ Creating phylogenetic trees

FrackyFrac features a builtin tree maker that uses the method suggested by <cit.>, Mash distance with UPGMA clustering, as a simple way to create phylogenetic trees. Since UniFrac requires a phylogenetic tree for its calculation, we provide this method in case no such tree is available.

§ CONCLUSIONS

In this work we introduce a standalone implementation of weighted and unweighted UniFrac. This evolved from a need to use UniFrac without having to go through software installation, format conversions, or having to switch programming languages or operating systems. FrackyFrac meets that need. Our implementation addresses mainly scientists who are not already versed in one of the existing frameworks for running UniFrac, but also scientists who have been using UniFrac and would like more freedom as to how and where they can use it. Faster alternatives exist, but for the sample sizes common in current studies, calculation speed has less impact than the human factor.

Some aspects are not covered by this work. In terms of scalability, it does not run on GPUs, making it less suitable for very large datasets. In terms of variation, it supports only weighted and unweighted UniFrac. As shown in the introduction, papers citing UniFrac tend to use the two original variations, so supporting only them covers the vast majority of use cases. Future improvements may include supporting additional flavors of UniFrac, and adding optional integration with third-party libraries for better scalability.

FrackyFrac has the potential to increase the use of UniFrac by reducing the amount of effort and expertise required to use it. It can expand the range of potential users to those with more diverse backgrounds and to more platforms.

§ IMPLEMENTATION

The FrackyFrac suite contains two executables: frcfrc which takes as input a TSV of abundances and a phylogenetic tree, and outputs pairwise UniFrac distances; and trtr which creates phylogenetic trees given input fasta sequences.

To support different platforms, FrackyFrac is cross-compiled into standalone binaries which contain all the required logic, with a native binary for every operating system. This removes the need for an installation in favor of downloading ready-to-use executables. There is no need to get FrackyFrac’s code in order to use it, so no particular programming language is required on the user’s side.

To support diverse pipelines, FrackyFrac uses the tab-separated-values (TSV) format for input and output, for which there exist encoders and decoders in most modern programming languages.

The algorithm for calculating UniFrac distances is similar to Fast UniFrac’s algorithm, where samples are represented as lists of tree-nodes and abundances. FrackyFrac uses the biostuff package suite <cit.> for processing biological data and the gostuff package suite <cit.> for algorithms and data structures.

§ AVAILABILITY AND REQUIREMENTS

Software: https://github.com/fluhus/frackyfrac. Operating systems: Linux, Windows, MacOS. Software requirements: none. License: MIT. Non-academic use is unrestricted.

§ DECLARATIONS

§.§ Availability of data and materials

The Earth Microbiome Project data that was used for demonstration can be downloaded from:

§.§ Competing interests

The authors declare that they have no competing interests.

§.§ Authors' contributions

AL: conceptualization, methodology, software, writing - original draft.
SS, AK: validation, data curation, writing - review & editing.
ES: supervision, writing - review & editing.

Missing 'biblatex' package
The bibliography requires the 'biblatex' package.

journaltitleEcological monographs
notePublisher: JSTOR
titleAn ordination of the upland forest communities of southern Wisconsin

journaltitleApplied and environmental microbiology
notePublisher: Am Soc Microbiol
titleUniFrac: a new phylogenetic method for comparing microbial communities

notePublisher: Oxford University Press
titleAssociating microbiome composition with environmental covariates using generalized UniFrac distances

journaltitleBMC bioinformatics
notePublisher: BioMed Central
titleVariance adjusted weighted UniFrac: a powerful beta diversity measure for comparing communities based on phylogeny

journaltitlePloS one
notePublisher: Public Library of Science San Francisco, CA USA
titleExpanding the UniFrac toolbox

journaltitleThe ISME journal
notePublisher: Nature Publishing Group
titleFast UniFrac: facilitating high-throughput phylogenetic analyses of microbial communities including analysis of pyrosequencing and PhyloChip data

journaltitleNature methods
notePublisher: Nature Publishing Group US New York
titleStriped UniFrac: enabling microbiome analysis at unprecedented scale

notePublisher: American Society for Microbiology 1752 N St., NW, Washington, DC
titleOptimizing UniFrac with OpenACC Yields Greater Than One Thousand Times Speed Increase

World Scientific

booktitleBiocomputing 2012
titlePhyloseq: a bioconductor package for handling and analysis of high-throughput phylogenetic sequence data

journaltitleNature biotechnology
notePublisher: Nature Publishing Group
titleReproducible, interactive, scalable and extensible microbiome data science using QIIME 2

journaltitleBMC biology
notePublisher: BioMed Central
titleThe Earth Microbiome project: successes and aspirations

journaltitleGenome biology
notePublisher: BioMed Central
titleMash: fast genome and metagenome distance estimation using MinHash

titlefluhus/biostuff: v0.1.19

titlefluhus/gostuff: v0.1.14

journaltitleScientific Reports
notePublisher: Nature Publishing Group UK London
titleWhite spot syndrome virus impact on the expression of immune genes and gut microbiome of black tiger shrimp Penaeus monodon

journaltitleScientific Reports
notePublisher: Nature Publishing Group UK London
titleObstructive sleep apnea is related to alterations in fecal microbiome and impaired intestinal barrier function

journaltitleScientific reports
notePublisher: Nature Publishing Group UK London
titleDiagnostic and prognostic potential of the microbiome in ovarian cancer treatment response

journaltitleNature Microbiology
notePublisher: Nature Publishing Group UK London
titleExacerbation of allergic rhinitis by the commensal bacterium Streptococcus salivarius

journaltitleScientific reports
notePublisher: Nature Publishing Group UK London
titleAltered gut microbiota in individuals with episodic and chronic migraine

journaltitleScientific Reports
notePublisher: Nature Publishing Group UK London
titleRectal swabs as a viable alternative to faecal sampling for the analysis of gut microbiota functionality and composition

journaltitleCommunications Biology
notePublisher: Nature Publishing Group UK London
titleHigh-fat diet and estrogen modulate the gut microbiota in a sex-dependent manner in mice

journaltitleTranslational Psychiatry
notePublisher: Nature Publishing Group UK London
titleInteractions between perceived stress and microbial-host immune components: two demographically and geographically distinct pregnancy cohorts

journaltitleScientific Reports
notePublisher: Nature Publishing Group UK London
titleSevere, short-term sleep restriction reduces gut microbiota community richness but does not alter intestinal permeability in healthy young men

journaltitleExperimental & Molecular Medicine
notePublisher: Nature Publishing Group UK London
titleStable colonization of Akkermansia muciniphila educates host intestinal microecology and immunity to battle against inflammatory intestinal diseases

journaltitleScientific Reports
notePublisher: Nature Publishing Group UK London
titleLong-term taxonomic and functional stability of the gut microbiome from human fecal samples

journaltitleScientific Reports
notePublisher: Nature Publishing Group UK London
titleThe rhizospheric bacterial diversity of Fritillaria taipaiensis under single planting pattern over five years

journaltitleScientific Reports
notePublisher: Nature Publishing Group UK London
titleMicrobial rewilding in the gut microbiomes of captive ring-tailed lemurs (Lemur catta) in Madagascar

journaltitleScientific Reports
notePublisher: Nature Publishing Group UK London
titleFecal level of butyric acid, a microbiome-derived metabolite, is increased in patients with severe carotid atherosclerosis

journaltitlenpj Biofilms and Microbiomes
notePublisher: Nature Publishing Group UK London
titleGut microbiota as an antioxidant system in centenarians associated with high antioxidant activities of gut-resident Lactobacillus

journaltitleMolecular Psychiatry
notePublisher: Nature Publishing Group UK London
titleHelicobacter hepaticus augmentation triggers Dopaminergic degeneration and motor disorders in mice with Parkinson’s disease

journaltitleThe ISME journal
notePublisher: Nature Publishing Group UK London
titlePlant microbiomes harbor potential to promote nutrient turnover in impoverished substrates of a Brazilian biodiversity hotspot

journaltitleScientific reports
notePublisher: Nature Publishing Group UK London
titleIncreased levels of oral Streptococcus-derived d-alanine in patients with chronic kidney disease and diabetes mellitus

journaltitleScientific reports
notePublisher: Nature Publishing Group UK London
titleDynamics of rumen microbiome in sika deer (Cervus nippon yakushimae) from unique subtropical ecosystem in Yakushima Island, Japan

journaltitleNature Communications
notePublisher: Nature Publishing Group UK London
titleThe gut microbiota and depressive symptoms across ethnic groups

§ SUPPLEMENTARY MATERIAL

Sample size     Type of UniFrac used     Paper title

12     W     White spot syndrome virus impact on the expression of immune genes and gut microbiome of black tiger shrimp Penaeus monodon <cit.>

48     U     Obstructive sleep apnea is related to alterations in fecal microbiome and impaired intestinal barrier function <cit.>

64     U,W,Generalized     Diagnostic and prognostic potential of the microbiome in ovarian cancer treatment response <cit.>

160     W     Exacerbation of allergic rhinitis by the commensal bacterium Streptococcus salivarius <cit.>

130     U,W     Altered gut microbiota in individuals with episodic and chronic migraine <cit.>

20     U     Rectal swabs as a viable alternative to faecal sampling for the analysis of gut microbiota functionality and composition <cit.>

64     U,W     High-fat diet and estrogen modulate the gut microbiota in a sex-dependent manner in mice <cit.>

84     Not specified     Interactions between perceived stress and microbial-host immune components: two demographically and geographically distinct pregnancy cohorts <cit.>

19     U,W     Severe, short-term sleep restriction reduces gut microbiota community richness but does not alter intestinal permeability in healthy young men <cit.>

21     W     Stable colonization of Akkermansia muciniphila educates host intestinal microecology and immunity to battle against inflammatory intestinal diseases <cit.>

10     U,W     Long-term taxonomic and functional stability of the gut microbiome from human fecal samples <cit.>

5     Not specified     The rhizospheric bacterial diversity of Fritillaria taipaiensis under single planting pattern over five years <cit.>

172     U,W     Microbial rewilding in the gut microbiomes of captive ring-tailed lemurs (Lemur catta) in Madagascar <cit.>

81     U     Fecal level of butyric acid, a microbiome-derived metabolite, is increased in patients with severe carotid atherosclerosis <cit.>

247     U     Gut microbiota as an antioxidant system in centenarians associated with high antioxidant activities of gut-resident Lactobacillus <cit.>

20     U     Helicobacter hepaticus augmentation triggers Dopaminergic degeneration and motor disorders in mice with Parkinson’s disease <cit.>

84     W     Plant microbiomes harbor potential to promote nutrient turnover in impoverished substrates of a Brazilian biodiversity hotspot <cit.>

17     U,W     Increased levels of oral Streptococcus-derived d-alanine in patients with chronic kidney disease and diabetes mellitus <cit.>

49     W     Dynamics of rumen microbiome in sika deer (Cervus nippon yakushimae) from unique subtropical ecosystem in Yakushima Island, Japan <cit.>

3211     W     The gut microbiota and depressive symptoms across ethnic groups <cit.>

Latest 20 papers mentioning UniFrac with source “nature” on Google Scholar (Jan 2023). UniFrac types are weighted (W) and unweighted (U).