#  A Simple Approach to Jointly Rank Passages and Select Relevant Sentences in
the OBQA Context

Man Luo
Arizona State University
<EMAIL_ADDRESS>
&Shuguang Chen
University of Houston
<EMAIL_ADDRESS>
&Chitta Baral
Arizona State University
<EMAIL_ADDRESS>
###### Abstract

In the open question answering (OBQA) task, how to select the relevant
information from a large corpus is a crucial problem for reasoning and
inference. Some datasets (e.g, HotpotQA) mainly focus on testing the model’s
reasoning ability at the sentence level. To overcome this challenge, many
existing frameworks use a deep learning model to select relevant passages and
then answer each question by matching a sentence in the corresponding passage.
However, such frameworks require long inference time and fail to take
advantage of the relationship between passages and sentences. In this work, we
present a simple yet effective framework to address these problems by jointly
ranking passages and selecting sentences. We propose consistency and
similarity constraints to promote the correlation and interaction between
passage ranking and sentence selection. In our experiments, we demonstrate
that our framework can achieve competitive results and outperform the baseline
by 28% in terms of exact matching of relevant sentences on the HotpotQA
dataset.

##  1 Introduction

Open book question answering (OBQA) task aims to answer questions with respect
to a given set of open book facts and common knowledge. This task has led to
many useful and practical Natural Language Processing (NLP) applications such
as search engines (Kwiatkowski et al., 2019) and dialogue systems (Reddy et
al., 2019; Choi et al., 2018). Among several OBQA datasets (Dhingra et al.,
2017; Mihaylov et al., 2018; Khot et al., 2020), HotpotQA Yang et al. (2018)
is more challenging because it requires 1) deeper reasoning with linguistic
understanding across multiple passages, and 2) a system to output supporting
facts to reach the answer.

Figure 1: An example from the HotpotQA dataset, where the question should be
answered by combining supporting facts(SP) from two passages. In the SP, the
first string refers to the title of passage, and the second integer means the
index of the sentence.

To address this challenge, previous work (nie2019revealing; Fang et al., 2020;
Tu et al., 2019; Groeneveld et al., 2020) mainly focuses on how to select the
most relevant sentences with a two-step pipeline: identify the most relevant
passage and then match each question with a single sentence in the
corresponding passage. This approach requires long inference time and neglects
the relationship between passages and sentences. Intuitively, the passage
ranking and selecting relevant sentences affect each other. If a passage has a
high score, then some sentences in the passage should be selected as relevant.
On the other hand, if a passage is scored low, then all sentences in the
passage should be classified as irrelevant.

In this paper, we describe a Two-in-One framework that aims to promote the
correlation and interaction between passages and sentences. We use both
passage ranking and sentence selection to jointly train the model. More
specifically, the model, based on the transformer (Vaswani et al., 2017)
structure, generates passage representations and sentence representation
simultaneously. Such representations are then fed to a passage ranker and
sentence classifier, respectively. Then we promote the interaction between
passage ranking and sentence classification using consistency and similarity
constraints. The consistency constraint is to enforce that the relevant
passage includes relevant sentences, while the similarity constraint enables
the model to represent relevant passages based on the relevant sentence
instead of irrelevant ones. By conducting experiments on the HotpotQA
datasets, we demonstrate that our framework, with a single model and fewer
parameters, can achieve competitive results and boost the performance of
baselines by 28%.

To summarize, we make the following contributions:

  1. 1.

We study the joint learning of passage ranking and sentence selection with
consistency and similarity constraints to model the relationship between
passages and sentences.

  2. 2.

We conduct experiments on the HotpotQA datasets and demonstrate the
effectiveness of our approach for deeper reasoning in the context of the OBQA
task.

##  2 Related Work

#### HotpotQA Systems

A straightforward way to solve the HotpotQA challenge is to build a
hierarchical system nie2019revealing, meaning a system first ranks relevant
passages and then identifies relevant sentences from the selected passages.
Such a hierarchical system involves multiple models thus requires long
inference time. More importantly, such a system only leverages the impact of
passage ranking on sentence selection but ignores the influence of the
sentence selection on the passage ranking. Our framework achieves these two
tasks by one model and facilitates the interaction by two constraints.
Groeneveld et al. (2020) proposes a pipeline based on three BERT models Devlin
et al. (2019) to solve the HotpotQA challenge. The system first selects
relevant sentences and then detects the answer span, finally, identifies the
relevant sentences according to the answer span. Though the pipeline is
strong, the way it solves the problem is opposite to human beings. We, humans,
identify the relevant sentences, and then give the answer span. Many existing
works demonstrate the effectiveness of graph neural networks(GNN) on HotpotQA
challenge (Fang et al., 2020; Tu et al., 2019). Since GNN is out of the scope
of this work, we do not compare it with these frameworks.

#### Joint Model for QA

Joint learning has been studied in Question Answering Tasks. Deng et al.
(2020) proposes a joint model to tackle community question answering such that
the model can simultaneously select the set of correct answers from candidates
and generate an abstractive summary for each selected answer. Sun et al.
(2019) proposes a generative collaborative network to answer questions and
generate questions. The main difference between our work and previous ones are
in two sense (1) our proposed model uses the shared encoder to tackle two
classification tasks (2) besides the loss function to optimize individual
tasks, we also propose two constraints that utilize the relation between these
two tasks.

##  3 HotpotQA Dataset

HotpotQA dataset (Yang et al., 2018) is designed for multi-hop reasoning
question answering tasks, i.e. to reason over multiple documents and answer
questions (see Figure 1). Particularly, HotpotQA challenge requires reasoning
over two passages. Furthermore, to guide the system to perform meaningful and
explainable reasoning, the dataset also provides supporting facts(SP) that
reach the answer to the question. HotpotQA provide two challenging settings:
in Fullwiki setting, a system needs to rank passage from the entire wiki
corpus; in Distractor setting, 10 distracting passages (including relevant
ones) are given for each question. In this work, we mainly focus on the latter
setting. From the training set, we find that 70.4% questions have exactly two
supporting facts (SP), and 60.0% of SP are the first sentence of passages.

##  4 Method

###  4.1 Task Formulation

We aim to jointly conduct two tasks, passage ranking and supporting facts
selection for HotpotQA. Given a question Q, the goal is to simultaneously rank
the set of candidates A = {$a_{1}$, …, $a_{i}$} and identify the supporting
facts for the TopK111The value of K depends on the task, and for HotpotQA, K
is 2. passages.

###  4.2 Model: Two-in-One Framework

We introduce the proposed joint model for passage ranking and support fact
selection, Two-in-One, which uses state-of-the-art transformer-based model
(Vaswani et al., 2017) to encode questions and contexts. In this work, we use
RoBERTa (Liu et al., 2019), however, any other variants like ELECTRA (Clark et
al., 2020) can be applied in this framework. The model architecture is given
in Figure 2. On top of the encoder, there are two MLP layers to score passages
and sentences, respectively. In details, given a question and a passage, we
firstly create an input to feed through RoBERTa Liu et al. (2019) by
concatenating the question and the passage as follows, $\langle s\rangle
Q\langle/s\rangle S_{1}\langle/s\rangle S_{2}...\langle/s\rangle
S_{k}\langle/s\rangle$ where $\langle s\rangle$ and $\langle/s\rangle$ are
special tokens in RoBERTa, $S_{i}$ is the $i^{th}$ sentence from a passage. We
take $\langle s\rangle$ as the contextual representation for passage ranking
and the $\langle/s\rangle$ in front of each sentence for sentence selection.
The passage ranker and the sentence classifier have identical structure but
different weights, which are two-layer Multiple-Layer Perceptron(MLP) with the
last layer predicting a logit.

Figure 2: The architecture of Two-in-One model for passage ranking and
relevant sentence selection. For HotpotQA dataset, K is two.

The model is jointly trained by passage loss and sentence loss. In detail,
during the training time, we assign the relevant passages and sentences with
ground truth score 1 while irrelevant passages and sentences with ground truth
score -1. Then, Mean Square Error(MSE) loss is applied to calculate the
passage and sentence loss as follows,

$\begin{split}&\mathcal{L}^{pass}=(\hat{y}-y)^{2},\\\
&\mathcal{L}^{sent}=\sum_{i=1}^{K}(\hat{x_{i}}-x_{i})^{2},\\\
&\mathcal{L}^{joint}=\mathcal{L}^{pass}+\mathcal{L}^{sent},\end{split}$ (1)

where $\hat{y}$ is the predicted passage score, $y$ is the ground truth score
of the passage, $\hat{x_{i}}$ and $x_{i}$ are the predicted sentence score and
ground truth score of $S_{i}$, respectively, and $K$ is the total number of
sentences in the passage. We simply sum up the passage loss and sentence loss
to jointly update model parameters.

During the inference time, passages are ranked based on the logits given by
the passage ranker. For the sentence classification, we take 0222The reason
for threshold “0” is that it is the middle value of 1 and -1, which are labels
for relevant and irrelevant sentences in the training time. as the threshold
to classify the relevance of each sentence: if the logit given by the sentence
classifier is larger than 0, then it is relevant; otherwise, irrelevant.

Next, we introduce two constraints to facilitate the interaction between these
two tasks.

###  4.3 Consistency Constraint

Intuitively, if a passage is relevant to the question, then there are some
sentences from the passages that are relevant; on the other hand, if a passage
is not relevant to the answer, then there should not be relevant sentences
inside the passage. Thus, we propose a consistency constraint over the passage
ranker and sentence classifier to minimize the gap between the passage score
and the maximum sentence score. The loss function is as follows:

$\mathcal{L}^{con}=(\hat{y}-max(\mathbf{x}))^{2},$ (2)

where $\mathbf{x}=[\hat{x}_{1}\dots\hat{x}_{n}]$ denotes a stack of predicted
sentence scores.

###  4.4 Similarity Constraint

As we have shown at the beginning of this section, the $\langle s\rangle$ is
used to get the passage score, and each $\langle/s\rangle$ is used to get the
sentence score. Intuitively, the similarity between $\langle s\rangle$ of a
relevant passage is more close to $\langle/s\rangle$ of a relevant sentence
than to $\langle/s\rangle$ of any irrelevant sentence. To enforce this
constraint, we use triplet as follows:

$\begin{split}\mathcal{L}^{sim}=\frac{1}{N\cdot
M}\sum_{i=1}^{N}\sum_{j=1}^{M}(&max\\{d(v^{p},v^{r}_{i})\\\
&-d(v^{p},v^{n}_{j})+m,0\\}),\end{split}$ (3)

where $d(\cdot,\cdot)$ is the Euclidian similarity, $N$ is the number of
relevant sentences, $M$ is the number of irrelevant sentences,
$v^{p},v^{r},v^{n}$ is the vector representation of the relevant passage,
relevant sentence, and irrelevant sentence respectively. In practice, we set
the margin $m$ at 1 and find optimum results. We train our model in an end-to-
end fashion by combining $\mathcal{L}^{joint}$, $\mathcal{L}^{con}$ and
$\mathcal{L}^{dis}$.

Model | # Parameters | SP Precision | SP Recall | SP F1 | SP EM | Passage EM
---|---|---|---|---|---|---
Sentence Selection Baseline |  $\sim$330M | 67.96 | 81.05 | 72.02 | 28.12 | 69.70
Passage Selection Baseline |  $\sim$330M | 66.43 | 56.55 | 60.20 | 27.30 | 90.44
Two-in-One + sim |  $\sim$330M | 88.06 | 85.68 | 85.82 | 59.17 | 91.11
QUARK |  $\sim$1020M∗ | N/A | N/A | 86.97 | 60.72 | N/A

Table 1:  The Results for two baselines and Two-in-One model with similarity
constraint on dev set of HotpotQA distracting dataset. The last raw, QUARK
Groeneveld et al. (2020), is the result of a framework with 3 BERT models. SP
stands for supporting facts and EM for Exact Match. ∗ refers to estimation.

##  5 Experiment

In this section, we first describe the training setup, and then introduce two
baselines. We evaluate the two baselines and our proposed joint model on the
HotpotQA dataset. Yang et al. (2018) provides two metrics for supporting facts
evaluation, exact matching (EM) and F1 score. We also present the precision
and recall of SP, and the exact matching of passages for detailed comparison.
Meanwhile, we compare our model with the QUARK system Groeneveld et al.
(2020). Lastly, we conduct an ablation study to show the effectiveness of the
proposed similarity loss and consistent loss.

###  5.1 Experiment Setup

We use Huggingface (Wolf et al., 2020) and Pytorch (Paszke et al., 2019)
libraries to implement each model. For all experiments, we use 4 TX1080 and
V100 NVIDIA, and all models are trained in 5 epochs with a learning rate of
1e-5, batch size of 32. We set the maximum input length in training to be 512.
We use 1 GPU to test models with a batch size of 16, and use the maximum
length as training to encode the input.

###  5.2 Baseline

To have comparable size of the model, two baselines have similar structure as
our Two-in-One model. Specifically, Two-in-One model has two classification
heads, whereas each of the baselines has only one classification head. One
baseline is to select relevant sentences, and the other one is to rank
passages.

#### Sentence Selection Baseline

The first baseline is to select relevant sentence, and particularly, we use a
RoBERTa-large with an additional MLP trained on question and a single
sentence: $\langle s\rangle Q\langle/s\rangle S\langle/s\rangle$, where $Q$ is
a question and $S$ is a sentence. Although this model can not predict the
relevant passage directly, based on the assumption that relevant passages
include relevant sentences, we pick up two relevant passages based on the top2
sentence scores. When the top1 and the top2 sentences are from the same
passage, we continue searching based on the ranking sentence scores until the
second document comes up. Then the supporting facts are those sentences from
the relevant documents with a score larger than 0.

#### Passage Selection Baseline

In the second baseline, again, we use RoBERTa-large but with the goal of
passage selection. The input to the model is a question and a passage:
$\langle s\rangle Q\langle/s\rangle P\langle/s\rangle$. Since such a model can
not predict sentence relevancy score, based on the statistic of HotpotQA that
majority of training set has two supporting facts and the most of them are the
first sentences in a paragraph (see Section 3), we select supporting facts by
the first sentence of the top1 and top2 passages.

###  5.3 Result

As we see from Table 1, Two-in-One framework outperforms two baselines with
large-margin improvement in all metrics, especially we see a significant
improvement on the EM of SP. Our framework outperforms the Sentence Selection
Baseline by 20% and 4.5% improvement on the precision and recall of SP,
respectively, which demonstrates that jointly learning is beneficial for
sentence classification. Also, jointly learning benefits for the passage
ranking by comparing Two-in-One with Passage Selection Baseline on the EM of
passage. Besides, we also compare Two-in-One with QUARK Groeneveld et al.
(2020), a framework involving three BERT models, (roughly three times larger
than ours). Two-in-One achieves comparable results in terms of F1 and EM of SP
regardless of much less parameters in our system. Notice that we do not have
the other three values because they are not presented in their original paper.

###  5.4 Ablation

To evaluate the impacts of the consistency constraint and the similarity
constraint, we conduct experiments with and without constraints. From Table 2,
we see that both consistency constraint and similarity constraint improve F1
and EM of SP and the similarity constraint also improves the EM of passages.
We found that without any constraint, though the model can rank the passages
well, it suffers from distinguishing between close sentences. The similarity
constraint addresses this issue in some sense by maximizing the distance
between relevant and irrelevant sentences.

Model | SP F1 | SP EM | Passage EM
---|---|---|---
Two-in-One | 85.52 | 58.67 | 90.93
Two-in-One + con | 85.55 | 58.98 | 90.29
Two-in-One + sim | 85.82 | 59.17 | 91.11
Two-in-One + con + sim | 85.63 | 58.74 | 90.78

Table 2: The results for Two-in-One model with or without consistency and
similarity constraints.

To better understand the impact of consistency constraint, we analyze the
consistency between the passage score and the sentence score. The prediction
of a model is consistent if the passage score agrees with the sentence scores
and the agreement can be measured by the gap between the passage score and the
maximum sentence score among all sentences in that passage. We observe that by
adding the consistency constraint, the gap between the passage score and the
sentence score is much smaller than without the consistency constraint as
shown in Table 3. It demonstrates that the constraint is beneficial for
consistent prediction.

Model | Mean Gap
---|---
Two-in-One | 0.11
Two-in-One + con | 0.03

Table 3: The comparison between the mean gap of passage score and sentence
score.

##  6 Discussion

In this work, we only evaluate the Two-in-One model on HotpotQA dataset since
to our knowledge, it is the only publicly available dataset that provides
supporting facts to a question. Nevertheless, evaluation on one dataset limits
the generalization capacity of our model. A potential solution to address the
dataset issue is selecting the supporting facts based on the answer.
Specifically, given a passage, a question, and an answer, if any sentence in
the passage contains the answer, then it is a supporting fact, otherwise, it
is not. With this assumption, we can evaluate Two-in-One models on more QA
datasets such as NaturalQuestions (Kwiatkowski et al., 2019).

##  7 Conclusion

In this work, we present a Two-in-One framework to answer questions with
factoid knowledge in the task of the OBQA. By jointly training with passage
ranking and sentence selection, the model is capable of capturing the
correlation between passages and sentences and reasoning over complex context
with a deeper linguistic understanding. We show the effectiveness of our
proposed framework by evaluating the model performance on the HotpotQA
datasets, concluding that jointly modeling passage ranking and sentence
selection is beneficial for the task of OBQA. Compared to the existing QA
systems, our framework, with fewer parameters, can achieve competitive
results. In the future, we will explore the relationship among passages,
supporting sentences, and answers in modeling and generalize our method on
more datasets.

## References

  * Choi et al. (2018) Eunsol Choi, He He, Mohit Iyyer, Mark Yatskar, Wen-tau Yih, Yejin Choi, Percy Liang, and Luke Zettlemoyer. 2018.  QuAC: Question answering in context.  In _Proceedings of the 2018 Conference on Empirical Methods in Natural Language Processing_ , pages 2174–2184, Brussels, Belgium. Association for Computational Linguistics.
  * Clark et al. (2020) Kevin Clark, Minh-Thang Luong, Quoc V. Le, and Christopher D. Manning. 2020.  ELECTRA: pre-training text encoders as discriminators rather than generators.  In _8th International Conference on Learning Representations, ICLR 2020, Addis Ababa, Ethiopia, April 26-30, 2020_. OpenReview.net.
  * Deng et al. (2020) Yang Deng, Wai Lam, Yuexiang Xie, Daoyuan Chen, Yaliang Li, Min Yang, and Ying Shen. 2020.  Joint learning of answer selection and answer summary generation in community question answering.  In _The Thirty-Fourth AAAI Conference on Artificial Intelligence, AAAI 2020, The Thirty-Second Innovative Applications of Artificial Intelligence Conference, IAAI 2020, The Tenth AAAI Symposium on Educational Advances in Artificial Intelligence, EAAI 2020, New York, NY, USA, February 7-12, 2020_ , pages 7651–7658. AAAI Press.
  * Devlin et al. (2019) Jacob Devlin, Ming-Wei Chang, Kenton Lee, and Kristina Toutanova. 2019.  BERT: Pre-training of deep bidirectional transformers for language understanding.  In _Proceedings of the 2019 Conference of the North American Chapter of the Association for Computational Linguistics: Human Language Technologies, Volume 1 (Long and Short Papers)_ , pages 4171–4186, Minneapolis, Minnesota. Association for Computational Linguistics.
  * Dhingra et al. (2017) Bhuwan Dhingra, Kathryn Mazaitis, and William W. Cohen. 2017.  Quasar: Datasets for question answering by search and reading.  _ArXiv preprint_ , abs/1707.03904.
  * Fang et al. (2020) Yuwei Fang, Siqi Sun, Zhe Gan, Rohit Pillai, Shuohang Wang, and Jingjing Liu. 2020\.  Hierarchical graph network for multi-hop question answering.  In _Proceedings of the 2020 Conference on Empirical Methods in Natural Language Processing (EMNLP)_ , pages 8823–8838, Online. Association for Computational Linguistics.
  * Groeneveld et al. (2020) Dirk Groeneveld, Tushar Khot, Mausam, and Ashish Sabharwal. 2020.  A simple yet strong pipeline for HotpotQA.  In _Proceedings of the 2020 Conference on Empirical Methods in Natural Language Processing (EMNLP)_ , pages 8839–8845, Online. Association for Computational Linguistics.
  * Khot et al. (2020) Tushar Khot, Peter Clark, Michal Guerquin, Peter Jansen, and Ashish Sabharwal. 2020\.  QASC: A dataset for question answering via sentence composition.  In _The Thirty-Fourth AAAI Conference on Artificial Intelligence, AAAI 2020, The Thirty-Second Innovative Applications of Artificial Intelligence Conference, IAAI 2020, The Tenth AAAI Symposium on Educational Advances in Artificial Intelligence, EAAI 2020, New York, NY, USA, February 7-12, 2020_ , pages 8082–8090. AAAI Press.
  * Kwiatkowski et al. (2019) Tom Kwiatkowski, Jennimaria Palomaki, Olivia Redfield, Michael Collins, Ankur Parikh, Chris Alberti, Danielle Epstein, Illia Polosukhin, Jacob Devlin, Kenton Lee, Kristina Toutanova, Llion Jones, Matthew Kelcey, Ming-Wei Chang, Andrew M. Dai, Jakob Uszkoreit, Quoc Le, and Slav Petrov. 2019.  Natural questions: A benchmark for question answering research.  _Transactions of the Association for Computational Linguistics_ , 7:452–466.
  * Liu et al. (2019) Yinhan Liu, Myle Ott, Naman Goyal, Jingfei Du, Mandar Joshi, Danqi Chen, Omer Levy, Mike Lewis, Luke Zettlemoyer, and Veselin Stoyanov. 2019.  Roberta: A robustly optimized bert pretraining approach.  _ArXiv preprint_ , abs/1907.11692.
  * Mihaylov et al. (2018) Todor Mihaylov, Peter Clark, Tushar Khot, and Ashish Sabharwal. 2018.  Can a suit of armor conduct electricity? a new dataset for open book question answering.  In _Proceedings of the 2018 Conference on Empirical Methods in Natural Language Processing_ , pages 2381–2391, Brussels, Belgium. Association for Computational Linguistics.
  * Paszke et al. (2019) Adam Paszke, Sam Gross, Francisco Massa, Adam Lerer, James Bradbury, Gregory Chanan, Trevor Killeen, Zeming Lin, Natalia Gimelshein, Luca Antiga, Alban Desmaison, Andreas Köpf, Edward Yang, Zachary DeVito, Martin Raison, Alykhan Tejani, Sasank Chilamkurthy, Benoit Steiner, Lu Fang, Junjie Bai, and Soumith Chintala. 2019.  Pytorch: An imperative style, high-performance deep learning library.  In _Advances in Neural Information Processing Systems 32: Annual Conference on Neural Information Processing Systems 2019, NeurIPS 2019, December 8-14, 2019, Vancouver, BC, Canada_ , pages 8024–8035.
  * Reddy et al. (2019) Siva Reddy, Danqi Chen, and Christopher D. Manning. 2019.  CoQA: A conversational question answering challenge.  _Transactions of the Association for Computational Linguistics_ , 7:249–266.
  * Sun et al. (2019) Yibo Sun, Duyu Tang, Nan Duan, Tao Qin, Shujie Liu, Zhao Yan, Ming Zhou, Yuanhua Lv, Wenpeng Yin, Xiaocheng Feng, et al. 2019.  Joint learning of question answering and question generation.  _IEEE Transactions on Knowledge and Data Engineering_ , 32(5):971–982.
  * Tu et al. (2019) Ming Tu, Kevin Huang, Guangtao Wang, Jing Huang, X. He, and Bowen Zhou. 2019.  Select, answer and explain: Interpretable multi-hop reading comprehension over multiple documents.  _ArXiv preprint_ , abs/1911.00484.
  * Vaswani et al. (2017) Ashish Vaswani, Noam Shazeer, Niki Parmar, Jakob Uszkoreit, Llion Jones, Aidan N. Gomez, Lukasz Kaiser, and Illia Polosukhin. 2017.  Attention is all you need.  In _Advances in Neural Information Processing Systems 30: Annual Conference on Neural Information Processing Systems 2017, December 4-9, 2017, Long Beach, CA, USA_ , pages 5998–6008.
  * Wolf et al. (2020) Thomas Wolf, Lysandre Debut, Victor Sanh, Julien Chaumond, Clement Delangue, Anthony Moi, Pierric Cistac, Tim Rault, Remi Louf, Morgan Funtowicz, Joe Davison, Sam Shleifer, Patrick von Platen, Clara Ma, Yacine Jernite, Julien Plu, Canwen Xu, Teven Le Scao, Sylvain Gugger, Mariama Drame, Quentin Lhoest, and Alexander Rush. 2020.  Transformers: State-of-the-art natural language processing.  In _Proceedings of the 2020 Conference on Empirical Methods in Natural Language Processing: System Demonstrations_ , pages 38–45, Online. Association for Computational Linguistics.
  * Yang et al. (2018) Zhilin Yang, Peng Qi, Saizheng Zhang, Yoshua Bengio, William Cohen, Ruslan Salakhutdinov, and Christopher D. Manning. 2018.  HotpotQA: A dataset for diverse, explainable multi-hop question answering.  In _Proceedings of the 2018 Conference on Empirical Methods in Natural Language Processing_ , pages 2369–2380, Brussels, Belgium. Association for Computational Linguistics.


Belle II Preprint 2024-012

KEK Preprint 2024-6

I. Adachi ,

# Search for lepton-flavor-violating
$\tau^{-}\rightarrow\mu^{-}\mu^{+}\mu^{-}$ decays at Belle II

0000-0003-2287-0173 L. Aggarwal ,  0000-0002-0909-7537 H. Aihara ,
0000-0002-1907-5964 N. Akopov  0000-0002-4425-2096 A. Aloisio
0000-0002-3883-6693 N. Althubiti  0000-0003-1513-0409 N. Anh Ky ,
0000-0003-0471-197X D. M. Asner  0000-0002-1586-5790 H. Atmacan
0000-0003-2435-501X V. Aushev ,  0000-0002-8588-5308 M. Aversano
0000-0001-9980-0953 R. Ayad ,  0000-0003-3466-9290 V. Babu ,
0000-0003-0419-6912 H. Bae ,  0000-0003-1393-8631 S. Bahinipati
0000-0002-3744-5332 P. Bambade ,  0000-0001-7378-4852 Sw. Banerjee
0000-0001-8852-2409 S. Bansal  0000-0003-1992-0336 M. Barrett
0000-0002-2095-603X J. Baudot ,  0000-0001-5585-0991 A. Baur ,
0000-0003-1360-3292 A. Beaubien  0000-0001-9438-089X F. Becherer
0000-0003-0562-4616 J. Becker  0000-0002-5082-5487 J. V. Bennett ,
0000-0002-5440-2668 F. U. Bernlochner ,  0000-0001-8153-2719 V. Bertacchi
0000-0001-9971-1176 M. Bertemes  0000-0001-5038-360X E. Bertholet
0000-0002-3792-2450 M. Bessner ,  0000-0003-1776-0439 S. Bettarini
0000-0001-7742-2998 F. Bianchi ,  0000-0002-1524-6236 L. Bierwirth
0009-0003-0192-9073 T. Bilka  0000-0003-1449-6986 D. Biswas ,
0000-0002-7543-3471 A. Bobrov ,  0000-0001-5735-8386 D. Bodrov
0000-0001-5279-4787 A. Bolz  0000-0002-4033-9223 J. Borah ,
0000-0003-2990-1913 A. Boschetti  0000-0001-6030-3087 A. Bozek
0000-0002-5915-1319 M. Bračko  0000-0002-2495-0524 P. Branchini
0000-0002-2270-9673 T. E. Browder ,  0000-0001-7357-9007 A. Budano ,
0000-0002-0856-1131 S. Bussino ,  0000-0002-3829-9592 Q. Campagna ,
0000-0002-3109-2046 M. Campajola  0000-0003-2518-7134 L. Cao
0000-0001-8332-5668 G. Casarosa  0000-0003-4137-938X C. Cecchi
0000-0002-2192-8233 J. Cerasoli  0000-0001-9777-881X M.-C. Chang
0000-0002-8650-6058 P. Chang  0000-0003-4064-388X R. Cheaib
0000-0001-5729-8926 P. Cheema  0000-0001-8472-5727 C. Chen ,
0000-0003-1589-9955 B. G. Cheon ,  0000-0002-8803-4429 K. Chilikin
0000-0001-7620-2053 K. Chirapatpimol  0000-0003-2099-7760 H.-E. Cho
0000-0002-7008-3759 K. Cho  0000-0003-1705-7399 S.-J. Cho  0000-0002-1673-5664
S.-K. Choi  0000-0003-2747-8277 S. Choudhury  0000-0001-9841-0216 J. Cochran ,
0000-0002-1492-914X L. Corona  0000-0002-2577-9909 J. X. Cui ,
0000-0002-2398-3754 S. Das ,  0000-0001-6857-966X F. Dattola
0000-0003-3316-8574 E. De La Cruz-Burelo  0000-0002-7469-6974 S. A. De La
Motte  0000-0003-3905-6805 G. De Nardo ,  0000-0002-2047-9675 M. De Nuccio
0000-0002-0972-9047 G. De Pietro  0000-0001-8442-107X R. de Sangro
0000-0002-3808-5455 M. Destefanis ,  0000-0003-1997-6751 S. Dey
0000-0003-2997-3829 R. Dhamija ,  0000-0001-7052-3163 A. Di Canto
0000-0003-1233-3876 F. Di Capua  0000-0001-9076-5936 J. Dingfelder
0000-0001-5767-2121 Z. Doležal  0000-0002-5662-3675 I. Domínguez Jiménez
0000-0001-6831-3159 T. V. Dong ,  0000-0003-3043-1939 M. Dorigo
0000-0002-0681-6946 D. Dorner  0000-0003-3628-9267 K. Dort
0000-0003-0849-8774 D. Dossett  0000-0002-5670-5582 S. Dreyer
0000-0002-6295-100X S. Dubey  0000-0002-1345-0970 K. Dugic ,
0009-0006-6056-546X G. Dujany  0000-0002-1345-8163 P. Ecker
0000-0002-6817-6868 M. Eliachevitch  0000-0003-2033-537X P. Feichtinger ,
0000-0003-3966-7497 T. Ferber  0000-0002-6849-0427 D. Ferlewicz
0000-0002-4374-1234 T. Fillinger  0000-0001-9795-7412 C. Finck
0000-0002-5068-5453 G. Finocchiaro  0000-0002-3936-2151 A. Fodor ,
0000-0002-2821-759X F. Forti  0000-0001-6535-7965 A. Frey  0000-0001-7470-3874
B. G. Fulsom ,  0000-0002-5862-9739 A. Gabrielli  0000-0001-7695-0537 E.
Ganiev ,  0000-0001-8346-8597 M. Garcia-Hernandez  0000-0003-2393-3367 R. Garg
0000-0002-7406-4707 G. Gaudino ,  0000-0001-5983-1552 V. Gaur
0000-0002-8880-6134 A. Gaz  0000-0001-6754-3315 A. Gellrich ,
0000-0003-0974-6231 G. Ghevondyan  0000-0003-0096-3555 D. Ghosh
0000-0002-3458-9824 H. Ghumaryan  0000-0001-6775-8893 G. Giakoustidis
0000-0001-5982-1784 R. Giordano  0000-0002-5496-7247 A. Giri
0000-0002-8895-0128 A. Glazov  0000-0002-8553-7338 B. Gobbo
0000-0002-3147-4562 R. Godang  0000-0002-8317-0579 O. Gogota
0000-0003-4108-7256 P. Goldenzweig  0000-0001-8785-847X W. Gradl ,
0000-0002-9974-8320 T. Grammatico ,  0000-0002-2818-9744 S. Granderath
0000-0002-9945-463X E. Graziani  0000-0001-8602-5652 D. Greenwald
0000-0001-6964-8399 Z. Gruberová  0000-0002-5691-1044 T. Gu
0000-0002-1470-6536 Y. Guan  0000-0002-5541-2278 K. Gudkova
0000-0002-5858-3187 S. Halder ,  0000-0002-6280-494X Y. Han
0000-0001-6775-5932 T. Hara ,  0000-0002-4321-0417 C. Harris
0000-0003-0448-4244 K. Hayasaka ,  0000-0002-6347-433X H. Hayashii
0000-0002-5138-5903 S. Hazra  0000-0001-6954-9593 C. Hearty
0000-0001-6568-0252 M. T. Hedges  0000-0001-6504-1872 A. Heidelbach
0000-0002-6663-5469 I. Heredia de la Cruz ,  0000-0002-8133-6467 M. Hernández
Villanueva  0000-0002-6322-5587 T. Higuchi ,  0000-0002-7761-3505 M. Hoek ,
0000-0002-1893-8764 M. Hohmann  0000-0001-5147-4781 P. Horak
0000-0001-9979-6501 C.-L. Hsu ,  0000-0002-1641-430X T. Humair ,
0000-0002-2922-9779 T. Iijima  0000-0002-4271-711X K. Inami
0000-0003-2765-7072 G. Inguglia  0000-0003-0331-8279 N. Ipsita
0000-0002-2927-3366 A. Ishikawa ,  0000-0002-3561-5633 R. Itoh ,
0000-0003-1590-0266 M. Iwasaki  0000-0002-9402-7559 W. W. Jacobs ,
0000-0002-9996-6336 D. E. Jaffe ,  0000-0003-3122-4384 E.-J. Jang
0000-0002-1935-9887 Q. P. Ji  0000-0003-2963-2565 S. Jia  0000-0001-8176-8545
Y. Jin  0000-0002-7323-0830 H. Junkerkalefeld ,  0000-0003-3987-9895 M. Kaleta
,  0000-0002-2863-5476 D. Kalita  0000-0003-3054-1222 A. B. Kaliyar
0000-0002-2211-619X J. Kandra  0000-0001-5635-1000 K. H. Kang
0000-0002-6816-0751 S. Kang  0000-0002-5320-7043 G. Karyan ,
0000-0001-5365-3716 T. Kawasaki ,  0000-0002-4089-5238 F. Keil
0000-0002-7278-2860 C. Kiesling ,  0000-0002-2209-535X C.-H. Kim ,
0000-0002-5743-7698 D. Y. Kim  0000-0001-8125-9070 K.-H. Kim ,
0000-0002-4659-1112 Y.-K. Kim ,  0000-0002-9695-8103 H. Kindo ,
0000-0002-6756-3591 K. Kinoshita  0000-0001-7175-4182 P. Kodyš ,
0000-0002-8644-2349 T. Koga  0000-0002-1644-2001 S. Kohani
0000-0003-3869-6552 K. Kojima  0000-0002-3638-0266 T. Konno
0000-0003-2487-8080 A. Korobov ,  0000-0001-5959-8172 S. Korpar
0000-0003-0971-0968 E. Kovalenko ,  0000-0001-8084-1931 R. Kowalewski
0000-0002-7314-0990 T. M. G. Kraetzschmar  0000-0001-8395-2928 P. Križan
0000-0002-4967-7675 P. Krokovny ,  0000-0002-1236-4667 T. Kuhr ,
0000-0001-6251-8049 Y. Kulii  0000-0001-6217-5162 J. Kumar
0000-0002-8465-433X M. Kumar  0000-0002-6627-9708 R. Kumar
0000-0002-6277-2626 K. Kumara  0000-0003-1572-5365 T. Kunigo ,
0000-0001-9613-2849 A. Kuzmin ,  0000-0002-7011-5044 Y.-J. Kwon ,
0000-0001-9448-5691 S. Lacaprara  0000-0002-0551-7696 K. Lalwani ,
0000-0002-7294-396X T. Lam  0000-0001-9128-6806 L. Lanceri
0000-0001-8220-3095 J. S. Lange  0000-0003-0234-0474 M. Laurenza
0000-0002-7400-6013 K. Lautenbach  0000-0003-3762-694X R. Leboucher ,
0000-0003-3097-6613 F. R. Le Diberder  0000-0002-9073-5689 M. J. Lee
0000-0003-4528-4601 P. Leo ,  0000-0003-3833-2900 C. Lemettais
0009-0008-5394-5100 D. Levit  0000-0001-5789-6205 P. M. Lewis
0000-0002-5991-622X L. K. Li ,  0000-0002-7366-1307 S. X. Li
0000-0003-4669-1495 Y. Li  0000-0002-4413-6247 Y. B. Li  0000-0002-9909-2851
J. Libby  0000-0002-1219-3247 M. H. Liu ,  0000-0002-9376-1487 Q. Y. Liu
0000-0002-7684-0415 Y. Liu  0000-0002-8374-3947 Z. Q. Liu ,
0000-0002-0290-3022 D. Liventsev  0000-0003-3416-0056 S. Longo
0000-0002-8124-8969 T. Lueck ,  0000-0003-3915-2506 C. Lyu ,
0000-0002-2275-0473 Y. Ma  0000-0001-8412-8308 M. Maggiora
0000-0003-4143-9127 S. P. Maharana  0000-0002-1746-4683 R. Maiti
0000-0001-5534-7149 S. Maity  0000-0003-3076-9243 G. Mancinelli
0000-0003-1144-3678 R. Manfredi  0000-0002-8552-6276 E. Manoni
0000-0002-9826-7947 M. Mantovano ,  0000-0002-5979-5050 D. Marcantonio
0000-0002-1315-8646 S. Marcello  0000-0003-4144-863X C. Marinas
0000-0003-1903-3251 C. Martellini ,  0000-0002-7189-8343 A. Martini
0000-0003-1161-4983 T. Martinov  0000-0001-7846-1913 L. Massaccesi
0000-0003-1762-4699 M. Masuda  0000-0002-7109-5583 K. Matsuoka ,
0000-0003-1706-9365 D. Matvienko ,  0000-0002-2698-5448 S. K. Maurya
0000-0002-7764-5777 J. A. McKenna  0000-0001-9871-9002 R. Mehta ,
0000-0001-8670-3409 F. Meier  0000-0002-6088-0412 M. Merola
0000-0002-7082-8108 F. Metzner  0000-0002-0128-264X C. Miller ,
0000-0003-2631-1790 M. Mirra  0000-0002-1190-2961 S. Mitra
0000-0002-1118-6344 K. Miyabayashi  0000-0003-4352-734X G. B. Mohanty ,
0000-0001-6850-7666 S. Mondal ,  0000-0002-3054-8400 S. Moneta
0000-0003-2184-7510 H.-G. Moser ,  0000-0003-3579-9951 M. Mrvar
0000-0001-6388-3005 R. Mussa ,  0000-0002-0294-9071 I. Nakamura
0000-0002-7640-5456 K. R. Nakamura ,  0000-0001-7012-7355 M. Nakao ,
0000-0001-8424-7075 Y. Nakazawa ,  0000-0002-6271-5808 A. Narimani Charan
0000-0002-5975-550X M. Naruki  0000-0003-1773-2999 D. Narwal
0000-0001-6585-7767 Z. Natkaniec  0000-0003-0486-9291 A. Natochii
0000-0002-1076-814X L. Nayak  0000-0002-7739-914X M. Nayak
0000-0002-2572-4692 G. Nazaryan  0000-0002-9434-6197 M. Neu
0000-0002-4564-8009 C. Niebuhr  0000-0002-4375-9741 J. Ninkovic ,
0000-0003-1523-3635 S. Nishida ,  0000-0001-6373-2346 S. Ogawa ,
0000-0002-7310-5079 Y. Onishchuk ,  0000-0002-8261-7543 H. Ono
0000-0003-4486-0064 F. Otani ,  0000-0001-6016-219X P. Pakhlov ,
0000-0001-7426-4824 G. Pakhlova  0000-0001-7518-3022 S. Pardi ,
0000-0001-7994-0537 K. Parham  0000-0001-9556-2433 H. Park
0000-0001-6087-2052 J. Park  0000-0001-6520-0028 S.-H. Park
0000-0001-6019-6218 B. Paschen  0000-0003-1546-4548 A. Passeri
0000-0003-4864-3411 S. Patra  0000-0002-4114-1091 S. Paul  0000-0002-8813-0437
T. K. Pedlar  0000-0001-9839-7373 R. Peschke ,  0000-0002-2529-8515 R.
Pestotnik  0000-0003-1804-9470 M. Piccolo ,  0000-0001-9750-0551 L. E.
Piilonen  0000-0001-6836-0748 G. Pinna Angioni  0000-0003-0808-8281 P. L. M.
Podesta-Lerma  0000-0002-8152-9605 T. Podobnik  0000-0002-6131-819X S.
Pokharel  0000-0002-3367-738X C. Praz ,  0000-0002-6154-885X S. Prell
0000-0002-0195-8005 E. Prencipe  0000-0002-9465-2493 M. T. Prim
0000-0002-1407-7450 I. Prudiev  0000-0002-0819-284X H. Purwar ,
0000-0002-3876-7069 P. Rados ,  0000-0003-0690-8100 G. Raeuber
0000-0003-2948-5155 S. Raiz  0000-0001-7010-8066 N. Rauls  0000-0002-6583-4888
M. Reif ,  0000-0002-0706-0247 S. Reiter ,  0000-0002-6542-9954 L. Reuter ,
0000-0002-5930-6237 I. Ripp-Baudot  0000-0002-1897-8272 G. Rizzo ,
0000-0003-1788-2866 S. H. Robertson ,  0000-0003-4096-8393 M. Roehrken ,
0000-0003-0654-2866 J. M. Roney  0000-0001-7802-4617 A. Rostomyan ,
0000-0003-1839-8152 N. Rout  0000-0002-4310-3638 G. Russo ,
0000-0001-5823-4393 D. A. Sanders ,  0000-0002-4902-966X S. Sandilya
0000-0002-4199-4369 L. Santelj ,  0000-0003-3904-2956 Y. Sato ,
0000-0003-3751-2803 V. Savinov ,  0000-0002-9184-2830 B. Scavino
0000-0003-1771-9161 S. Schneider ,  0009-0002-5899-0353 C. Schwanda ,
0000-0003-4844-5028 M. Schwickardi ,  0000-0003-2033-6700 Y. Seino
0000-0002-8378-4255 A. Selce  0000-0001-8228-9781 K. Senyo
0000-0002-1615-9118 J. Serrano  0000-0003-2489-7812 M. E. Sevior
0000-0002-4824-101X C. Sfienti  0000-0002-5921-8819 W. Shan
0000-0003-2811-2218 C. Sharma  0000-0002-1312-0429 C. P. Shen ,
0000-0002-9012-4618 X. D. Shi  0000-0002-7006-6107 T. Shillington ,
0000-0003-3862-4380 T. Shimasaki  0000-0003-3291-9532 J.-G. Shiu ,
0000-0002-8478-5639 D. Shtol  0000-0002-0622-6065 A. Sibidanov ,
0000-0001-8805-4895 F. Simon  0000-0002-5978-0289 J. B. Singh
0000-0001-9029-2462 J. Skorupa  0000-0002-8566-621X K. Smith
0000-0003-0446-9474 R. J. Sobie  0000-0001-7430-7599 M. Sobotzik
0000-0002-1773-5455 A. Soffer  0000-0002-0749-2146 A. Sokolov
0000-0002-9420-0091 E. Solovieva ,  0000-0002-5735-4059 S. Spataro
0000-0001-9601-405X B. Spruck  0000-0002-3060-2729 M. Starič
0000-0001-8751-5944 P. Stavroulakis  0000-0001-9914-7261 S. Stefkova
0000-0003-2628-530X R. Stroili ,  0000-0002-3453-142X Y. Sue ,
0000-0003-2430-8707 M. Sumihama ,  0000-0002-8954-0585 K. Sumisawa ,
0000-0001-7003-7210 W. Sutcliffe  0000-0002-9795-3582 N. Suwonjandee
0009-0000-2819-5020 H. Svidras ,  0000-0003-4198-2517 M. Takahashi ,
0000-0003-1171-5960 M. Takizawa  0000-0001-8225-3973 U. Tamponi
0000-0001-6651-0706 S. Tanaka  0000-0002-6029-6216 K. Tanida
0000-0002-8255-3746 F. Tenchini ,  0000-0003-3469-9377 A. Thaller
0000-0003-4171-6219 O. Tittel  0000-0001-9128-6240 R. Tiwary
0000-0002-5887-1883 D. Tonelli  0000-0002-1494-7882 E. Torassa
0000-0003-2321-0599 K. Trabelsi ,  0000-0001-6567-3036 I. Tsaklidis
0000-0003-3584-4484 I. Ueda ,  0000-0002-6833-4344 T. Uglov ,
0000-0002-4944-1830 K. Unger  0000-0001-7378-6671 Y. Unno  0000-0003-3355-765X
K. Uno  0000-0002-2209-8198 S. Uno  0000-0002-3401-0480 P. Urquijo
0000-0002-0887-7953 Y. Ushiroda  0000-0003-3174-403X S. E. Vahsen ,
0000-0003-1685-9824 R. van Tonder  0000-0002-7448-4816 K. E. Varvell ,
0000-0003-1017-1295 M. Veronesi  0000-0002-1916-3884 A. Vinokurova
0000-0003-4220-8056 V. S. Vismaya  0000-0002-1606-5349 L. Vitale
0000-0003-3354-2300 V. Vobbilisetti  0000-0002-4399-5082 R. Volpe
0000-0003-1782-2978 A. Vossen ,  0000-0003-0983-4936 B. Wach
0000-0003-3533-7669 M. Wakai ,  0000-0003-2818-3155 S. Wallner ,
0000-0002-9105-1625 E. Wang ,  0000-0001-6391-5118 M.-Z. Wang
0000-0002-0979-8341 X. L. Wang  0000-0001-5805-1255 Z. Wang
0000-0002-3536-4950 A. Warburton  0000-0002-2298-7315 M. Watanabe
0000-0001-6917-6694 S. Watanuki  0000-0002-5241-6628 C. Wessel ,
0000-0003-0959-4784 X. P. Xu ,  0000-0001-5096-1182 B. D. Yabsley
0000-0002-2680-0474 S. Yamada  0000-0002-8858-9336 W. Yan  0000-0003-0713-0871
S. B. Yang ,  0000-0002-9543-7971 J. Yelton  0000-0001-8840-3346 J. H. Yin
0000-0002-1479-9349 Y. M. Yook  0000-0002-4912-048X K. Yoshihara
0000-0002-3656-2326 C. Z. Yuan  0000-0002-1652-6686 L. Zani ,
0000-0003-4957-805X F. Zeng  0009-0003-6474-3508 B. Zhang  0000-0002-5065-8762
Y. Zhang ,  0000-0003-2961-2820 V. Zhilich ,  0000-0002-0907-5565 Q. D. Zhou ,
0000-0001-5968-6359 X. Y. Zhou  0000-0002-0299-4657 V. I. Zhukova
0000-0002-8253-641X R. Žlebčík  0000-0003-1644-8523 coll-
<EMAIL_ADDRESS>
###### Abstract

We present the result of a search for the charged-lepton-flavor violating
decay $\tau^{-}\rightarrow\mu^{-}\mu^{+}\mu^{-}$ using a 424$\mbox{\,fb}^{-1}$
sample of data recorded by the Belle II experiment at the SuperKEKB
$e^{+}e^{-}$ collider. The selection of
$e^{+}e^{-}\rightarrow\tau^{+}\tau^{-}$ events is based on an inclusive
reconstruction of the non-signal tau decay, and on a boosted decision tree to
suppress background. We observe one signal candidate, which is compatible with
the expectation from background processes. We set a $90\%$ confidence level
upper limit of $1.9\times 10^{-8}$ on the branching fraction of the
$\tau^{-}\rightarrow\mu^{-}\mu^{+}\mu^{-}$ decay, which is the most stringent
bound to date.

##  1 Introduction

In the standard model (SM) with massless neutrino hypotheses, the charged-
lepton flavor is accidentally conserved. However, this symmetry is broken at
loop-level when taking into account neutrino mixing, which implies the
existence of charged-lepton-flavor violation, and thus the existence of
processes such as $\mu\rightarrow e$, $\tau\rightarrow e$ and
$\tau\rightarrow\mu$ conversions. In the simplest SM extension that allows for
massive neutrinos, all charged-lepton-flavor-violating (LFV) rates are
proportional to the square of neutrino masses. This results in predicted decay
rates of $10^{-50}$ ref:SM_cLFV ; ref:tau_lfv ; Blackstone:2019njl , well
below the sensitivities of any experiment. The observation of LFV decays
would, therefore, provide indisputable evidence of physics beyond the SM. Over
the past four decades, the CLEO experiment at CESR, and the first generation
$B$-factory experiments BaBar at SLAC and Belle at KEK, have searched for LFV
in $\tau$-lepton decays. In total, 52 LFV $\tau$ decays with neutrinoless two-
body or three-body final states have been investigated. Among these,
$\tau^{-}\rightarrow\ell^{-}\ell^{+\prime}\ell^{-\prime\prime}$
decays,111Charge conjugation is implied throughout this paper. where
$\ell^{(\prime,\prime\prime)}=e,\mu$, and in particular
$\tau^{-}\rightarrow\mu^{-}\mu^{+}\mu^{-}$, have garnered significant
attention in recent years. This is due to the potential enhancement of the
branching fraction up to a value of $10^{-8}$ in scenarios beyond the SM
Abada:2021zcm ; Raidal:2008jk ; Teixeira:2016ecr ; Bhattacharya:2016mcc ;
Feruglio:2017rjo ; Greljo:2015mma . The most stringent upper limit for this
decay was set by the Belle collaboration at the level of $2.1\times 10^{-8}$
at the 90% confidence level (C.L.) using electron-positron data corresponding
to an integrated luminosity of 782$\mbox{\,fb}^{-1}$ Hayasaka:2010np . The
upper limits on the $\tau^{-}\rightarrow\mu^{-}\mu^{+}\mu^{-}$ branching
fraction set by the CLEO CLEO:1997aqy  and BaBar BaBar:2010axs  collaborations
are 1.9$\times 10^{-6}$ and $3.3\times 10^{-8}$, using 4.79$\mbox{\,fb}^{-1}$
and 468$\mbox{\,fb}^{-1}$ of data, respectively. In addition, experiments at
the Large Hadron Collider have also contributed to this search. The LHCb
LHCb:2014kws , ATLAS ATLAS:2016jts , and CMS CMS:2020kwy  collaborations
reported upper limits at the 90% C.L. of $4.6\times 10^{-8}$, $3.8\times
10^{-7}$ and $8.0\times 10^{-8}$ using proton-proton collision data
corresponding to integrated luminosities of 3$\mbox{\,fb}^{-1}$,
20.3$\mbox{\,fb}^{-1}$, and 33.2$\mbox{\,fb}^{-1}$ respectively. Recently, CMS
presented a new result incorporating 2017 and 2018 data, amounting to a total
integrated luminosity of 131$\mbox{\,fb}^{-1}$, placing an upper limit at 2.9
$\times 10^{-8}$ at 90% C.L CMS:2023iqy . We report the results of a search
for the LFV $\tau^{-}\rightarrow\mu^{-}\mu^{+}\mu^{-}$ decay using an untagged
selection, in which, in contrast to the previous searches, only the signal tau
decay mode is explicitly reconstructed. We use a sample of 389 million
$e^{+}e^{-}\rightarrow\tau^{+}\tau^{-}$ events recorded with the Belle II
detector Abe:2010gxa  at the asymmetric-energy $e^{+}e^{-}$ SuperKEKB collider
Akai:2018mbz . The data, collected between 2019 and 2022, corresponds to an
integrated luminosity of 424$\mbox{\,fb}^{-1}$. Section 2 gives an overview of
the Belle II detector and the data samples used. Section 3 presents the
overall strategy for this analysis and the details of the candidate
reconstruction and selection. Section 4 discusses the systematic
uncertainties, and Section 5 presents the branching fraction measurement and
limit computation. Finally, a summary is given in Section 6.

##  2 The Belle II detector, simulation and data samples

The Belle II detector consists of several subdetectors arranged in a
cylindrical structure around the $e^{+}e^{-}$ interaction point Abe:2010gxa .
Charged-particle trajectories (tracks) are reconstructed using a two-layer
silicon-pixel detector, surrounded by a four-layer double-sided silicon-strip
detector and a central drift chamber (CDC). Only 15% of the second pixel layer
was installed when the data were collected. Outside the CDC, which also
provides $d{E}$/$d{x}$ energy-loss measurements, particle identification using
Cherenkov radiation is provided by a time-of-propagation (TOP) detector and an
aerogel ring-imaging Cherenkov (ARICH) detector which cover the barrel and
forward endcap regions, respectively. An electromagnetic calorimeter (ECL),
divided into the forward endcap, barrel, and backward endcap regions, fills
the remaining volume inside a 1.5 T superconducting solenoid and is used to
reconstruct photons and electrons. A $K_{L}^{0}$ and muon detection system
(KLM) is installed in the iron flux return of the solenoid. The $z$ axis of
the laboratory frame is defined as the detector solenoid axis, with the
positive direction along the electron beam. The polar angle $\theta$ and the
transverse plane are defined relative to this axis. The search presented here
is based on $e^{+}e^{-}$ collisions at the center of mass energy of 10.58 GeV
corresponding to the mass of the $\mathchar 28935\relax{(4S)}$
(362$\mbox{\,fb}^{-1}$), 60 $\mathrm{\,Me\kern-1.00006ptV}$ below it
(42$\mbox{\,fb}^{-1}$), and around 10.75 GeV (19$\mbox{\,fb}^{-1}$). The
corresponding cross sections for $e^{+}e^{-}\rightarrow\tau^{+}\tau^{-}$
production are $\sigma_{\tau\tau}$=0.919 $\rm\,nb$ at the $\mathchar
28935\relax{(4S)}$ energy, 0.929 $\rm\,nb$ for off-resonance data, and 0.891
$\rm\,nb$ around the $\mathchar 28935\relax{(5S)}$ energy, leading to a data
sample of 389 million $\tau$-pairs Banerjee:2007is . Monte-Carlo simulated
events are used to optimize the selection and background rejection and to
measure the signal efficiency. To study the signal process, we use 10 million
$e^{+}e^{-}\rightarrow\tau^{+}\tau^{-}$ events, in which one $\tau$ decays to
three muons and the other to SM-allowed decays. The potential background
processes studied using simulation include $e^{+}e^{-}\rightarrow
q\overline{}q$ events, where $q$ indicates a $u$, $d$, $c$, or $s$ quark;
$e^{+}e^{-}\rightarrow b\bar{b}$ events;
$e^{+}e^{-}\rightarrow\ell^{+}\ell^{-}(\gamma)$, where $\ell=e,\mu$;
$e^{+}e^{-}\rightarrow e^{+}e^{-}h^{+}h^{-}$ events, where $h$ indicates a
pion, kaon, or proton; and four-lepton processes: $e^{+}e^{-}\rightarrow
e^{+}e^{-}e^{+}e^{-},\mu^{+}\mu^{-}\mu^{+}\mu^{-},\mu^{+}\mu^{-}e^{+}e^{-},e^{+}e^{-}\tau^{+}\tau^{-},\mu^{+}\mu^{-}\tau^{+}\tau^{-}$.
The $e^{+}e^{-}\rightarrow\tau^{+}\tau^{-}$ process is generated using the
KKMC generator Jadach:1999vf . The $\tau$ decays are simulated by the TAUOLA
generator Jadach:1990mz  and their final state radiated photons by PHOTOS
Barberio:1990ms . We use KKMC to simulate $\mu^{+}\mu^{-}(\gamma)$ and
$q\overline{}q$ production; PYTHIA Sjostrand:2014zea  for the fragmentation of
the $q\overline{}q$ pair; PYTHIA interfaced with EvtGen Lange:2001uf  for the
production and decay of $e^{+}e^{-}\rightarrow b\bar{b}$ events; BabaYaga@NLO
Balossini:2006wc ; Balossini:2008xr ; CarloniCalame:2003yt ;
CarloniCalame:2001ny ; CarloniCalame:2000pz  for $e^{+}e^{-}\rightarrow
e^{+}e^{-}(\gamma)$ events; and AAFH BERENDS1985421 ; BERENDS1985441 ;
BERENDS1986285  and TREPS Uehara:1996bgt  for the production of non-radiative
four-leptons and $e^{+}e^{-}h^{+}h^{-}$ final states. The size of the
simulated samples for $e^{+}e^{-}\rightarrow\tau^{+}\tau^{-}$ and
$e^{+}e^{-}\rightarrow q\overline{}q$ events is equivalent to an integrated
luminosity of 8 ab-1, while it ranges between 100$\mbox{\,fb}^{-1}$ and 2 ab-1
for the other processes. The online event selection (hardware trigger) is
based on the energy deposits (clusters) and their topologies in the ECL, or on
independent trigger selection based on the number of charged particles
reconstructed in the CDC. Most of the events are selected requiring a ECL
total energy larger than 1 GeV and a topology incompatible with Bhabha events.
The Belle II analysis software Kuhr:2018lps ; basf2-zenodo  uses the GEANT4
Agostinelli:2002hh  package to simulate the response of the detector to the
passage of the particles and also provides a simulation of the triggers.

##  3 Event reconstruction and background rejection

###  3.1 Overall strategy

In the $e^{+}e^{-}$ center-of-mass (c.m.) frame, the $\tau$ leptons are
produced in opposite directions, with the decay products of one $\tau$
isolated from those of the other $\tau$ and contained in opposite hemispheres.
The boundary between the hemispheres is experimentally defined by the plane
perpendicular to the vector $\mathbf{\hat{n}}_{T}$ that maximizes the thrust
value ($T$):

$T=\max_{\mathbf{\hat{n}}_{T}}\left(\dfrac{\sum_{i}\left|\mathbf{p^{*}}_{i}\cdot\mathbf{\hat{n}}_{T}\right|}{\sum_{i}\left|\mathbf{p^{*}}_{i}\right|}\right),$
(1)

where $\mathbf{p}_{i}^{*}$ is the momentum of final state particle $i$ in the
$e^{+}e^{-}$ c.m. frame Brandt:1964sa ; Farhi:1977sg , including both charged
and neutral particles. Here, and throughout the paper, quantities in the
$e^{+}e^{-}$ c.m. frame are indicated by an asterisk. We define the signal
hemisphere as the one that contains the
$\tau^{-}\rightarrow\mu^{-}\mu^{+}\mu^{-}$ decay candidate, reconstructed by
combining three charged particles. The trajectories of those particles are
required to be displaced from the average interaction point by less than 3
$\rm\,cm$ along the $z$ axis and less than 1 $\rm\,cm$ in the transverse
plane. Identification of muons of momentum greater than 0.7
${\mathrm{\,Ge\kern-1.00006ptV\\!/}c}$ relies mostly on their penetration
depth in the KLM whereas for lower momenta information from the CDC and ECL
dominate. Muons are identified using the discriminator
$\mathcal{P}_{\mu}={\cal L}_{\mu}/({\cal L}_{e}+{\cal L}_{\mu}+{\cal
L}_{\pi}+{\cal L}_{K}+{\cal L}_{p}+{\cal L}_{d})$ where the likelihoods ${\cal
L}_{i}$ for each charged-particle hypothesis
$i=e,\mu,\pi,K,\text{proton}(p),\text{deuteron}(d)$ combine particle-
identification information from CDC, TOP, ARICH, ECL, and KLM. The muon
identification efficiency for the requirement $\mathcal{P}_{\mu}>0.5(0.95)$ is
92(89)% with a pion misidentification probability of 5(3)%, for particle
momenta between 1.0 and 1.5${\mathrm{\,Ge\kern-1.00006ptV\\!/}c}$. Since the
decay $\tau^{-}\rightarrow\mu^{-}\mu^{+}\mu^{-}$ is a neutrinoless process,
the invariant mass $M_{3\mu}$ of the reconstructed muons should be consistent
with the mass of the $\tau$ lepton, except for decays affected by final state
radiation (FSR) from the $\tau$ or its decay products. The energy
$E_{\tau}^{*}$ in the c.m. system should be half of the $e^{+}e^{-}$ c.m.
energy $\sqrt{s}/2$ except for corrections from initial state radiation (ISR)
from the $e^{\pm}$ beams and FSR. Thus, the energy difference $\Delta
E_{3\mu}=E_{\tau}^{*}-\sqrt{s}/2$ should be close to zero. These
characteristic features are used to define the signal region, which is hidden
until the finalization of the selection procedure to avoid any experimenter’s
bias, and to optimize the selection criteria, which maximize the signal
efficiency and suppress the contribution of background events. The
distribution of signal in the ($M_{3\mu}$, $\Delta E_{3\mu}$) plane is
broadened by detector resolution and radiative effects. The radiation of
photons from the initial state leads to a tail at low values of $\Delta
E_{3\mu}$. Instead, FSR produces a tail at high values of $M_{3\mu}$ and
$\Delta E_{3\mu}$. Previous $B$-factory experiments relied on the
reconstruction of $e^{+}e^{-}\rightarrow\tau^{+}\tau^{-}$ events where one
$\tau$ decays to $\mu^{-}\mu^{+}\mu^{-}$ and the other $\tau$, named the _tag_
, decays into a final state with a single charged particle. The background
rejection strategy relied on a cut-based selection. Here, we expand the search
by using inclusive tagging to include a wider range of decays for the tag
$\tau$ lepton. We also optimize the signal efficiency and reduce the
background by using boosted decision trees (BDT). To compare the sensitivity
of the inclusive tagging to the traditional one-prong tagging, we also perform
a validation using a reconstruction method similar to the one used previously
by Belle Hayasaka:2010np  and BaBar BaBar:2010axs .

###  3.2 Inclusive tagging selection

Signal candidate $\tau^{-}\rightarrow\mu^{-}\mu^{+}\mu^{-}$ decays are
obtained by combining three muons originating from the interaction point with
a total charge equal to $\pm 1$, belonging to the same hemisphere and having
$\mathcal{P}_{\mu}>0.5$. Instead of explicitly reconstructing the tagged
$\tau$ lepton, we use the inclusive properties of all other particles in the
event. Photons are reconstructed from ECL clusters within the CDC acceptance
($17^{\circ}<\theta<150^{\circ}$) and not associated with any tracks. Photons
used for $\pi^{0}$ reconstruction must have an energy deposit of at least 0.1
$\mathrm{\,Ge\kern-1.00006ptV}$. Neutral pions are then identified as photon
pairs with invariant masses within $0.115<M_{\gamma\gamma}<0.152$
${\mathrm{\,Ge\kern-1.00006ptV\\!/}c^{2}}$, which corresponds to a range of
approximately $\pm 2.5$ units of resolution about the known $\pi^{0}$ mass
ParticleDataGroup:2022pth . All photons with energies greater than 0.2
$\mathrm{\,Ge\kern-1.00006ptV}$, along with photons part of reconstructed
$\pi^{0}$ candidates, are used to define variables related to the kinematic
properties of the event such as the missing momentum, defined as the
difference between the momenta of the initial $e^{+}e^{-}$ and that of all
reconstructed tracks and photons in the event, its mass and the thrust axis.
This high energy threshold for photon candidates aims to reduce photons from
beam background. All tracks and clusters that are not used in the signal
reconstruction form the rest of the event (ROE), whose kinematic properties
are exploited for further background suppression. In order to suppress beam-
related backgrounds, the following criteria must be satisfied for particles to
be retained as part of the ROE: photons are required to have energies of at
least 0.2 $\mathrm{\,Ge\kern-1.00006ptV}$, and tracks must have transverse
momenta higher than 0.075 ${\mathrm{\,Ge\kern-1.00006ptV\\!/}c}$, and be
displaced from the average interaction point by less than 3 $\rm\,cm$ along
the $z$ axis and less than 1 $\rm\,cm$ in the transverse plane. In addition,
tracks and photons should lie within the CDC angular acceptance. We define
rectangular regions in the two-dimensional plane ($M_{3\mu}$, $\Delta
E_{3\mu}$) centered at the expected position of the signal peak ($\bar{\mu}$),
with side lengths proportional to the expected resolutions ($\delta$). For
each variable, $\delta$ values are estimated as the widths of a bifurcated
Gaussian function fitted to the simulated signal distribution. The expected
peak position $\bar{\mu}$, and the high- and low-mass resolutions
$\delta_{high(low)}$, are reported in Table 1.

Table 1: Fitted central values and resolutions for $M_{3\mu}$ and $\Delta E_{3\mu}$. Variable |
$\bar{\mu}$ | $\delta_{low}$ | $\delta_{high}$
---|---|---|---
$M_{3\mu}$ (${\mathrm{\,Me\kern-1.00006ptV\\!/}c^{2}}$) | 1777.35$\pm$0.07 | 4.80$\pm$0.07 | 4.44$\pm$0.06
$\Delta E_{3\mu}$ ($\mathrm{\,Me\kern-1.00006ptV}$) | 0.7$\pm$0.3 | 14.9$\pm$0.3 | 10.0$\pm$0.5

A $\pm 20\,\delta$-wide rectangular region in the ($M_{3\mu}$, $\Delta
E_{3\mu}$) plane is used for the background-rejection optimization using
simulated data. The final yield extraction is performed in a 5 $\delta$ semi-
axis wide, asymmetric elliptical signal region (SR). The rotation angle of the
axis of the ellipse is obtained from a linear fit to the profile distribution
of $M_{3\mu}$ versus $\Delta E_{3\mu}$. The sideband region (SB), defined as
the area covering $\pm 20\,\delta_{M}$, $\pm 10\,\delta_{\Delta E}$ subtracted
of the SR, is used for the validation of the background rejection. The
background originates from radiative dilepton and four-lepton final states
(low-multiplicity backgrounds) with potential electrons misidentified as
muons, incorrectly reconstructed SM $e^{+}e^{-}\rightarrow\tau^{+}\tau^{-}$
events, and continuum hadronization processes from $e^{+}e^{-}\rightarrow
q\overline{}q$ events, where pions are misidentified as muons. The other
simulated processes are found to be negligible. To suppress events with pions
or other particles misidentified as muons, two of the muon candidates are
required to have $\mathcal{P}_{\mu}>0.95$. Low-multiplicity backgrounds are
concentrated at very high thrust values and have a missing momentum pointing
at the boundaries of the polar angular acceptance in the c.m. frame. In Fig.
1, there are disagreements between data and simulation in these regions,
consistent with unsimulated low-multiplicity backgrounds such as four-muon
final state processes with initial and final state radiations. To suppress
such backgrounds, we place requirements on the polar angle of the missing
momentum, $0.3<\theta^{*}_{\rm miss}<2.7\rm\,rad$, and the thrust,
$0.89<T<0.97$. The overall signal efficiency after those requirements,
including event reconstruction and selection, is $30.7\%$, while analysis of
the simulated background events predicts about 300 events in the $\pm
20\,\delta$ region. They are mostly due to $e^{+}e^{-}\rightarrow
q\overline{}q$ processes in which one or two of the final state muons are
misidentified.

Figure 1: Comparison between sideband data (black points with error bars) and
simulation after the inclusive tagging reconstruction for the polar angle of
missing momentum (left) and thrust (right) distributions. The various
simulated background processes are shown as a stack of color-filled
histograms, with statistical uncertainties displayed as hatched areas, while
the signal is shown as a red histogram with an arbitrary scale. The vertical
dashed green lines indicate the boundaries of the selection criteria.

###  3.3 BDT-based background rejection

To suppress the remaining background events, a boosted decision tree
classifier (BDT) is trained using the XGBoost library XGBoostPaper . The BDT
uses 32 variables related to three distinct categories. The first category
consists of variables associated with the signal $\tau$, such as the muon
energies, polar angles and ordered transverse momenta in the c.m. frame; the
flight time of the $\tau$ divided by its uncertainty; the cosine of the angle
between the $\tau$ momentum and its direction reconstructed from the origin
and decay vertices; the angle between the $\tau$ momentum and the closest
track in the event (including final-state muons); and the $\chi^{2}$ of the
kinematic fit of the $\tau$ decay chain. The second category involves
variables related to the ROE properties such as its mass, defined as the mass
of the four-vector resulting from the sum of all reconstructed objects forming
the ROE, its $\Delta E$, computed as the difference between the total energy
of all ROE objects and half the center-of-mass energy; the numbers of muon,
pion and electron candidates in the ROE, which are identified according to the
highest $\mathcal{P}_{\mu,\pi,e}$ value; the product of the ROE charge (which
can be anything) and $\tau$ signal charge (which is equal to $\pm 1$); the
thrust and the following quantities related to the ROE thrust axis, where the
ROE thrust axis is computed in a manner similar to the event thrust axis, but
using only the ROE detector objects: the cosine of the angle between the ROE
thrust axis and the signal thrust axis, and the cosine of the angle between
the ROE thrust axis and the $z$ axis. The third category of input variables
comprises the thrust, the numbers of tracks and photons in the event, the
total photon energy in the c.m. frame, and variables related to the missing
momentum of the event: its transverse component in the c.m. frame, its polar
angle in the c.m. frame, the cosine of the angle between the missing momentum
and each of the three muons, and the invariant missing mass squared. The BDT
is trained on 176000 simulated signal events and 3400 simulated background
$e^{+}e^{-}\rightarrow\tau^{+}\tau^{-}$ and $e^{+}e^{-}\rightarrow
q\overline{}q$ events in the $\pm 20\,\delta$ region, corresponding to a data
set equivalent to an integrated luminosity of 4$\mbox{\,ab}^{-1}$. Given that
the background retention is already low before the BDT training, a cross-
validation (_k-folding_) algorithm is used to reduce the impact of statistical
fluctuations on the training results 10.1145/307400.307439 . The BDT
parameters are optimized using the Optuna library akiba_optuna_2019 , which
minimizes the logarithmic loss function. The requirement on the BDT output
($p^{\rm BDT}$) is chosen to maximize the Punzi figure-of-merit Punzi:2003bu
defined as $\varepsilon_{3\mu}/(\alpha/2+\sqrt{B})$, where
$\varepsilon_{3\mu}$ is the signal efficiency, $B$ the number of background
events in the signal region and $\alpha=3$ (i.e. we optimise the search for a
three-standard-deviation signal). The size of the elliptical signal region is
also optimized, and the initial width of 5 $\delta$ is found to be the best
one, together with a requirement on the $p^{\rm BDT}$ to be larger than 0.9.
Finally, a data-driven requirement is introduced only to keep events in which
the sum of the charged-particles charges is zero. The total signal efficiency
evaluated on an independent simulation sample is $(20.42\pm 0.06)\%$, where
the uncertainty is due to the limited size of the sample, for a background
reduction of 98.2%. These values include corrections that account for the
differences between data and simulation in particle identification
efficiencies. After all background rejection steps, 88% of the selected signal
candidates in the simulation have one track in the ROE, while 12% contain
three tracks.

###  3.4 Simulation validation and expected background events

The agreement between data and simulation for the BDT output is checked using
the SB region as shown in Fig. 2. After the BDT selection, the number of
background events expected in the SB region from the simulation is
$2.0^{+0.7}_{-0.5}$, while in data we observe $3$.

Figure 2: Comparison between data (black points with error bars) and
simulation for the BDT output distribution for events in the sideband region
after the preselection. The various simulated background processes are shown
as a stack of color-filled histograms, with statistical uncertainties
displayed as hatched areas, while the signal is shown as a red histogram with
an arbitrary scale.

The number of expected events in the signal region after the BDT selection is
obtained from a data-driven method using three other regions of the plane
consisting of the BDT output and the distance to the signal peak: region $A$
is outside the SR with $0.2<p^{\rm BDT}<0.85$; region $B$ is inside the SR
with $0.2<p^{\rm BDT}<0.85$; and region $C$ is outside the SR with $p^{\rm
BDT}>0.9$. The distributions of events in data and simulation for those two
variables are shown in Fig. 3.

Figure 3: Comparison between data (black points with error bars) and
simulation for the BTD output (left) and distance to the signal peak (right).
The various simulated background processes are shown as a stack of color-
filled histograms, with statistical uncertainties displayed as hatched areas,
while the signal is shown as a red histogram with an arbitrary scale. The
green dashed line in the right-hand plot corresponds to the boundary of the
SR.

Events in region $B$ correspond to a subset of the originally hidden $\pm
5\,\delta$ SR. However, due to the BDT requirement, the signal efficiency in
this region is 50 times lower than in the search region, and the potential
signal is negligible compared with the background. The number of expected
events in the $\pm 5\,\delta$ SR region is given by

$N_{exp}=N_{C}\times R_{B/A},$ (2)

where $R_{B/A}=N_{B}/N_{A}$ is the transfer factor between the region $A$ and
$B$. The measured event yields are $N_{A}=4$, $N_{B}=2$, and $N_{C}=1$,
resulting in $N_{exp}=0.7_{-0.5}^{+0.6}$, obtained by assuming Poisson
distributions with means corresponding to the yields in regions A, B and C.

##  4 Systematic uncertainties

One category of systematic uncertainties arises from differences between
experimental data and simulation due to possible mismodeling in the generation
and reconstruction of the simulated samples, and affects the signal efficiency
uncertainty. We take into account the systematic uncertainty associated with
the corrections to the simulated muon-identification efficiencies, derived
from auxiliary measurements in data using
${J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}\rightarrow\mu^{+}\mu^{-}$,
$e^{+}e^{-}\rightarrow\mu^{+}\mu^{-}\gamma$, and $e^{+}e^{-}\rightarrow
e^{+}e^{-}\mu^{+}\mu^{-}$ events. These corrections are obtained as functions
of momentum, polar angle and charge, and applied to events reconstructed from
simulation. The systematic uncertainty is obtained by varying the corrections
within their statistical and systematic uncertainties and estimating the
impact of these variations on the selection efficiency. Adding the statistical
and systematic variations in quadrature, the result is a relative uncertainty
in the signal efficiency of ${}^{+2.4}_{-2.1}$%. The difference between data
and simulation in track-reconstruction efficiency is measured in
$e^{+}e^{-}\rightarrow\tau^{+}\tau^{-}$ events, selecting $\tau^{-}\rightarrow
e^{-}\overline{\nu}_{e}\nu_{\tau}$ and
$\tau^{-}\rightarrow\pi^{-}\pi^{+}\pi^{-}\nu_{\tau}$ decays. A discrepancy of
0.24% per track is observed, resulting in a systematic uncertainty of $1.0\%$.
The agreement between data and simulation for the trigger efficiency is
evaluated using the $\tau^{-}\rightarrow\pi^{-}\pi^{+}\pi^{-}\nu_{\tau}$
control sample. In data, the trigger efficiency is computed using independent
trigger selections: the efficiency of the ECL-based trigger selection is
obtained using events triggered by the CDC, while the efficiency of the CDC-
based trigger selection is evaluated using events passing the ECL trigger
requirements. The agreement between data and simulation efficiencies is 0.5%
for the ECL trigger selection and is 4.3% for the CDC trigger selection. Given
that the efficiency of trigger selections based on the ECL only is 88%, the
weighted average of the discrepancies is computed to be 0.7% and this is used
as the systematic uncertainty of the trigger efficiency. The possible bias
coming from the use of independent trigger selections to evaluate the trigger
efficiency is tested on simulated events. A difference of 0.5% is found with
respect to the absolute efficiency, and added in quadrature to the previous
uncertainty. The $\tau^{-}\rightarrow\pi^{-}\pi^{+}\pi^{-}\nu_{\tau}$ control
sample is also used to obtain a systematic uncertainty on the BDT selection.
The same BDT that was trained for the signal decay is applied to these events
and a requirement is chosen on the output in order to have the same efficiency
on $\tau^{-}\rightarrow\pi^{-}\pi^{+}\pi^{-}\nu_{\tau}$ events as on the
signal decay. The difference between the efficiency in data and simulation for
this BDT requirement, $1.5\%$, is considered as a systematic uncertainty. In
order to take into account possible mismodeling of the ISR, FSR, and
resolution of the signal extraction variables, we change the signal region
definition by $\pm 1\delta$. This results in a ${}^{+2.9}_{-3.9}\%$ variation
of the signal efficiency, that we take as a systematic uncertainty. The number
of expected background events is affected by an uncertainty that originates
from imperfections in the magnetic field description used in the event
reconstruction and material mismodeling. To correct for those effects, the
momenta of charged particles are scaled with a factor of 0.99987, evaluated by
measuring the mass-peak position of high-yield samples of $D^{0}\rightarrow
K^{-}\pi^{+}$ decays reconstructed in data and comparing it to the world
average value of the $D^{0}$ mass ParticleDataGroup:2022pth . The
corresponding systematic uncertainty is obtained by varying the correction
factor according to its uncertainty (${}^{+3.8}_{-5.7}\times 10^{-4}$), which
leads to different data yields in the sideband region and, thus, different
numbers of expected events $N_{\mathrm{exp}}$. The resulting systematic
uncertainty, taken as the difference from the nominal value, is $16$%, and is
mainly due to the small number of events. The systematic uncertainty in the
integrated luminosity $\mathcal{L}$, measured with samples of Bhabha and
diphoton events Belle-II:2019usr , is evaluated from the difference observed
between the results from the two methods and amounts to 0.6% relative
uncertainty. Finally, we also assign an uncertainty of 0.003 nb on the
$\tau$-pair production cross section, as evaluated in Ref. Banerjee:2007is . A
summary of the systematic uncertainties is given in Table 2.

Table 2: Summary of relative systematic uncertainties. |  | Uncertainty (%)
---|---|---
Quantity | Source | Low | High
$\varepsilon_{\mathrm{3\mu}}$ | PID | 2.1 | 2.4
Tracking | 1.0 | 1.0
Trigger | 0.9 | 0.9
BDT | 1.5 | 1.5
Signal region | 3.9 | 2.9
$N_{exp}$ | Momentum Scale | 16 | 16
$\mathcal{L}$ |  | 0.6 | 0.6
$\sigma_{\tau\tau}$ |  | 0.3 | 0.3

In order to perform a validation of the inclusive reconstruction and BDT
selection, we carry out an alternative analysis using one-prong tagging
reconstruction similar to that done previously by Belle and Babar. We select
events containing exactly four charged particles with zero total charge.
Within the tag side, we categorize the charged particle as either a lepton
(leptonic tag) or a hadron (hadronic tag). The contamination from
$e^{+}e^{-}\rightarrow e^{+}e^{-}nh$ events is mitigated using the data-driven
requirements $0.4<\theta_{{\rm{miss}}}^{*}<2.7$$\rm\,rad$ for leptonic and
$0.3<\theta_{{\rm{miss}}}^{*}<2.8$$\rm\,rad$ for hadronic tags, along with a
requirement at most one signal track points to the endcaps. These requirements
also suppress the $e^{+}e^{-}\rightarrow\ell^{+}\ell^{-}\ell^{+}\ell^{-}$ and
$e^{+}e^{-}\rightarrow e^{+}e^{-}h^{+}h^{-}$ background components, which are
further rejected by requiring the difference between the tag tau energy and
half the center-of-mass energy to be less than
$-0.8$$\mathrm{\,Ge\kern-1.00006ptV}$, and $0.90<T<0.97$. The
$e^{+}e^{-}\rightarrow q\bar{q}$ background is suppressed by requiring the
visible c.m. energy of all reconstructed particles in the event ${E_{\rm
vis}^{*}}$ to be smaller than 10.2$\mathrm{\,Ge\kern-1.00006ptV}$ and the
missing momentum $p_{\rm{miss}}^{*}$ to be larger than 0.4
${\mathrm{\,Ge\kern-1.00006ptV\\!/}c}$. To optimize the muon identification
selection, charged particles are separated according to their momentum into
three categories corresponding to tracks not reaching, partially crossing, or
fully crossing the KLM. Selection criteria are optimized independently for the
three categories using the Punzi figure-of-merit. The signal efficiency is
14.9%, which is nearly twice that achieved at the first generation $B$-factory
experiments. This increase is due to a more optimal selection, which includes
low momentum muons on the signal side and
$\tau^{-}\rightarrow\mu^{-}\nu\bar{\nu}$ decays in the tag side, as well as
improved muon identification. The number of background events expected from
simulation is 0.43. The inclusive-tagging selection has a 37% higher
efficiency than that of the one-prong tagging method for a similar level of
background and is thus used as the final result.

##  5 Result

The distribution of events in the ($M_{3\mu}$, $\Delta E_{3\mu}$) plane is
shown in Fig. 4 for data and simulated signal events. We observe one event in
the signal region, $N_{obs}=1$. Using a signal efficiency of
$\varepsilon_{\mathrm{3\mu}}=(20.42\pm 0.06^{+1.02}_{-0.84})\%$, where the
first uncertainty is statistical and the second are systematic, and a number
of expected background events $N_{exp}=0.7_{-0.5-0.1}^{+0.6+0.1}$ in the
formula given in Eq. 3,

$\mathcal{B}(\tau^{-}\rightarrow\mu^{-}\mu^{+}\mu^{-})=\frac{N_{obs}-N_{exp}}{\mathcal{L}\times
2\sigma_{\tau\tau}\times\varepsilon_{\mathrm{3\mu}}},$ (3)

we compute a branching fraction of
$\mathcal{B}(\tau^{-}\rightarrow\mu^{-}\mu^{+}\mu^{-})=(2.1^{+5.1}_{-2.4}\pm
0.4)\times 10^{-9}$, where the first uncertainty is statistical and the second
is systematic. The other inputs to Eq. 3 are the integrated luminosity
$\mathcal{L}=424\pm 3\leavevmode\nobreak\ \mbox{\,fb}^{-1}$ of the analyzed
data sample; the $\tau$-pair production cross section
$\sigma_{\tau\tau}=0.919\pm 0.003\leavevmode\nobreak\ \rm\,nb$, where a
weighted average of the cross sections at the different data taking energies
is used.

Figure 4: Scatter plot of selected events in the ($M_{3\mu}$, $\Delta
E_{3\mu}$) plane for data (orange crosses) and simulated signal (color-filled
area). The SB region is shown as the red rectangle. The yellow ellipse
represents the signal region. Events outside the black rectangle, representing
the $\pm 20\delta$ region, are discarded.

The distributions of all selected events in the $\pm 20\,\delta$ region for
data and simulation for the $M_{3\mu}$ and $\Delta E_{3\mu}$ variables are
shown in Fig. 5.

Figure 5: Comparison between data (black points with error bars) and
simulation for the $M_{3\mu}$ and $\Delta E_{3\mu}$ distributions. Events are
those populating the entire black rectangle shown in Fig. 4. The various
simulated background processes are shown as a stack of color-filled
histograms, while the signal is shown as a red histogram with an arbitrary
scale. The statistical uncertainties are displayed as hatched areas.

As no signal is found, we compute a 90% confidence level (C.L.) upper limit on
the $\tau^{-}\rightarrow\mu^{-}\mu^{+}\mu^{-}$ branching fraction. We estimate
the upper limit using the modified frequentist $\mathrm{CL}_{s}$ Junk:1999kv ;
Read:2002hq  method implemented in the RooStat framework. We generate $5\times
10^{4}$ pseudo-experiments at 40 points distributed uniformly in the branching
fraction range $(0-5)\times 10^{-8}$. The total statistical and systematic
uncertainties affecting each experimental input, discussed in Sect. 4, are
combined in quadrature. Figure 6 shows the $\mathrm{CL}_{s}$ curves computed
as a function of the upper limit on the branching fractions for the inclusive
tagging analysis. The dashed black line shows the expected $\mathrm{CL}_{s}$
and the green and yellow bands give the $\pm 1\sigma$ and $\pm 2\sigma$
contours, respectively. The expected limit assuming an observed number of
events equal to 0.7, as expected from background estimation, is $1.8\times
10^{-8}$ at 90% C.L, while the observed limit on the branching fraction of
$\tau^{-}\rightarrow\mu^{-}\mu^{+}\mu^{-}$ is $1.9\times 10^{-8}$ at 90% C.L.

Figure 6: Observed (solid red curve with dots) and expected (dashed black
curve) $\mathrm{CL}_{s}$ as a function of the assumed branching fractions of
$\tau^{-}\rightarrow\mu^{-}\mu^{+}\mu^{-}$. The red line corresponds to the
90% C.L upper limit.

The limit is also evaluated using a Bayesian approach CALDWELL20092197
yielding the same result. Using the alternative one-prong tagging analysis,
zero events are observed in the signal region. This corresponds to a branching
fraction of
$\mathcal{B}(\tau^{-}\rightarrow\mu^{-}\mu^{+}\mu^{-})=(-3.7^{+1.9}_{-0.2}\pm
0.1)\times 10^{-9}$, where the first uncertainty is statistical and the second
is systematic. The small statistical uncertainty is coming from the fact that
the number of expected background events is derived from simulation. The
corresponding limit is $2.0\times 10^{-8}$ at 90% C.L.

##  6 Summary

We present a search for the LFV decay
$\tau^{-}\rightarrow\mu^{-}\mu^{+}\mu^{-}$ using a 424$\mbox{\,fb}^{-1}$
sample of electron-positron collision data collected by the Belle II
experiment. Using a novel inclusive-tagging reconstruction followed by a BDT-
based selection, the efficiency is higher by a factor of 2.5 than the
efficiency in the latest Belle analysis Hayasaka:2010np  and 37% higher than a
one-prong tagging reconstruction performed on the same Belle II dataset, for
an expected number of background events compatible with zero. We observe one
event in the signal region, which corresponds to a branching fraction of
$\mathcal{B}(\tau^{-}\rightarrow\mu^{-}\mu^{+}\mu^{-})=(2.1^{+5.1}_{-2.4}\pm
0.4)\times 10^{-9}$, where the first uncertainty is statistical and the second
one systematic. The observed (expected) limit at 90% C.L. computed in a
frequentist approach is 1.9 (1.8) $\times 10^{-8}$, which is more restrictive
than the previous lowest limit. This work, based on data collected using the
Belle II detector, which was built and commissioned prior to March 2019, was
supported by Higher Education and Science Committee of the Republic of Armenia
Grant No. 23LCG-1C011; Australian Research Council and Research Grants No.
DP200101792, No. DP210101900, No. DP210102831, No. DE220100462, No.
LE210100098, and No. LE230100085; Austrian Federal Ministry of Education,
Science and Research, Austrian Science Fund No. P 34529, No. J 4731, No. J
4625, and No. M 3153, and Horizon 2020 ERC Starting Grant No. 947006
“InterLeptons”; Natural Sciences and Engineering Research Council of Canada,
Compute Canada and CANARIE; National Key R&D Program of China under Contract
No. 2022YFA1601903, National Natural Science Foundation of China and Research
Grants No. 11575017, No. 11761141009, No. 11705209, No. 11975076, No.
12135005, No. 12150004, No. 12161141008, and No. 12175041, and Shandong
Provincial Natural Science Foundation Project ZR2022JQ02; the Czech Science
Foundation Grant No. 22-18469S; European Research Council, Seventh Framework
PIEF-GA-2013-622527, Horizon 2020 ERC-Advanced Grants No. 267104 and No.
884719, Horizon 2020 ERC-Consolidator Grant No. 819127, Horizon 2020 Marie
Sklodowska-Curie Grant Agreement No. 700525 “NIOBE” and No. 101026516, and
Horizon 2020 Marie Sklodowska-Curie RISE project JENNIFER2 Grant Agreement No.
822070 (European grants); L’Institut National de Physique Nucléaire et de
Physique des Particules (IN2P3) du CNRS and L’Agence Nationale de la Recherche
(ANR) under grant ANR-21-CE31-0009 (France); BMBF, DFG, HGF, MPG, and AvH
Foundation (Germany); Department of Atomic Energy under Project Identification
No. RTI 4002, Department of Science and Technology, and UPES SEED funding
programs No. UPES/R&D-SEED-INFRA/17052023/01 and No. UPES/R&D-SOE/20062022/06
(India); Israel Science Foundation Grant No. 2476/17, U.S.-Israel Binational
Science Foundation Grant No. 2016113, and Israel Ministry of Science Grant No.
3-16543; Istituto Nazionale di Fisica Nucleare and the Research Grants BELLE2;
Japan Society for the Promotion of Science, Grant-in-Aid for Scientific
Research Grants No. 16H03968, No. 16H03993, No. 16H06492, No. 16K05323, No.
17H01133, No. 17H05405, No. 18K03621, No. 18H03710, No. 18H05226, No.
19H00682, No. 20H05850, No. 20H05858, No. 22H00144, No. 22K14056, No.
22K21347, No. 23H05433, No. 26220706, and No. 26400255, and the Ministry of
Education, Culture, Sports, Science, and Technology (MEXT) of Japan; National
Research Foundation (NRF) of Korea Grants No. 2016R1D1A1B02012900, No.
2018R1A2B3003643, No. 2018R1A6A1A06024970, No. 2019R1I1A3A01058933, No.
2021R1A6A1A03043957, No. 2021R1F1A1060423, No. 2021R1F1A1064008, No.
2022R1A2C1003993, and No. RS-2022-00197659, Radiation Science Research
Institute, Foreign Large-Size Research Facility Application Supporting
project, the Global Science Experimental Data Hub Center of the Korea
Institute of Science and Technology Information and KREONET/GLORIAD;
Universiti Malaya RU grant, Akademi Sains Malaysia, and Ministry of Education
Malaysia; Frontiers of Science Program Contracts No. FOINS-296, No. CB-221329,
No. CB-236394, No. CB-254409, and No. CB-180023, and SEP-CINVESTAV Research
Grant No. 237 (Mexico); the Polish Ministry of Science and Higher Education
and the National Science Center; the Ministry of Science and Higher Education
of the Russian Federation and the HSE University Basic Research Program,
Moscow; University of Tabuk Research Grants No. S-0256-1438 and No.
S-0280-1439 (Saudi Arabia); Slovenian Research Agency and Research Grants No.
J1-9124 and No. P1-0135; Agencia Estatal de Investigacion, Spain Grant No.
RYC2020-029875-I and Generalitat Valenciana, Spain Grant No.
CIDEGENT/2018/020; National Science and Technology Council, and Ministry of
Education (Taiwan); Thailand Center of Excellence in Physics; TUBITAK ULAKBIM
(Turkey); National Research Foundation of Ukraine, Project No. 2020.02/0257,
and Ministry of Education and Science of Ukraine; the U.S. National Science
Foundation and Research Grants No. PHY-1913789 and No. PHY-2111604, and the
U.S. Department of Energy and Research Awards No. DE-AC06-76RLO1830, No. DE-
SC0007983, No. DE-SC0009824, No. DE-SC0009973, No. DE-SC0010007, No. DE-
SC0010073, No. DE-SC0010118, No. DE-SC0010504, No. DE-SC0011784, No. DE-
SC0012704, No. DE-SC0019230, No. DE-SC0021274, No. DE-SC0021616, No. DE-
SC0022350, No. DE-SC0023470; and the Vietnam Academy of Science and Technology
(VAST) under Grants No. NVCC.05.12/22-23 and No. DL0000.02/24-25. These
acknowledgements are not to be interpreted as an endorsement of any statement
made by any of our institutes, funding agencies, governments, or their
representatives. We thank the SuperKEKB team for delivering high-luminosity
collisions; the KEK cryogenics group for the efficient operation of the
detector solenoid magnet; the KEK Computer Research Center for on-site
computing support; the NII for SINET6 network support; and the raw-data
centers hosted by BNL, DESY, GridKa, IN2P3, INFN, and the University of
Victoria.

## References

  * (1) T. Li, M.A. Schmidt, C.-Y. Yao and M. Yuan, _Charged lepton flavor violation in light of the muon magnetic moment anomaly and colliders_ , _Eur. Phys. J. C_ 81 (2021) 811 [2104.04494].
  * (2) G. Hernández-Tomé, G. López Castro and P. Roig, _Flavor violating leptonic decays of $\tau$ and $\mu$ leptons in the Standard Model with massive neutrinos_, _Eur. Phys. J. C_ 79 (2019) 84 [1807.06050].
  * (3) P. Blackstone, M. Fael and E. Passemar, _$\tau\rightarrow\mu\mu\mu$ at a rate of one out of $10^{14}$ tau decays?_, _Eur. Phys. J. C_ 80 (2020) 506 [1912.09862].
  * (4) A. Abada, J. Kriewald and A.M. Teixeira, _On the role of leptonic CPV phases in cLFV observables_ , _Eur. Phys. J. C_ 81 (2021) 1016 [2107.06313].
  * (5) M. Raidal et al., _Flavour physics of leptons and dipole moments_ , _Eur. Phys. J. C_ 57 (2008) 13 [0801.1826].
  * (6) A.M. Teixeira, _Theoretical aspects of charged Lepton Flavour Violation_ , _J. Phys. Conf. Ser._ 888 (2017) 012029 [1612.05561].
  * (7) B. Bhattacharya, A. Datta, J.-P. Guévin, D. London and R. Watanabe, _Simultaneous Explanation of the $R_{K}$ and $R_{D^{(*)}}$ Puzzles: a Model Analysis_, _JHEP_ 01 (2017) 015 [1609.09078].
  * (8) F. Feruglio, P. Paradisi and A. Pattori, _On the Importance of Electroweak Corrections for B Anomalies_ , _JHEP_ 09 (2017) 061 [1705.00929].
  * (9) A. Greljo, G. Isidori and D. Marzocca, _On the breaking of Lepton Flavor Universality in B decays_ , _JHEP_ 07 (2015) 142 [1506.01705].
  * (10) Belle collaboration, _Search for Lepton Flavor Violating Tau Decays into Three Leptons with 719 Million Produced $\tau^{+}\tau^{-}$ Pairs_, _Phys. Lett. B_ 687 (2010) 139 [1001.3221].
  * (11) CLEO collaboration, _New limits for neutrinoless tau decays_ , _Phys. Rev. D_ 57 (1998) 5903 [hep-ex/9712010].
  * (12) BaBar collaboration, _Limits on tau Lepton-Flavor Violating Decays in three charged leptons_ , _Phys. Rev. D_ 81 (2010) 111101 [1002.4550].
  * (13) LHCb collaboration, _Search for the lepton flavour violating decay $\tau$- → $\mu$- $\mu$+ $\mu$-_, _JHEP_ 02 (2015) 121 [1409.8548].
  * (14) ATLAS collaboration, _Probing lepton flavour violation via neutrinoless $\tau\rightarrow 3\mu$ decays with the ATLAS detector_, _Eur. Phys. J. C_ 76 (2016) 232 [1601.03567].
  * (15) CMS collaboration, _Search for the lepton flavor violating decay $\tau$ $\rightarrow$ 3$\mu$ in proton-proton collisions at $\sqrt{s}=$ 13 TeV_, _JHEP_ 01 (2021) 163 [2007.05658].
  * (16) CMS collaboration, _Search for the lepton flavor violating $\tau$ → 3$\mu$ decay in proton-proton collisions at $\sqrt{s}$=13TeV_, _Phys. Lett. B_ 853 (2024) 138633 [2312.02371].
  * (17) Belle II collaboration, _Belle II technical design report_ , 1011.0352.
  * (18) K. Akai, K. Furukawa and H. Koiso, _SuperKEKB collider_ , _Nucl. Instrum. Meth._ A907 (2018) 188 [1809.01958].
  * (19) S. Banerjee, B. Pietrzyk, J.M. Roney and Z. Wa̧s, _Tau and muon pair production cross-sections in electron-positron annihilations at $\sqrt{s}=10.58$ GeV_, _Phys. Rev. D_ 77 (2008) 054012 [0706.3235].
  * (20) S. Jadach, B.F.L. Ward and Z. Wa̧s, _The precision Monte Carlo event generator KK for two-fermion final states in $e^{+}e^{-}$ collisions_, _Comput. Phys. Commun._ 130 (2000) 260 [hep-ph/9912214].
  * (21) S. Jadach, J.H. Kuhn and Z. Wa̧s, _TAUOLA: A library of Monte Carlo programs to simulate decays of polarized tau leptons_ , _Comput. Phys. Commun._ 64 (1990) 275.
  * (22) E. Barberio, B. van Eijk and Z. Wa̧s, _PHOTOS: A universal Monte Carlo for QED radiative corrections in decays_ , _Comput. Phys. Commun._ 66 (1991) 115.
  * (23) T. Sjöstrand, S. Ask, J.R. Christiansen, R. Corke, N. Desai, P. Ilten et al., _An Introduction to PYTHIA 8.2_ , _Comput. Phys. Commun._ 191 (2015) 159 [1410.3012].
  * (24) D.J. Lange, _The EvtGen particle decay simulation package_ , _Nucl. Instrum. Meth._ A462 (2001) 152.
  * (25) G. Balossini, C.M. Carloni Calame, G. Montagna, O. Nicrosini and F. Piccinini, _Matching perturbative and parton shower corrections to Bhabha process at flavour factories_ , _Nucl. Phys._ B 758 (2006) 227.
  * (26) G. Balossini, C. Bignamini, C.M.C. Calame, G. Montagna, O. Nicrosini and F. Piccinini, _Photon pair production at flavour factories with per mille accuracy_ , _Phys. Lett._ B 663 (2008) 209.
  * (27) C.M. Carloni Calame, G. Montagna, O. Nicrosini and F. Piccinini, _The BABAYAGA event generator_ , _Nucl. Phys. B Proc. Suppl._ 131 (2004) 48.
  * (28) C.M. Carloni Calame, _An improved parton shower algorithm in QED_ , _Phys. Lett._ B 520 (2001) 16.
  * (29) C.M. Carloni Calame, C. Lunardini, G. Montagna, O. Nicrosini and F. Piccinini, _Large angle Bhabha scattering and luminosity at flavor factories_ , _Nucl. Phys._ B 584 (2000) 459.
  * (30) F. Berends, P. Daverveldt and R. Kleiss, _Radiative corrections to the process $e^{+}e^{-}\rightarrow e^{+}e^{-}\mu^{+}\mu^{-}$_, _Nucl. Phys._ B 253 (1985) 421.
  * (31) F. Berends, P. Daverveldt and R. Kleiss, _Complete lowest-order calculations for four-lepton final states in electron-positron collisions_ , _Nucl. Phys._ B 253 (1985) 441.
  * (32) F. Berends, P. Daverveldt and R. Kleiss, _Monte Carlo simulation of two-photon processes: II: Complete lowest order calculations for four-lepton production processes in electron-positron collisions_ , _Comp. Phys. Commun._ 40 (1986) 285.
  * (33) S. Uehara, _TREPS: A Monte-Carlo Event Generator for Two-photon Processes at $e^{+}e^{-}$ Colliders using an Equivalent Photon Approximation_, 1310.0157.
  * (34) Belle II Framework Software Group, _The Belle II Core Software_ , _Comput. Softw. Big Sci._ 3 (2019) 1 [1809.04299].
  * (35) Belle II collaboration, “Belle II Analysis Software Framework (basf2).” https://doi.org/10.5281/zenodo.5574115.
  * (36) GEANT4 collaboration, _GEANT4: A simulation toolkit_ , _Nucl.Instrum.Meth._ A506 (2003) 250.
  * (37) S. Brandt, C. Peyrou, R. Sosnowski and A. Wroblewski, _The Principal axis of jets. An Attempt to analyze high-energy collisions as two-body processes_ , _Phys. Lett._ 12 (1964) 57.
  * (38) E. Farhi, _A QCD Test for Jets_ , _Phys. Rev. Lett._ 39 (1977) 1587.
  * (39) Particle Data Group, _Review of Particle Physics_ , _PTEP_ 2022 (2022) 083C01.
  * (40) 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, ACM, Aug., 2016, DOI.
  * (41) A. Blum, A. Kalai and J. Langford, _Beating the hold-out: Bounds for k-fold and progressive cross-validation_ , in _Proceedings of the Twelfth Annual Conference on Computational Learning Theory_ , COLT ’99, (New York, NY, USA), p. 203–208, Association for Computing Machinery, 1999, DOI.
  * (42) T. Akiba, S. Sano, T. Yanase, T. Ohta and M. Koyama, _Optuna: A Next-generation Hyperparameter Optimization Framework_ , 1907.10902.
  * (43) G. Punzi, _Sensitivity of searches for new signals and its optimization_ , _eConf_ C030908 (2003) MODT002 [physics/0308063].
  * (44) Belle II collaboration, _Measurement of the integrated luminosity of the Phase 2 data of the Belle II experiment_ , _Chin. Phys. C_ 44 (2020) 021001 [1910.05365].
  * (45) T. Junk, _Confidence level computation for combining searches with small statistics_ , _Nucl. Instrum. Meth. A_ 434 (1999) 435 [hep-ex/9902006].
  * (46) A.L. Read, _Presentation of search results: The CL(s) technique_ , _J. Phys. G_ 28 (2002) 2693.
  * (47) A. Caldwell, D. Kollár and K. Kröninger, _BAT – The Bayesian analysis toolkit_ , _Computer Physics Communications_ 180 (2009) 2197.

## Additional Material

The signal efficiency as function of the two-dimensional plane defined by the
mass squared of the opposite charge muons is provided in Fig. 7.

Figure 7: Signal efficiency distribution in the two-dimensional plane defined
by the mass squared of the opposite charge muons, for candidates passing the
full selection.

Comparison between data and simulation for the most discriminating variables
used in the BDT are shown in Fig. 8.

Figure 8: Comparison between data (black points with error bars) and
simulation for the most discriminating variables used in the BDT: mass of the
rest of the event (upper-left), difference of energy between the rest of event
and the beam (upper-right), transverse momentum of the second highest momentum
muon (bottom-left) and transverse momentum of the lowest momentum muon
(bottom-right). Events in the sideband and signal regions are used. The
various simulated background processes are shown as a stack of color-filled
histograms, while the signal is shown as a red histogram with an arbitrary
scale. The statistical uncertainties are displayed as the grey-hatched areas.