# Hybrid Model For Word Prediction Using Naive Bayes and Latent Information

Henrique X. Goulart, Mauro D. L. Tosi, Daniel Soares-Gonçalves, Rodrigo F.
Maia and Guilherme Wachs-Lopes

###### Abstract

Historically, the Natural Language Processing area has been given too much
attention by many researchers. One of the main motivation beyond this interest
is related to the word prediction problem, which states that given a set words
in a sentence, one can recommend the next word. In literature, this problem is
solved by methods based on syntactic or semantic analysis. Solely, each of
these analysis cannot achieve practical results for end-user applications. For
instance, the Latent Semantic Analysis can handle semantic features of text,
but cannot suggest words considering syntactical rules [1]. On the other hand,
there are models that treat both methods together and achieve state-of-the-art
results, e.g. Deep Learning. These models can demand high computational
effort, which can make the model infeasible for certain types of applications.
With the advance of the technology and mathematical models, it is possible to
develop faster systems with more accuracy. This work proposes a hybrid word
suggestion model, based on Naive Bayes and Latent Semantic Analysis,
considering neighbouring words around unfilled gaps. Results show that this
model could achieve 44.2% of accuracy in the MSR Sentence Completion
Challenge.

###### Index Terms:

Naive Bayes, Latent Semantic Analysis, Sentence Completion, Word Prediction.

##  1 Introduction

In present days, there is a increasing demand by impaired people to use
computer programs, mainly for social interaction. Most of these interactions
are made by textual messages, which makes a difficult task for physical
impaired people to communicate with others. Thereby, one technique which
treats this problem that got attention by researchers is called word
prediction.

Word prediction is a word processing feature that aims reduce the number of
keystrokes necessary for typing words [2]. Usually, these models predict the
next word given a set of words based on a context. Because of that, the
Natural Language Processing (NLP) area, which performs tasks such word
prediction through understanding and interpretation of texts and speeches,
became popular.

One of the first NLP approaches in the computer science field is the n-gram
model. In this model, each gram is a word in a known set of given words (also
known as history). It performs a training in a text database in order to
extract information about its language and create a set of features [3]. This
approach uses the joint probability table, which increases the hit rate, but
that demands too much data knownledge and, many times, it uses unfeasable
computer resources [4]. For that reason, in order to use less resources,
approaches as the Naive Bayes has been developed.

The Naive Bayes is an approach that assumes the conditional independence of
its variables. It requires less memory and processing than approaches that
consider this dependence. However, the Naive Bayes accuracy tends to be worst
than those because of its knowledge loss during the independence assumption
[5].

The NLP area still has researches for an accurate method that could be
implemented and predict in an applicable time. In this context, the latent
semantic analysis (LSA) was created. This technique is used to semantically
analyze texts through the relationship between the words in different text
levels, as phrases, paragraphs, among others [6]. Besides the high accuracy of
the LSA technique, its inferences still are based only on the text frequency,
which means that it does not consider word orders and consequently the text
syntax.

With the advances of technology, in terms of memory space and processing time,
other methods could be developed and consequently implemented, such as the
Deep Learning [1]. Basically, it consists in a set of artificial neurons which
is trained from a database through a multi-layer neural network with the aim
of minimize an error function [7]. The Deep Learning process achieves a great
precision in sentences completion challenges, however, it has some issues. The
computational cost and time to train the amount of data needed is higher
compared to other approaches. Beside this, since the learned knowledge is in
the weights of neurons connections, it is not possible to interpret what was
learned, which means that this is a black box model.

There are many Word prediction models in the NLP area, however, as far as we
know none of them can accomplish the task of predict a word with precision in
an applicable time, as a human being would do. Therefore, the area still
demands development, and a hybrid model might help researches to explore new
approaches.

This paper proposes a new hybrid model to predict words using two well-known
models: the Naive Bayes and the LSA. In addition, to this we optimize
parameters used to improve the prediction precision through the Gradient
Descent technique.

This paper is organized as following. In the next section, will be presented
the relevant areas backgrounds. Then, in the third section, the methodology
used to the model development and its optimization is explained. Next, in the
fourth section, the tests performed to qualify the results that prove the
model precision are presented. Then, in the fifth section, the benefits and
issues are discussed in the paper conclusion. At last, in the sixth section,
some issues that can be improved in future works are presented.

##  2 Background

In NLP area, many models were proposed to predict a word, the Table 2 shows a
list with some related models.

Table I: Related models. Work | Naive Bayes | Other n-gram | LSA | Other model | Year
---|---|---|---|---|---
Hunnicutt and Carlberger [8] | X |  |  |  | 2001
Al-Mubaid [9] | X |  |  | X | 2003
Al-Mubaid [10] | X |  |  | X | 2007
Aliprandi et al. [2] |  |  |  | X | 2008
Zweig et al. [11] |  | X | X |  | 2012
Koutný[12] |  |  |  | X | 2012
Mikolov et al.[1] |  |  |  | X | 2013
Kleinman et al. [13] |  |  |  | X | 2015
Spiccia et al. [14] |  |  | X |  | 2015
Spiccia et al. [15] |  | X | X |  | 2015
Luke and Christianson [16] |  |  | X |  | 2016
Cavalieri et al.[17] | X |  |  |  | 2016

###  2.1 Naive Bayes

The Naive Bayes is a probabilistic model used in Natural Language Processing
(_NLP_) as a n-gram, developed through Bayesian networks that are vastly used
in the machine learning area, as stated in [5].

Usually, full n-gram models are not feasible to be implemented due to its
complexity, because they require a joint probability table. Those tables grow
exponentially with the insertion of new variables.

This model states that its variables can be divided in to _cause_ and _effect_
behaviors; thereby, it can be assumed that the _effects_ are conditionally
independent between themselves, which reduces the computational cost of the
model. This model is mathematically represented in Equation (1), in which, the
effect and the cause are represented by $e$ and $c_{i}$ respectively. Thus,
using this model it is not necessary to build a joint probability table.

$P(e|c_{1},\cdots,c_{n})=\frac{P(e)\prod_{i=1}^{n}P(c_{i}|e)}{\gamma}$ (1)

The normalization factor $\gamma$ is calculated by the Equation (2), in which
$k$ represents the number of words learned in the train.

$\gamma=\sum_{j=1}^{k}P(e_{j}|c_{1},...,c_{n})$ (2)

###  2.2 Latent Semantic Analysis

The latent semantic analysis (LSA) has been vastly used in systems that
analyze textual contents [18], since it performs a comparison between words to
infer.

The LSA can use any cohesive textual level, such as phrases, paragraphs,
entire documents in its traning set [19]. The information (text) of this
textual levels contains the semantic relationship between the words in it. To
store these relations, a table is constructed (_relationship table_)
containing the frequency ($f_{i,j}$) that any word $i$ appeared in a textual
level $j$.

With the _relationship table_ , comparisons between words can be performed to
discover similar words. Therefore, it is necessary to compute the distance
between the _relationship table_ lines using some metrics [20].

Thereby, the LSA is a technique that is used to extract and infer the words
contextual usage relationship in a vector space, as stated in [6].

A common problem of the LSA is the dimension of the _relationship table_ ,
that is usually sparse and consumes unnecessary memory space. A solution used
in the area for this issue [14][19] is the usage of the Singular Value
Decomposition (SVD) technique. With it, it is possible to reduce the
_relationship table_ size without losing all the original information, keeping
it dense. In this work, the reduced table is noted as _Semantic Reduced Table_
($SRT$).

###  2.3 Gradient Descent

The gradient descent is a mathematical optimization algorithm used to modify
variables values based on their contribution to minimize an error function
[21]. Often, this technique is used in prediction models such as: neural
networks, deep learning, bayesian learning, among many others.

To perform the optimization, the error function has to be differentiable,
because the gradient descent algorithm uses derivatives. The most used error
function is the squared-error [21], modeled as Equation (3), in which $f(x)$
is the model prediction, $y$ is the ideal value of the prediction (usually 0
or 1) and the number $\frac{1}{2}$ is used to simplify the further derivatives
calculus, as it will be dropped in the next processes.

$E=\frac{1}{2}(y-f(x))^{2}$ (3)

The variables values changes are proportional to the error function gradient
[22]. Thereby, to minimize the error, the variables have to be updated
proportionally in the opposite direction of gradient function. To smooth the
optimization updates, a factor $\eta$ is often used to control and keep the
current knowledge of the model, compared with the recent information acquired.

##  3 Proposed Model

This paper proposes a hybrid word prediction model, that performs inferences
based on Naive Bayes and Latent Semantic Analysis (LSA) theories. The
methodology used to develop the proposed model is divided into three stages:
Training, Optimization and Inferences, illustrated in Figure 1. Those steps
are described and analyzed in this section.

Figure 1: Proposed model

###  3.1 Training

To perform probabilistic inferences, it is necessary to create and train the
LSA and Naive Bayes sub-models. These sub-models are trained in order to
represent two types of information. The former, Naive Bayes, is used to store
co-occurrences patterns of words; the later, LSA, is used to model the
language semantics.

In order to train this model, a textual database is required. Therefore, the
Project Gutenberg database [23] was used to accomplish this pre-requisite as
it is used vastly in literature such as in [14][24][11]. A set composed by 522
19th Century literature books was used in this paper.

####  3.1.1 Naive Bayes Graph

The co-occurrences patterns which represents the Naive Bayes network can be
stored in a set of graphs $G=(g_{0},g_{1},...,g_{d-1})$, where nodes represent
words and edges represent the number of times that each pair of words co-
occurred in a same textual level. Each graph $g_{i}$ with $0<=i<={n-1}$
represents the connections of the words that appear in the text with $i$ words
between them. This distances are illustrated in the Figures 2 and 3, taking
the phrases ”The sky is blue” and ”The blue is a color” as example.

Figure 2: Words Distance ($d$) for phrase: ”The sky is blue”. The sky is blue
$d$ = 2$d$ = 1$d$ = 0$d$ = 1$d$ = 0$d$ = 0 Figure 3: Words Distance ($d$) for
phrase: ”The blue is a color”. The blue is a color $d$ = 3$d$ = 2$d$ = 1$d$ =
0$d$ = 1$d$ = 0$d$ = 2$d$ = 0$d$ = 1$d$ = 0 Figure 4: Graph 0 with weights
($w$).theskyisblueacolor$w^{0}_{the,sky}$ = 1$w^{0}_{sky,is}$ =
1$w^{0}_{the,blue}$ = 1$w^{0}_{is,a}$ = 1$w^{0}_{a,color}$ =
1$w^{0}_{is,blue}$ = 1$w^{0}_{blue,is}$ = 1

The $G$ set is constructed from many different graphs $g_{i}$ for all $i$.
Each $g_{i}$ is constructed taking every co-occurrences of distance $i$ in a
text. The weight of edges will be higher as many co-occurrences appear in it,
this is represented by $w^{d}_{i,j}$, where $d$ is the graph and $i,j$
represents the edge which connects words $i$ and $j$; the value of these edges
is the number of times that words $i$ and $j$ occurred in the database,
establishing this approach as a frequentist model. Take as example the Figures
4 and 5 that represents the $g_{0}$ and $g_{1}$ respectively which was
constructed using the same sentences used previously: ”The sky is blue” and
”The blue is a color”.

Figure 5: Graph 1 with weights ($w$).theskyisblueacolor$w^{1}_{the,is}$ =
2$w^{1}_{sky,blue}$ = 1$w^{1}_{blue,a}$ = 1$w^{1}_{is,color}$ =
1$w^{1}_{a,color}$ = 1

####  3.1.2 Latent Semantic Analysis (LSA)

This step consists of building the _relationship table_. Firstly, it is
necessary to adopt the textual level to be used in the LSA. In this paper,
only sentences with more than 4 nonstop-words are used. This measure was
established to gather only the most semantic relevant phrases 111The source-
code used as part of implementation of this sub-model can be found in
following link: https://github.com/chiawen/sentence-completion.

To construct the _relationship table_ , each different word in the database is
represented as a row and each textual level as a column. The value of each
table cell $c_{word,t}$ is defined by the number of times that the word is
present in the textual level $t$.

As stated before in Section 2.2, the $STR$ is obtained from _relationship
table_ , where the inferences can be performed using acceptable computer
resources.

###  3.2 Inferences

Based on the $n$ previously used words in a text, inferences can be performed
in the constructed models to predict a word usage probability in the analyzed
context. This section will describe these inferences using the notation
$prev_{0:{n-1}}$ to represent all the $n$ previous words.

####  3.2.1 Naive Bayes inference

The bayesian inferences are computed to establish the
$P(sugg_{j}|prev_{0},...,prev_{n-1})$, in which $sugg_{j}$ is a possible word
recommendation (suggestion) from $prev_{0},...,prev_{n-1}$ history.

Thereby, the Equation (4) represents the probability of the $sugg_{j}$ be used
just after the other $prev_{0:n-1}$ words in the same sentence. This
probability is represented by the notation $B_{j}$. Furthermore, each
$P^{i}(prev_{i}|sugg_{j})$ is a probability created by a normalization of
$w^{d}_{prev_{i},sugg_{j}}$ (explained in Section 3.1.1) illustrated in
Equation (5). The prior probability $P(sugg_{j})$ was gathered directly from
database, which is represented in Equation (6).

$B_{j}=\frac{P(sugg_{j})\prod_{i=0}^{n-1}P^{i}(prev_{i}|sugg_{j})}{\gamma}$
(4)

$P^{d}(i|j)=\frac{w^{d}_{i,j}}{\sum_{j}w^{d}_{i,j}}$ (5)

$P(i)=\frac{\textnormal{number of occurrences of }i}{\textnormal{total of word
occurrences in database}}$ (6)

Therefore, all the previously used words have the same relevance to the model,
as the Naive Bayes considers the events independently, which decreases the
model precision, because as explained previously, the Naive Bayes not
considers all probabilities involved in inference. In order to minimize this
issue, variables responsible to weight this relevance were inserted in the
model, noted as $\lambda_{i}$, in which $i$ refers to the distance between the
weighted and the analyzed words.

The equation $B_{j}$ now represents the weighted model, Figure 7, in which
$P(prev_{0}|sugg_{j})$ is inversely weighted by the $\lambda_{1,...,n-1}$.
Thereby, the value of $\lambda_{1,...,n-1}$ are normalized between themselves,
and the equation energy is preserved.

$B_{j}=\frac{P(sugg_{j})P(prev_{0}|sugg_{j})^{\frac{1}{\lambda_{1}*...*\lambda_{n-1}}}\prod_{i=1}^{n-1}P^{i}(prev_{i}|sugg_{j})^{\lambda_{i}}}{\gamma}$
(7)

####  3.2.2 Latent Information Inference

In order to use the knowledge of the $prev_{0:n-1}$ words, $n$ inferences have
to be computed in the $SRT$, which directly impacts the computational cost of
the model. Those are executed by calculating the semantic similarity of each
$prev_{i}$ word, $0<=i<n-1$, with all $SRT$ words in vocabulary. Those
similarities are calculated by the second norm distance, in which, each vector
represents a word.

Therefore, each of the $n$ $SRT$ cells represents the distance between each
candidate word $j$ and $prev_{i}$. Thus, to obtain the similarity of $j$ to
all $prev$ words, the inverse of the distance vector values are summed, as the
distance is inversely proportional to the semantic similarity of the words.
Thus, the semantic distance $d_{j}$ between $j$ and $prev_{0:n-1}$ is
represented by Equation (8).

$d_{j}=\sum_{i=0}^{n-1}\frac{1}{\|j-i\|_{2}+1}$ (8)

After calculating $d_{j}$, it is necessary to normalize the results based in
$n$, as in Equation (9). Therefore, $L_{j}$ represents the normalized semantic
similarity of the word $j$ based on the $n$ previous known words.

$L_{j}=\frac{1}{n}*d_{j}s$ (9)

####  3.2.3 Hybrid Inference

The trained networks achieved through the Naive Bayes and LSA models, and
their respective inferences, Equations (7) and (9), output different pattern
results, which affects the hybrid model development. The probability variance
of the Naive Bayes inferences is much larger than the LSA ones.

To successfully merge the models, it is necessary to establish an output
pattern to equalized them. Therefore, the probability values of the next word
to be inferred, from both models, are sorted in ascending order crescent in
two vectors. Thereby, those probabilities are replaced by their vector index
and normalized by the sum of all indexes, as illustrated below.

Probabilities | 0.5 | 0.15 | 0.3 | 0.05
---|---|---|---|---
$\Downarrow$
---
Probalities Indexes | 4 | 2 | 3 | 1
---|---|---|---|---
$\Downarrow$
---
Equalized Probabilities | $\frac{4}{10}$ | $\frac{2}{10}$ | $\frac{3}{10}$ | $\frac{1}{10}$
---|---|---|---|---

Thus, the Naive Bayes and LSA inferences are weighted, which creates the
hybrid model, that is represented in Equation (10).

$\theta_{j}=\alpha B_{j}+(1-\alpha)L_{j}$ (10)

The proposed model uses a weighting constant $alpha$, that aims to improve the
inferences precision. This constant represents the Naive Bayes percentage
relevance in comparison with the LSA. For example, if $\alpha=1$ all the
inference will be done through the Naive Bayes else, if $\alpha=0$, through
the LSA.

###  3.3 Optimization

Once Naive Bayes and LSA networks are trained, the inferences can be
performed. However, those models provide a different set of probabilities
($B_{j}$ and $L_{j}$). Considering the hybrid model of Equation (10), there is
a set parameters that can be optimized in order to provide an ideal weight
between the models, improving the recommendation precision. Moreover, the
values of $\lambda$ parameters can also be optimized, achieving a better
weight between the partial probabilities of Naive Bayes.

Therefore, optimization through supervised training was performed. The set of
optimized parameters was composed by $\alpha$, that represents the relevant
percentage of the Naive Bayes compared to the LSA; and
$\lambda_{1},...,\lambda_{n-1}$, that express the weight/relevance relation
between the words distances in the Naive Bayes inference.

As the inference goal is to predict the correct word with the maximum
precision, it is desired that the correct inferred value be as close as
possible to 100%. Thus, the Equation (11) can be determined, in which $E$
represents the system error. The optimization proposed in this work is
performed using only the probability output of the correct suggestion
$sugg_{c}$ (thus, $B_{c}$ and $L_{c}$ are used for optimization).

$E=\frac{1}{2}(1-\theta_{j=c})^{2}\leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ ,\leavevmode\nobreak\ \mbox{for}\leavevmode\nobreak\
c\leavevmode\nobreak\ =\leavevmode\nobreak\ \mbox{``Correct suggestion''}$
(11)

After expanding the equation of error, the Equation (12) is obtained.

$E=\frac{1}{2}\left[1-(\alpha
B_{c}+(1-\alpha)L_{c})\right]^{2}\leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ ,\leavevmode\nobreak\ \mbox{for}\leavevmode\nobreak\
c\leavevmode\nobreak\ =\leavevmode\nobreak\ \mbox{``Correct suggestion''}$
(12)

####  3.3.1 Alpha ($\alpha$)

The gradient descent optimization technique was used to find the ideal value
of $\alpha$. In this case, its value was updated using the error function
through an iterative process. Using Equation (12), the $\alpha$ update is
performed using its derivative over the $\alpha$ itself represented by
Equation (13). Therefore, the error value will tend to decrease.

$\frac{\partial E}{\partial\alpha}=-(1-\theta_{c})*(B_{c}-L_{c})$ (13)

Equations (14) and (15) add a new term to the optimization, in which
$\eta_{\alpha}$ represents the learning rate of each iteration, regarding the
current knowledge of the model.

$\alpha^{\prime}=\alpha-\left(\frac{\partial
E}{\partial\alpha}\right)*\eta_{\alpha}$ (14)

$\alpha^{\prime}=\alpha+(1-\theta_{c})*(B_{c}-L_{c})*\eta_{\alpha}$ (15)

####  3.3.2 Lambda ($\lambda$)

The $\lambda_{1},...,\lambda_{n-1}$ values were optimized similarly to the
$\alpha$ one. A supervised training was performed with thousands of iterations
in which, each $\lambda_{1},...,\lambda_{n-1}$ was updated regarding it
respective error in each iteration.

The error, Equation (12), can be expanded to Equation (16), in which the
equalized $P(prev_{i}|sugg_{j})$, as explained in Section 3.2.3, and the
$P(sugg_{j})$ have been replaced by $t_{i}$ and $t$ respectively.

$E=\frac{1}{2}\left(1-\left(\frac{\alpha*t*...*t_{2}^{\lambda_{2}}*t_{1}^{\lambda_{1}}*t_{0}^{\frac{1}{\lambda_{1}*...\lambda_{n-1}}}}{\gamma}+(1-\alpha)L_{j}\right)\right)^{2}$
(16)

Therefore, with the derivative of the $E$ over each $\lambda$, an individual
error rate value is obtained. The Equations (17) and (18) represent the
derivative of the $E$ regarding $\lambda_{1}$ and $\lambda_{2}$ respectively,
in a model that considers only three previous words.

$\frac{\partial
E}{\partial\lambda_{1}}=-\alpha(1-\theta)\left[\frac{t*t_{2}^{\lambda_{2}}*t_{1}^{\lambda_{1}}*t_{0}^{\frac{1}{\lambda_{1}*\lambda_{2}}}}{\gamma}\left(ln(t_{1})-\frac{ln(t_{0})}{\lambda_{1}^{2}*\lambda_{2}}\right)\right]$
(17)

$\frac{\partial
E}{\partial\lambda_{2}}=-\alpha(1-\theta)\left[\frac{t*t_{2}^{\lambda_{2}}*t_{1}^{\lambda_{1}}*t_{0}^{\frac{1}{\lambda_{1}*\lambda_{2}}}}{\gamma}\left(ln(t_{2})-\frac{ln(t_{0})}{\lambda_{1}*\lambda_{2}^{2}}\right)\right]$
(18)

With the derivative of the $E$ over each $\lambda_{x}$, $1<=x<=n-1$, the
$\lambda_{x}$ value is updated over each iteration, Equation (19), in which
$\eta_{\lambda}$ represents the learning rate of each iteration, regarding the
current knowledge of the model.

$\lambda_{x}^{\prime}=\lambda_{x}-\frac{\partial
E}{\partial\lambda_{x}}*\eta_{\lambda}$ (19)

##  4 Results

In order to validate the proposed model, it was applied the experiment MSR
Sentence Completion Challenge [25], this experiment consists of the 1040
sentences with a missing word. Moreover, it is provided 5 options to choose a
word for complete the sentence.

In the MSR Sentence Completion Challenge, it is possible that the missing word
is in the middle of sentence. For this case, the bayesian model also considers
the posterior words from the gap. This case is illustrated in Figure 6 using
as example the following sentence: _Lorem ipsum dolor ___ amet vitae elit_.
Thereby, the LSA also consider all words it was considered by Naive Bayes.

Figure 6: Example of the bayesian inference for MSR Sentence Completion
Challenge.?LoremipsumdolorDistance 2Distance 1Distance 0ametvitaeelitDistance
-2Distance -1Distance 0

Thus, the Naive Bayes inference was described by Equation (20).

$B_{j}=\frac{t*b_{2}^{\lambda_{b2}}*b_{1}^{\lambda_{b1}}*b_{0}^{\frac{1}{\lambda_{b1}*\lambda_{b2}}}*a_{0}^{\frac{1}{\lambda_{a-1}*\lambda_{a-2}}}*a_{-1}^{\lambda_{a-1}}*a_{-2}^{\lambda_{a-2}}}{\gamma}$
(20)

To guarantee cross-validation on experiments, the MSR Sentence Completion
Challenge was separated into 5 groups. Therefore, one of these groups was used
for test and the others for optimization, this is illustrated in Figure 7.
Furthermore, is possible to have 5 different and independent configurations
for experiments.

Figure 7: Experiment configurations.Group 1Group 2Grupo 3Group 4Group 5Config.
1Group 1Group 2Group 3Group 4Group 5Config. 2Group 1Group 2Group 3Group 4Group
5Config. 3Group 1Group 2Group 3Group 4Group 5Config. 4Group 1Group 2Group
3Group 4Group 5Config. 5OptimizationTest

###  4.1 Error Variation

In order to validate the proposed optimization, it was created a scenario
where the parameters $\lambda$ and $\alpha$ were initialized randomly in the
range 0 and 1 for each configuration. Furthermore, the optimization groups
were disposed in a cyclic loop (i.e. epochs), thus extending the optimization.
In this scenario, it was used a 3-gram history.

Figure 8: Variation of error throughout the epochs.

The Figure 8 shows the error variation in each configuration. It is possible
to observe that the error is minimized for all configurations. This error
metric is obtained through the sum of all individual sentences errors.

###  4.2 Alpha Variation

Figure 9: Variation (interpolated) of alpha.

Following the same scenario, the Figure 9 illustrates the variation of alpha
value over epochs. Note that $\alpha$ lies between $\sim$0.2 and $\sim$0.4.
Therefore it can be stated that LSA is more relevant than Naive Bayes, since
this value is lower than 0.5.

###  4.3 Lambda Optimization

In the scenario described above, it is possible to visualize the values of
parameters $\lambda$ in each configuration over the optimization sentences, as
illustrated in Figures 10, 11, 12, 13, 14 and 15.

Figure 10: Variation (interpolated) of previous word $\lambda$ about distance
0. Figure 11: Variation (interpolated) of previous word $\lambda$ about
distance 1. Figure 12: Variation (interpolated) of previous word $\lambda$
about distance 2. Figure 13: Variation (interpolated) of posterior word
$\lambda$ about distance 0. Figure 14: Variation (interpolated) of posterior
word $\lambda$ about distance -1. Figure 15: Variation (interpolated) of
posterior word $\lambda$ about distance -2.

Results show that $\lambda$ values converge to the same region even with the
initial random values and configurations. Therefore, it is feasible to state
that in this scenario the word about distance 0 is the most relevant for
inference and the words about distance 2 and -2 are the least important for
inference.

###  4.4 Size of History Influence

Figure 16: Error variation depending on the size of history.

In order to understand how the size of history can influence the final result,
it was created a scenario where the history was tuned between 2-gram and
15-gram using configuration 1. Furthermore, the proposed model was tested
using only Naive Bayes, only LSA and weighing both models (i.e. $\alpha$
optimized), as illustrated in Figure 16.

As observed, up to 4-gram there are a weighing between both models that
provide better results than these models individually. From that point, LSA
provides better results alone. This probably happens because the bayesian
model tries to detect patterns which were not present in training.

###  4.5 Influence of Each Word for the Naive Bayes

Through the Naive Bayes it is possible to determine which distances are more
or less important for inference, observing the parameters $\lambda$. Thereby,
it was created a scenario where a 15-gram history is used and only Naive Bayes
is considered (i.e. $\alpha$ = 0) as illustrated in Figure 17.

Figure 17: Value of each $\lambda$ in a 15-gram history.

Thereby, both sides of inference have a similar pattern, where the distance 0
words are more relevant.

###  4.6 Tests and Comparative

Figure 18: Results of MSR Sentence Completion Challenge and comparative with
(A) RNNLM + Skip-gram [1], (B) LSA [25] and (C) N-gram [24].

As explained previously, the proposed model in this work was tested in each
configuration using a 3-gram history. The Figure 18 shows the performance of
model optimized and non-optimized. Furthermore, it is compared with the state-
of-the-art [1] and related models.

As observed, the proposed model in this work is superior to another n-gram
model which considers history between 5-gram; but it is worse than LSA model,
which considers all words in the sentence for inference; also it is worse than
neural network model. However, the proposed model provides recommendations
considering syntactic and semantic features using less training data than a
neural network model.

##  5 Conclusion

In this work, it was developed a hybrid model to complete sentences employing
Naive Bayes and LSA models. This hybrid model uses a set words to infer the
next word considering each sub-model individually.

In order to validate the proposed hybrid model, it was conducted an experiment
in which it is composed by sentences with a missing word. Therefore, there
were proposed the parameters $\lambda_{1},...,\lambda_{n-1}$ to optimize the
inferences of Naive Bayes and another parameter $\alpha$ to improve their
integration with the LSA. To optimize the values of these parameters, it was
applied the Gradient Descent technique, where these values are adjusted based
on an error function.

Results show that lambda values converge to same values. That is, for each
lambda $\lambda$, different experiments result to similar final values. This
is a strong evidence that the optimization technique is working. Furthermore,
results show that the most relevant history information for the suggestion are
the previous and posterior neighboring words, as illustrated in Figure 17.

Within the optimization experiment, the $\alpha$ parameter was also trained.
Results show that $\alpha$ tends to give more importance to LSA model than
Naive Bayes, since its value converged to the region between $\sim$0.2 and
$\sim$0.4. Another useful result is that the size of history words affects the
convergence value of $\alpha$. That is, if the size of history is less than 4,
$alpha$ lies between $0$ and $1$. However, if the size of history is higher
than 4, $alpha$ tends to $0$, which means that LSA model is the unique model
used to suggest the next word. Maybe, this is a evidence that the Naive Bayes
model could not find suitable conditional probabilities in the Naive Bayes
graph. This behavior could happen if the train database is small and does not
contains some relations between words.

One important regard about the experiments is that it was used the cross-
validation in order to validate the convergence values. It was found that, for
different sets of train samples, the $\lambda$ and $\alpha$ parameters did not
change the final optimized value.

During tests, the proposed model proved to be relevant to other related
models. Even not providing better results than actual state-of-the-art, the
proposed model brought relevant results, considering syntactic/semantic rules
among words and consuming a fraction of training required by a model based on
Deep Learning.

Thereby, the proposed model in this work achieved the objectives providing a
great time $\times$ accuracy relation which can be improved, thus creating new
perspectives in NLP as also different applications such: a better word
prediction model for mobile platforms, improve the conversation ability of
robots, enhance search engines and improve the social interaction of impaired
people.

##  6 Future Work

The proposed model has several features which were not investigated. One of
these features is to use a bigger database for optimization since the MSR
Sentence Completion Challenge has a few supervised examples to do this task.

In this work it was used a fixed number for learning rate $\eta_{\alpha}$ and
$\eta_{\lambda}$. In case of those number were changed over the epochs, could
provide better results.

It was proposed a minimal conditional probability for the Naive Bayes
inferences whenever a co-occurrence is not present in the words graph. This
gives a minimal chance for a word to be recommended, instead of zeroed product
in the bayesian inference. In case of this value is optimized, maybe the
bayesian model provides greater results.

The database from Project Gutenberg gives a set of full books (where are
included meta-information such as: chapter indicator, paragraphs indicators,
publisher notes, etc). In this work, we used a simple regular expression rule
to extract useful information from this database. However, much information
was lost during this process. It would be better to use text mining algorithms
to extract more information in order to get a bigger dataset.

Over the optimization, only the correct option was considered for error
minimization. Also, it is possible to consider wrong words to compute the
error function. This could bring more precision to the inference process.

To fairly compare the Naive Bayes and LSA, it was proposed a equalization
procedure of output probabilities for each model. This technique worsens the
optimization of parameters $\lambda$ when giving normalized values instead of
the original. Thereby, another technique that avoids this problem could
improve the proposed model.

In Naive Bayes, it was added many parameters $\lambda$ to enhance the bayesian
performance. However, the prior probability does not have one of these
parameters. In case of implementing this new $\lambda$ term, it would be
possible to observe not only the influence of previous words than posterior
but also the prior probability term.

## Acknowledgments

The authors would like to thank Telefônica VIVO partnership with Centro
Universitário FEI and IIoT Lab Telefônica/FEI (Intelligent Internet of Things
Lab) by the funding and hardware resources provided.

## References

  * [1] T. Mikolov, K. Chen, G. Corrado, and J. Dean, “Efficient estimation of word representations in vector space,” _arXiv preprint arXiv:1301.3781_ , 2013\.
  * [2] C. Aliprandi, N. Carmignani, N. Deha, P. Mancarella, and M. Rubino, “Advances in nlp applied to word prediction,” 2008.
  * [3] W. B. Cavnar, J. M. Trenkle _et al._ , “N-gram-based text categorization,” _Ann Arbor MI_ , vol. 48113, no. 2, pp. 161–175, 1994.
  * [4] T. Brants, “Natural language processing in information retrieval.” in _CLIN_ , 2003.
  * [5] S. Russell and P. Norvig, _Artificial Intelligence: A Modern Approach_ , 3rd ed. Upper Saddle River, NJ, USA: Prentice Hall Press, 2009.
  * [6] T. K. Landauer, P. W. Foltz, and D. Laham, “An introduction to latent semantic analysis,” _Discourse processes_ , vol. 25, no. 2-3, pp. 259–284, 1998.
  * [7] I. Goodfellow, Y. Bengio, and A. Courville, _Deep Learning_. MIT Press, 2016, http://www.deeplearningbook.org.
  * [8] S. Hunnicutt and J. Carlberger, “Improving word prediction using markov models and heuristic methods,” _Augmentative and Alternative Communication_ , vol. 17, no. 4, pp. 255–264, 2001.
  * [9] H. Al-Mubaid, “Context-based word prediction and classification,” in _Proceedings of the 18th International Conference on Computers and their Applications CATA’2003_ , 2003, pp. 384–388.
  * [10] ——, “A learning-classification based approach for word prediction.” _Int. Arab J. Inf. Technol._ , vol. 4, no. 3, pp. 264–271, 2007.
  * [11] G. Zweig, J. C. Platt, C. Meek, C. J. Burges, A. Yessenalina, and Q. Liu, “Computational approaches to sentence completion,” in _Proceedings of the 50th Annual Meeting of the Association for Computational Linguistics: Long Papers-Volume 1_. Association for Computational Linguistics, 2012, pp. 601–610.
  * [12] M. Koutný, “Word prediction using language models,” 2012.
  * [13] D. Kleinman, E. Runnqvist, and V. S. Ferreira, “Single-word predictions of upcoming language during comprehension: Evidence from the cumulative semantic interference task,” _Cognitive psychology_ , vol. 79, pp. 68–101, 2015.
  * [14] C. Spiccia, A. Augello, G. Pilato, and G. Vassallo, “A word prediction methodology for automatic sentence completion,” in _Semantic Computing (ICSC), 2015 IEEE International Conference on_. IEEE, 2015, pp. 240–243.
  * [15] C. Spiccia, A. Augello, and G. Pilato, “A word prediction methodology based on posgrams,” in _International Joint Conference on Knowledge Discovery, Knowledge Engineering, and Knowledge Management_. Springer, 2015, pp. 139–154.
  * [16] S. G. Luke and K. Christianson, “Limits on lexical prediction during reading,” _Cognitive Psychology_ , vol. 88, pp. 22–60, 2016.
  * [17] D. C. Cavalieri, S. E. Palazuelos-Cagigas, T. F. Bastos-Filho, and M. Sarcinelli-Filho, “Combination of language models for word prediction: an exponential approach,” _IEEE/ACM Transactions on Audio, Speech, and Language Processing_ , vol. 24, no. 9, pp. 1481–1494, 2016.
  * [18] K. Zupanc and Z. Bosnić, “Automated essay evaluation with semantic analysis,” _Knowledge-Based Systems_ , 2017.
  * [19] N. Coccaro and D. Jurafsky, “Towards better integration of semantic predictors in statistical language modeling.” in _ICSLP_ , 1998, pp. 2403–2406.
  * [20] A. Zhila, W.-t. Yih, C. Meek, G. Zweig, and T. Mikolov, “Combining heterogeneous models for measuring relational similarity.” in _HLT-NAACL_ , 2013, pp. 1000–1009.
  * [21] I. H. Witten, E. Frank, M. A. Hall, and C. J. Pal, _Data Mining: Practical machine learning tools and techniques_. Morgan Kaufmann, 2016.
  * [22] R. Z. Homod, K. S. M. Sahari, H. A. Almurib, and F. H. Nagi, “Gradient auto-tuned takagi–sugeno fuzzy forward control of a hvac system using predicted mean vote index,” _Energy and Buildings_ , vol. 49, pp. 254–267, 2012.
  * [23] (2017, nov) Project gutenberg. [Online]. Available: www.gutenberg.org
  * [24] J. Gubbins and A. Vlachos, “Dependency language models for sentence completion.” in _EMNLP_ , vol. 13, 2013, pp. 1405–1410.
  * [25] G. Zweig and C. J. Burges, “The microsoft research sentence completion challenge,” Technical Report MSR-TR-2011-129, Microsoft, Tech. Rep., 2011.


# Reanalysis of the top-quark pair production via the $e^{+}e^{-}$
annihilation near the threshold region up to N3LO QCD corrections

Jiang Yan<EMAIL_ADDRESS>Xing-Gang Wu<EMAIL_ADDRESS>Zhi-Fei Wu
<EMAIL_ADDRESS>Jing-Hao Shan<EMAIL_ADDRESS>Department of Physics, Chongqing
Key Laboratory for Strongly Coupled Physics, Chongqing University, Chongqing
401331, P.R. China  Hua Zhou<EMAIL_ADDRESS>School of Science, Southwest
University of Science and Technology, Mianyang 621010, P.R. China

###### Abstract

In this paper, we present an improved analysis of the top-quark pair
production via the process $e^{+}e^{-}\to\gamma^{*}\to t\bar{t}$ near the
threshold region up to next-to-next-to-next-to-leading order (N3LO) QCD
corrections. Near the threshold region, the top-quark velocity $v$ tends to
zero, leading to Coulomb singularity. To achieve a reasonable prediction in
the threshold region, we reconstruct the analytical expression for the
Coulomb-terms up to N3LO accuracy by using the PSLQ algorithm, whose numerical
values agree well with the previous N3LO-level calculations. It is found that
the N3LO series still has sizable renormalization scale dependence, and to
improve the precision of the series, we apply the Principle of Maximum
Conformality to eliminate such scale dependence. After that, the Coulomb part
is resummed into a Sommerfeld-Gamow-Sakharov factor, which finally leads to a
much more reasonable behavior near the threshold region.

##  I Introduction

The production of top-quark pair in electron-positron annihilation [1, 2, 3,
4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14] holds substantial relevance for both the
validation of the Standard Model (SM) [15] and the pursuit of new physics
beyond it. As the heaviest species within the SM, a comprehensive probe into
the diverse properties of the top quark is anticipated via the exploration of
its threshold production. For instance, pinpoint theoretical estimations of
the threshold production cross section for $e^{+}e^{-}\to t\bar{t}$ facilitate
the ascertainment of the top quark’s exact mass, which also serves as a vital
observable for future high-energy electron-positron colliders such as the
International Linear Collider (ILC) [16], the Circular Electron Positron
Collider (CEPC) [17] and the Future Circular Collider [18], and etc.

Recently, the full next-to-next-to-next-to-leading order (N3LO) QCD
corrections to the heavy quark pair production $e^{+}e^{-}\to\gamma^{*}\to
Q\bar{Q}$ has been completed in Ref.[19], which still exhibits significant
renormalization scale dependence. Especially, the Coulomb-like interaction
between the quark and the anti-quark triggers a pronounced enhancement of the
cross section close to the production threshold, which has been incorporated
to all orders in perturbation theory. To refine the current theoretical
calculation, it is essential to address the mentioned renormalization scale
dependence and the Coulomb singularity near the threshold. This can be
accomplished by investigating renormalization scale-setting problem,
estimating unknown higher-order (UHO) contributions, and resumming all of the
Coulomb corrections. Such an endeavor will not only lead to a more robust
theoretical calculation but also provide a more precise perturbative QCD
(pQCD) prediction on the top-quark pair production cross section near the
threshold region.

The fixed-order pQCD predictions for observables are usually assumed to suffer
from an uncertainty in fixing the renormalization scale. This uncertainty in
making fixed-order predictions occurs because one usually assumes an arbitrary
renormalization scale together with an arbitrary range to ascertain its
uncertainty. This simple treatment of the pQCD series and the renormalization
scale causes the coefficients of the QCD running coupling ($\alpha_{s}$) at
each order to be strongly dependent on the choice of both the renormalization
scale and the renormalization scheme. The $\alpha_{s}$-running behavior is
governed by the renormalization group equation (RGE),

$\displaystyle\beta(\alpha_{s})=\frac{{\rm d}\alpha_{s}(\mu^{2})}{{\rm
d}\ln\mu^{2}}=-\sum_{i=0}^{\infty}\beta_{i}\alpha_{s}^{i+2}(\mu^{2}),$ (1)

where the $\\{\beta_{i}\\}$-functions have been computed up to five-loop level
in modified minimal-subtraction scheme ($\overline{\rm MS}$-scheme) [20, 21,
22, 23, 24, 25, 26, 27, 28, 29]. Consequently, the $\\{\beta_{i}\\}$-terms
emerged in pQCD series can be reabsorbed into $\alpha_{s}$ to accurately fix
the correct magnitude of $\alpha_{s}$. In respect of this fact, the Principle
of Maximum Conformality (PMC) [30, 32, 33, 31, 34] has been proposed to
address the scale-setting problem in fixed-order pQCD calculation. It
eliminates conventional renormalization scale ambiguity by using the above RGE
recursively so as to set the correct $\alpha_{s}$ value and to achieve a well-
matched (maximumly-approaching) conformal series between $\alpha_{s}$ and its
expansion coefficients, yielding naturally scheme-independent theoretical
predictions at any fixed-order, cf. the recent reviews [35, 36, 37, 38]. This
treatment is a kind of resummation of all the known type of RG-involved
$\\{\beta_{i}\\}$-terms. More specifically, the PMC single-scale-setting
approach (PMCs) [39, 40] determines a global effective coupling
$\alpha_{s}(Q_{*}^{2})$ ($Q_{*}$ is called as the PMC scale which can be
viewed as the effective momentum flow of the process) for any fixed-order
prediction by using the RGE, which also provides a solid way to extend the
well-known Brodsky-Lepage-Mackenzie (BLM) approach [41] to all orders 111A
short review of the development of PMC from BLM and replies to some wrong
comments on the PMC can be found in Ref.[42].. Using the PMCs, it has been
demonstrated that the PMC prediction is independent to any choice of
renormalization scheme and scale [43], being consistent with the fundamental
renormalization group approaches [44, 78, 74, 47] and also the self-
consistency requirements of the renormalization group such as reflectivity,
symmetry and transitivity [48].

It has been found that Coulomb-type corrections for the heavy-quark pair
($Q\bar{Q}$) production via the electron-positron annihilation have sizable
contributions in the threshold region [49, 50], which similar to the hadronic
production case [51] shall also be enhanced by factors of $\pi$ and the
resultant PMC scale can be relatively small when the heavy quark velocity
$v\to 0$ [52]. Thus the terms which are proportional to various powers of
$(\pi/v)$ need to be treated carefully. In Ref.[52], the total cross section
of $e^{+}e^{-}\to\gamma^{*}\to Q\bar{Q}$ near the threshold region has been
investigated up to next-to-next-to leading order (N2LO) level. At the N3LO
level and higher, the question will be much more involved due to the much more
complex Coulomb structures. In this paper, we will begin by partitioning the
N3LO-level total cross section of $e^{+}e^{-}\to\gamma^{*}\to Q\bar{Q}$ into
two components, Coulomb and non-Coulomb corrections. Subsequently, we will
employ the PMCs approach to deal with each components individually. Following
that, we will perform a resummation of the divergent Coulomb corrections using
the well-established Sommerfeld-Gamow-Sakharov (SGS) factor [53]. Then as an
explicit example, we will give the numerical results for the case of top-quark
pair production. Finally, we will utilize the Bayesian analysis (BA) [54, 55,
56, 57] to estimate the contributions from the UHO-terms.

The remaining parts of paper are organized as follows. In Sec.II, we present
the total cross section of $e^{+}e^{-}\to\gamma^{*}\to Q\bar{Q}$ near the
threshold region up to N3LO-level, and explain how to deal with its Coulomb
and non-Coulomb components separately. And then we derive the analytical
expression of the Coulomb part by using the integer relation finding
algorithm, e.g. the so-called PSLQ algorithm [58, 59]. The results for the
N3LO-level Coulomb coefficients will be given in the Appendix. In Sec.III, we
give a concise introduction to the PMCs method. Numerical results and
discussions for the top-quark pair production are given in Sec.IV. Sec.V is
reserved for a summary.

##  II The heavy quark pair production via the $e^{+}e^{-}$ annihilation near
the threshold region

Up to N3LO level QCD corrections, the total cross section of
$e^{+}e^{-}\to\gamma^{*}\to Q\bar{Q}$ can be written as

$\displaystyle\sigma(s)=\sigma_{0}\left(1+r_{1}\alpha_{s}+r_{2}\alpha_{s}^{2}+r_{3}\alpha_{s}^{3}+{\cal
O}(\alpha_{s}^{4})\right),$ (2)

where $\alpha_{s}\equiv\alpha_{s}^{(n_{l})}$ is QCD running coupling in
$\overline{\rm MS}$-scheme with $n_{l}$ represents the active number of light
flavors, and $\sigma_{0}$ is the LO cross section which takes the form

$\displaystyle\sigma_{0}=N_{C}\frac{4\pi\alpha^{2}}{3s}\frac{v(3-v^{2})}{2}e_{Q}^{2}.$
(3)

Here $N_{C}=3$ for $SU(3)_{C}$ color group, $\alpha$ is electromagnetic
coupling constant, $v=\sqrt{1-4m_{Q}^{2}/s}$ is heavy quark velocity, $m_{Q}$
is the on-shell mass of the heavy quark $Q$, $s$ is the squared $e^{+}e^{-}$
collision energy, and $e_{Q}$ is the charge of the heavy quark $Q$ in units of
the electric charge, which equals to $+2/3$ for $Q=t$. It is convenient to
express the perturbative coefficients $r_{1}$, $r_{2}$ and $r_{3}$ as

$\displaystyle r_{1}$ $\displaystyle=\frac{1}{v}r_{1,v}+r_{1,+},$ (4)
$\displaystyle r_{2}$
$\displaystyle=\frac{1}{v^{2}}r_{2,v^{2}}+\frac{1}{v}r_{2,v}+r_{2,+},$ (5)
$\displaystyle r_{3}$
$\displaystyle=\frac{1}{v^{2}}r_{3,v^{2}}+\frac{1}{v}r_{3,v}+r_{3,+},$ (6)

where $r_{i,+}$ represents the terms of non-negative powers of $v$ in $r_{i}$
with $i=(1,2,3)$, respectively. Clearly, in the vicinity of the two-particle
threshold, i.e., as the limit $v\to 0$, the cross section exhibits the well-
known Coulomb singularity, which is dominated by the negative power term of
$v$ in Eqs.(4,5,6).

The analytical results of $r_{1}$, $r_{2}$ and all $n_{l}$-terms in $r_{3}$
can be found in Ref.[11, 12], and the semi-numerical expression of $r_{3}$ can
be found in Ref.[13, 14]. In Ref.[19], numerical results up to N3LO are
expressed as asymptotic expansions of $x=4m_{Q}^{2}/s$, which have been
expanded up to $40_{\rm th}$-orders. Using those numerical results, we are
able to reconstruct the analytical results of $r_{3}$ by using the PSLQ
algorithm [58, 59]. The PSLQ algorithm can be used to reconstruct analytical
expression for given high-precision numerical values. As an explanation, the
inputs to the PSLQ algorithm are high-precision numerical value denoted as
$y$, and a vector $x=(x_{1},\cdots,x_{n})$ consisting of real numbers that are
analytically well-defined, such as $1$, $\ln 2$, $\pi$, ${\rm Li}_{n}(1/2)$,
$\zeta_{n}$ and so on. Subsequently, the PSLQ algorithm determines a set of
integer coefficients $m_{0},m_{1},\cdots,m_{n}$, such that the linear
combination with these integer coefficients is close to zero within the
required numerical precision, i.e., $m_{0}y+\sum_{i=1}^{n}m_{i}x_{i}=0$.
Therefore, the analytical expression of $y$ can be written as
$y=-\sum_{i=1}^{n}\frac{m_{i}}{m_{0}}x_{i}$. The PSLQ algorithm can be
implemented using the Mathematica package PolyLogTools [60]. Our results of
those coefficients at the scale $\mu=\sqrt{s}$ are given in the Appendix.

Note that the coefficients in Eq.(2) are parameterized in terms of
$\alpha_{s}^{(n_{l})}$, which can be transformed from the results
parameterized in terms of $\alpha_{s}^{(n_{l}+1)}$ by relating the coupling
for the theory with $n_{l}+1$ flavors to the theory with $n_{l}$ flavors via
the following decoupling relationship,

$\displaystyle\frac{\alpha_{s}^{(n_{l}+1)}(\mu^{2})}{\alpha_{s}^{(n_{l})}(\mu^{2})}=1+\frac{1}{3}T_{F}\ln\frac{\mu^{2}}{m_{h}^{2}}\frac{\alpha_{s}^{(n_{l})}(\mu^{2})}{\pi}+T_{F}\left[-\frac{2}{9}C_{A}+\frac{15}{16}C_{F}+\left(\frac{5}{12}C_{A}+\frac{1}{4}C_{F}\right)\ln\frac{\mu^{2}}{m_{h}^{2}}+\frac{1}{9}T_{F}\ln^{2}\frac{\mu^{2}}{m_{h}^{2}}\right]\left(\frac{\alpha_{s}^{(n_{l})}(\mu^{2})}{\pi}\right)^{2},$
(7)

where $m_{h}$ is the on-shell mass of the decoupled heavy quark, $C_{A}=3$,
$C_{F}=4/3$ and $T_{F}=1/2$.

Since the Coulomb singularity plays a dominant role for the total cross
section near the threshold region ($v\to 0$), we can factorize the total cross
section (2) into the Coulomb and non-Coulomb parts, respectively, i.e.,

$\displaystyle\sigma=\sigma_{0}\times\mathcal{R}_{\rm
NC}\times\mathcal{R}_{\rm C},$ (8)

where $\mathcal{R}_{\rm C}$ and $\mathcal{R}_{\rm NC}$ represent the Coulomb
part and the non-Coulomb part of the total cross section, respectively, which
can be expressed as

$\displaystyle\mathcal{R}_{\rm C}$ $\displaystyle=1+r_{1}^{\rm
C}\alpha_{s}^{\rm V}(sv^{2})+r_{2}^{\rm C}\alpha_{s}^{\rm
V,2}(sv^{2})+r_{3}^{\rm C}\alpha_{s}^{\rm V,3}(sv^{2}),$ (9)
$\displaystyle\mathcal{R}_{\rm NC}$ $\displaystyle=1+r_{1}^{\rm
NC}\alpha_{s}(s)+r_{2}^{\rm NC}\alpha_{s}^{2}(s)+r_{3}^{\rm
NC}\alpha_{s}^{3}(s).$ (10)

Here the QCD coupling $\alpha_{s}^{\rm V}(\vec{q}^{2})$ has been introduced
for describing the interaction of the non-relativistic heavy quark-antiquark
pair, which is defined as the effective charge in the following Coulomb-like
potential [61, 62, 63, 64, 65, 66, 67, 68]:

$\displaystyle V(\vec{q}^{2})=-4\pi C_{F}\frac{\alpha_{s}^{\rm
V}(\vec{q}^{2})}{\vec{q}^{2}},$ (11)

where $\alpha_{s}^{\rm V}(\vec{q}^{2})$ absorbs all the higher-order QCD
corrections, which is related to the $\overline{\rm MS}$-scheme coupling via
the following way

$\displaystyle\frac{\alpha_{s}^{\rm
V}(\mu^{2})}{\alpha_{s}(\mu^{2})}=1+a_{1}\alpha_{s}(\mu^{2})+a_{2}\alpha_{s}^{2}(\mu^{2})+\mathcal{O}(\alpha_{s}^{3}),$
(12)

where the expansion coefficients $a_{1}$ and $a_{2}$ have been calculated in
Ref.[68], and their specific forms are as follows:

$\displaystyle a_{1}=$
$\displaystyle\,\frac{1}{4\pi}\left(\frac{31}{9}C_{A}-\frac{20}{9}T_{F}n_{l}\right),$
(13) $\displaystyle a_{2}=$
$\displaystyle\,\frac{1}{(4\pi)^{2}}\bigg{[}\left(\frac{4343}{162}+4\pi^{2}-\frac{\pi^{2}}{4}+\frac{22}{3}\zeta_{3}\right)C_{A}^{2}$
$\displaystyle\,-\left(\frac{1798}{81}+\frac{56}{3}\zeta_{3}\right)C_{A}T_{F}n_{l}-\left(\frac{55}{3}-16\zeta_{3}\right)C_{F}T_{F}n_{l}$
$\displaystyle\,+\left(\frac{20}{9}T_{F}n_{l}\right)^{2}\bigg{]},$ (14)

where $\zeta_{n}=\zeta(n)$ is Riemann’s Zeta function. Then we can obtain

$\displaystyle r_{1}^{\rm C}$ $\displaystyle=C_{F}\frac{\pi}{2v},$ (15)
$\displaystyle r_{2}^{\rm C}$
$\displaystyle=C_{F}^{2}\frac{\pi^{2}}{12v^{2}},$ (16) $\displaystyle
r_{3}^{\rm C}$
$\displaystyle=C_{F}\left(\frac{\pi^{3}}{3v}\beta_{0}^{2}(n_{l})-C_{F}\frac{2\zeta_{3}}{v^{2}}\beta_{0}(n_{l})\right),$
(17)

and the coefficients $r_{i}^{\rm NC}$ can be obtained by expanding
$\sigma(s)/(\sigma_{0}\mathcal{R}_{\rm C})$ over $\alpha_{s}$ and $v$.

##  III Simple introduction of principle of maximum conformality single-scale-
setting method

In the previous literature [39, 40], it has been found that the PMCs approach
provides an all-orders single-scale-setting method to fix the conventional
renormalization scale ambiguity, which adopts the standard RGE and the general
QCD degeneracy relations among different orders [69] to transform the RG-
involved $n_{f}$-series into the $\\{\beta_{i}\\}$-series. In this subsection,
we give a simple introduction of PMCs for self-consistency.

If the pQCD approximant $\rho(Q)$ starts and ends at the $n_{\rm th}$ and
$(n+p-1)_{\rm th}$ order of $\alpha_{s}$, one has

$\displaystyle\rho(Q)=\sum_{k=1}^{p}c_{k}\alpha_{s}^{n+k-1}(\mu_{R}^{2}),$
(18)

where $\mu_{R}$ is an arbitrary initial choice of renormalization scale, $Q$
represents the kinetic scale at which the observable is measured or the
typical momentum flow of the reaction. The perturbative coefficients $c_{k}$
can be rewritten as the following forms by using the general QCD degeneracy
relations:

$\displaystyle c_{1}$ $\displaystyle=c_{1,0},$ (19) $\displaystyle c_{2}$
$\displaystyle=c_{2,0}+n\beta_{0}c_{2,1},$ (20) $\displaystyle c_{3}$
$\displaystyle=c_{3,0}+(n+1)\beta_{0}c_{3,1}+\frac{n(n+1)}{2}\beta_{0}^{2}c_{3,2}+n\beta_{1}c_{2,1},$
(21) $\displaystyle\cdots.$

Specifically, the fixed-order pQCD series $\rho$ can be further rewritten as
the following compact form:

$\displaystyle\rho(Q)=\sum_{k=1}^{p}\hat{c}_{k,0}\alpha_{s}^{n+k-1}(\mu_{R}^{2})+\sum_{k=1}^{p-1}\left[(n+k-1)\beta(\alpha_{s})\alpha_{s}^{n+k-2}(\mu_{R}^{2})\sum_{i=1}^{p-k}(-1)^{i}\Delta_{n,k}^{(i-1)}\sum_{j=0}^{i}C_{i}^{j}\hat{c}_{k+i-j,i-j}L_{\mu_{R}}^{j}\right],$
(22)

where $L_{\mu_{R}}\equiv\ln(\mu_{R}^{2}/Q^{2})$, and the combinatorial
coefficients are $C_{i}^{j}=j!/i!(j-i)!$. $\hat{c}_{i,0}=c_{i,0}$ is the
conformal coefficients, $\hat{c}_{i,j\geq 1}=c_{i,j}|_{\mu_{R}=Q}$ is non-
conformal coefficients, $c_{i,j}$ and $\hat{c}_{i,j}$ can be connected as the
relation $c_{i,j}=\sum_{k=0}^{j}C_{j}^{k}\hat{c}_{i-k,j-k}L_{\mu_{R}}^{k}$.
$\Delta_{n,k}^{(i-1)}$ is related to the running of the RGE, which is defined
as

$\displaystyle\Delta_{n,k}^{(i-1)}(\mu_{R})=\frac{1}{(n+k-1)\beta\alpha_{s}^{n+k-2}}\frac{1}{i!}\frac{{\rm
d}^{i}\alpha_{s}^{n+k-1}}{({\rm d}\ln\mu^{2})^{i}}\bigg{|}_{\mu=\mu_{R}},$
(23)

we thus have

$\displaystyle\Delta_{n,k}^{(0)}$ $\displaystyle=1,$ (24)
$\displaystyle\Delta_{n,k}^{(1)}$
$\displaystyle=-\frac{1}{2!}\sum_{i=0}(n+k+i)\beta_{i}\alpha_{s}^{i+1},$ (25)
$\displaystyle\Delta_{n,k}^{(2)}$
$\displaystyle=\frac{1}{3!}\sum_{j=0}\sum_{i=0}^{j}(n+k+i)(n+k+j+1)\beta_{i}\beta_{j-i}\alpha_{s}^{j+2},$
(26) $\displaystyle\cdots.$

According to the PMCs procedure, all non-conformal terms in Eq.(22) should be
absorbed into $\alpha_{s}(Q_{*}^{2})$, where $Q_{*}$ is a global effective
scale and its value can be determined by requiring all non-conformal terms
vanish in Eq.(22), i.e., one can obtain PMC scale $Q_{*}$ by solving the
following equation:

$\displaystyle\sum_{k=1}^{\ell+1}$
$\displaystyle\,\bigg{[}(n+k-1)\alpha_{s}^{k-1}(Q_{*}^{2})\sum_{i=1}^{\ell-k+2}(-1)^{i}\Delta_{n,k}^{(i-1)}(Q_{*})$
$\displaystyle\,\times\sum_{j=0}^{i}C_{i}^{j}\hat{c}_{k+i-j,i-j}L_{Q_{*}}^{j}\bigg{]}=0,$
(27)

where $\ell$ indicates that the PMC scale $Q_{*}$ can be fixed up to NℓLL
accuracy. Since the $\\{\beta_{i}\\}$-terms in the perturbative series (22)
have been resummed to determine the accurate value of $\alpha_{s}$, and the
resultant PMC series is free of divergent renormalon terms and becomes scheme-
and scale-independent pQCD conformal series [43], e.g.

$\displaystyle\rho|_{\rm
PMC}=\sum_{k=1}^{p}\hat{c}_{k,0}\alpha_{s}^{n+k-1}(Q_{*}^{2}).$ (28)

Practically, to achieve the analytic expression, the PMC scale $Q_{*}$ is
expressed as the form of $\ln\left(Q_{*}^{2}/Q^{2}\right)$ by expanding
Eq.(III) in terms of $\alpha_{s}(Q_{*}^{2})$ or $\alpha_{s}(Q^{2})$ up to the
required NℓLL accuracy [39, 40]. In this paper, we will directly solve
Eq.(III) numerically so as to derive more accurate PMC scales $Q_{*,\rm C}$
and $Q_{*,\rm NC}$ 222Due to double suppression from both the
$\alpha_{s}$-suppression and the exponential suppression, the usual analytic
treatment is accurate enough. For the present case, the NℓLL-level expansion
series of $\ln\left(Q_{*}^{2}/Q^{2}\right)$ does not convergent enough in the
threshold region; It is thus helpful to keep all known divergent $1/v$-terms
so as to achieve a more accurate prediction..

##  IV Numerical Results

We are now ready to deal with the Coulomb and non-Coulomb corrections within
the PMC framework.

For the Coulomb part $\mathcal{R}_{\rm C}$ (9), after applying the PMC, we
obtain the resultant PMC series, i.e.

$\displaystyle\mathcal{R}_{\rm C}|_{\rm PMC}$ $\displaystyle=$ $\displaystyle
1+\sum_{i}\hat{r}_{i,0}^{\rm C}\alpha_{s}^{{\rm V},i}(Q_{*,\rm C}^{2}),$ (29)
$\displaystyle=$ $\displaystyle 1+\frac{\pi}{2v}C_{F}\alpha_{s}^{\rm
V}(Q_{*,\rm C}^{2})+\frac{\pi^{2}}{12v^{2}}C_{F}^{2}\alpha_{s}^{\rm
V,2}(Q_{*,\rm C}^{2})$ $\displaystyle+{\cal O}(\alpha_{s}^{\rm V,4}),$

where the first three conformal coefficients $\hat{r}_{i,0}^{\rm C}$ can be
read from the coefficients $r_{i}^{\rm C}$ given in Eqs.(15,16,17) by applying
the degeneracy relations (19)-(21) with the help of the replacement $c_{i}\to
r_{i}^{\rm C}$. It is interesting to find that since $\hat{r}_{3,0}^{\rm
C}=0$, the N3LO Coulomb coefficient $r_{3}^{\rm C}$ (17) is exactly non-
conformal, which shall be resummed through the RGE and determines the
effective magnitude of $\alpha_{s}$, giving the effective momentum flow near
the threshold region (corresponding to the PMC scale $Q_{*,C}$); And the
resultant Coulomb series becomes a N2LO series, which can be further resummed
into a Sommerfeld-Gamow-Sakharov (SGS) factor so as to achieve a precise
prediction for the Coulomb part [11]:

$\mathcal{R}_{\rm C}|_{\rm PMC}^{\rm
resum}=\frac{X}{1-e^{-X}},\;\;X=\frac{\pi}{v}C_{F}\alpha_{s}^{\rm V}(Q_{*,\rm
C}^{2}).$ (30)

The $1/v$ in the numerator will be canceled with the $v$ in $\sigma_{0}$, thus
the Coulomb singularity near the threshold can be eliminated.

Similarly, for the non-Coulomb part (9), its PMC conformal series up to N3LO
can be achieved and be written in the following form:

$\displaystyle\mathcal{R}_{\rm NC}|_{\rm
PMC}=1+\sum_{i=1}^{3}\hat{r}_{i,0}^{\rm NC}\alpha_{s}^{i}(Q_{*,\rm NC}),$ (31)

where conformal coefficients $\hat{r}_{i,0}^{\rm NC}$ can be obtained in the
same way as that of $\hat{r}_{i,0}^{\rm C}$.

It is worth mentioning that, even though the Coulomb divergence has been
eliminated, $\ln v$-terms still remain in the $v^{0}$ terms of
$\mathcal{R}_{\rm NC}$. However, when compared to the terms with negative
powers of $v$, these terms are relatively more benign, and therefore, further
treatment of those terms is not addressed in this article.

###  IV.1 The cross section of $t\bar{t}$ production

To do the numerical calculation, we take $m_{c}=1.5~{}{\rm GeV}$,
$m_{b}=4.8~{}{\rm GeV}$, $m_{t}=173.34~{}{\rm GeV}$, $\alpha=1/132.2$ and
$\alpha_{s}(M_{Z})=0.1179$ as their central values.

Figure 1: (Color online) The PMC scales for the Coulomb and non-Coulomb
corrections. The blue dashed, red dashed, blue solid and red solid lines
represent the LL- and NLL- non-Coulomb and Coulomb scales, respectively.

According to the procedure mentioned in the previous section, we can
numerically obtain the PMC scales corresponding the Coulomb and non-Coulomb
corrections, as shown in Fig.1, the Coulomb PMC scale $Q_{*,\rm C}$ is more
sensitive to changes in $\sqrt{s}$ (or $v$).

Figure 2: Total cross section $\sigma_{t\bar{t}}$ with different QCD
corrections under conventional (dashed line) and PMCs (solid line) scale-
setting approaches, respectively.

We then present the total cross section before and after applying the PMCs in
Fig.2, where the middle lines correspond to the choice of $\mu_{R}=\sqrt{s}$
and the wide band correspond to the choice of
$\mu_{R}\in[\sqrt{s}/2,2\sqrt{s}]$ for conventional results. The application
of PMC single-scale setting approach removes the renormalization scale
dependence. Further, by further doing the resummation of $1/v$-power terms,
the divergency near the top-quark pair production threshold due to Coulomb-
like interaction has been successfully suppressed.

###  IV.2 Estimation of the N4LO QCD contributions using the Bayesian analysis

The resulting PMC series offers not only precise predictions for the known
fixed-order pQCD series but also provides a reliable foundation for assessing
the possible contributions from the UHO-terms, thus greatly enhancing the
predictive power of perturbation theory.

In the following, we will employ the Bayesian approach (BA) [54, 55, 56, 57]
to estimate the magnitude of the uncalculated N4LO-terms for the total cross
section $\sigma(s)$ by utilizing the known N3LO-level conventional and PMC
series, respectively.

The BA quantifies the contributions of the UHO-terms in terms of the
probability distribution. This method is powerful in constructing probability
distributions, wherein Bayes’ theorem is applied to iteratively update the
probabilities as new information becomes available. A comprehensive
introduction to the BA and its combination with the PMC approach can be found
in Ref.[70]. Generally, the conditional probability density function (p.d.f.)
$f_{r}(r_{p+1}|r_{1},\cdots,r_{p})$ for the uncalculated higher-order
coefficient $r_{p+1}$ with known coefficients $\\{r_{1},r_{2},\cdots,r_{p}\\}$
is given by

$\displaystyle f_{r}(r_{p+1}|r_{1},\cdots,r_{p})=\int
h_{0}(r_{p+1}|\bar{r})f_{\bar{r}}(r_{p+1}|r_{1},\cdots,r_{p})\,{\rm
d}\bar{r},$ (32)

where $f_{\bar{r}}(r_{p+1}|r_{1},\cdots,r_{p})$ is the conditional p.d.f. for
$\bar{r}$ with known coefficients $\\{r_{1},r_{2},\cdots,r_{p}\\}$ and
$h_{0}(r_{p+1}|\bar{r})$ is the conditional p.d.f of $r_{p+1}$ given
$\bar{r}$. Applying Bayes’ theorem, we have

$\displaystyle
f_{\bar{r}}(r_{p+1}|r_{1},\cdots,r_{p})=\frac{h(r_{1},\cdots,r_{p}|\bar{r})g_{0}(\bar{r})}{\int
h(r_{1},\cdots,r_{p}|\bar{r})g_{0}(\bar{r})\,{\rm d}\bar{r}},$ (33)

where $h(r_{1},\cdots,r_{p}|\bar{r})=\prod_{k=1}^{p}h_{0}(r_{k}|\bar{r})$ and
$g_{0}(\bar{r})$ are the likelihood function and the prior p.d.f. for
$\bar{r}$ corresponding to the CH model [54], respectively. We can then
calculate the degree-of-belief (DoB) that the value of $r_{p+1}$ belongs to a
certain credible interval (CI) by using following formula

$\displaystyle{\rm DoB}=\int_{-r_{p+1}^{(\rm DoB)}}^{r_{p+1}^{(\rm
DoB)}}f_{r}(r_{p+1}|r_{1},\cdots,r_{p})\,{\rm d}r_{p+1}.$ (34)

As a result, we obtain the expression for the boundary $r_{p+1}^{(\rm DoB)}$
of the symmetric smallest CI with fixed DoB for $r_{p+1}$ as

$\displaystyle r_{p+1}^{(\rm DoB)}=\left\\{\begin{aligned}
&\bar{r}_{(p)}\frac{p+1}{p}{\rm DoB},&{\rm DoB}&\leq\frac{p}{p+1}\\\
&\bar{r}_{(p)}\left[(p+1)(1-{\rm DoB})\right]^{-1/p},&{\rm
DoB}&\geq\frac{p}{p+1}\end{aligned}\right.$ (35)

where $\bar{r}_{(p)}={\rm max}\\{|r_{1}|,|r_{2}|,\cdots,|r_{p}|\\}$. Finally,
the estimated $r_{p+1}$ for fixed DoB is

$\displaystyle r_{p+1}\in\left[-r_{p+1}^{(\rm DoB)},r_{p+1}^{(\rm
DoB)}\right].$ (36)

In the subsequent calculation, ${\rm DoB}=95.5\%$ is taken to estimate the
contributions from the UHO-terms.

Figure 3: Total cross section $\sigma(s)$ under conventional (Dashed line) and
PMCs (Solid line) scale-setting approaches, respectively.

Since the Coulomb corrections have been resummed into a SGS factor, we only
need to estimate the contributions from the UHO-terms for the non-Coulomb QCD
series. Because the coefficients of the conventional pQCD series are known to
be scale-dependent at every order, the BA can only be applied after one
specifies the choices for the renormalization scale, thereby introducing extra
uncertainties for the BA. Conversely, the PMC conformal series is scale-
independent, providing a more reliable foundation for constraining predictions
regarding the UHO contributions. By taking into account the conventional scale
dependence for $\mu_{R}\in[\sqrt{s}/2,2\sqrt{s}]$ and the predicted magnitudes
of the UHO-terms, we present the total uncertainties caused by the UHO-terms
using the BA in Fig.3. The results confirm the importance of the PMC scale-
setting approach.

##  V Summary

In the paper, we have presented an improved analysis for the total cross-
section of $e^{+}e^{-}\to\gamma^{*}\to t\bar{t}$ near the threshold region up
to N3LO QCD corrections. Analytical expression for the N3LO-level Coulomb-
terms has been given by employing the PSLQ algorithm, which agree well with
the recently given numerical values in Ref.[19]. The PMCs approach has been
applied to eliminate the large renormalization scale uncertainties for the
Coulomb and non-Coulomb terms separately. After that, the Coulomb-terms have
been resummed as a Sommerfeld-Gamow-Sakharov factor, which great improves the
reliability and precision of total cross-section near the threshold region.
Our results confirm the correctness and the importance of the PMC procedures
for achieving scheme-and-scale invariant pQCD predictions, which can also be
ensured by the commensurate scale relations among different schemes [71, 72]
and agrees with the renormalization group invariance [73, 74, 78, 75, 76, 77].
The elimination of the uncertainty in setting the renormalization scale for
pQCD predictions using the PMCs, together with the reliable estimate for the
UHO contributions obtained using the Bayesian analysis, greatly increases the
precision of collider tests of the Standard Model and thus the sensitivity to
new phenomena.

Acknowledgements: This work was supported in part by the Chongqing Graduate
Research and Innovation Foundation under Grant No.CYB23011 and No.ydstd1912,
and by the Natural Science Foundation of China under Grant No.12175025 and
No.12347101.

## Appendix A

For the sake of completeness and convenience, we give the NLO and N2LO
coefficients $r_{1}$ and $r_{2}$ at the scale $\mu=\sqrt{s}$ in the following,

$\displaystyle r_{1}=$
$\displaystyle\,\frac{C_{F}}{\pi}\bigg{\\{}\frac{1+v^{2}}{v^{2}}\left[4{\rm
Li}_{2}\left(\frac{1-v}{1+v}\right)+2{\rm
Li}_{2}\left(-\frac{1-v}{1+v}\right)+\ln\left(\frac{1-v}{1+v}\right)\left(\ln\left(\frac{2}{1+v}\right)+2\ln\left(\frac{2v}{1+v}\right)\right)\right]$
$\displaystyle\qquad-2\left(\ln\left(\frac{2}{1+v}\right)+2\ln\left(\frac{2v}{1+v}\right)\right)-\frac{(1-v)(33-39v-17v^{2}+7v^{3})}{8v(3-v^{2})}\ln\left(\frac{1-v}{1+v}\right)+\frac{3(5-3v^{2})}{4(3-v^{2})}\bigg{\\}}$
(37) $\displaystyle=$
$\displaystyle\,\frac{C_{F}}{\pi}\left[\frac{\pi^{2}}{2v}-4+\frac{\pi^{2}}{2}v+\left(-\frac{86}{9}+8\ln
2+\frac{16}{3}\ln v\right)v^{2}+\mathcal{O}(v^{4})\right],$ (38)
$\displaystyle r_{2}=$
$\displaystyle\,C_{F}\bigg{\\{}\frac{\pi^{2}C_{F}}{12v^{2}}+\frac{1}{v}\left[-2C_{F}+C_{A}\left(\frac{31}{72}-\frac{11}{12}\ln
v\right)+T_{F}n_{l}\left(-\frac{5}{18}+\frac{1}{3}\ln v\right)\right]$
$\displaystyle\qquad+C_{F}\left(-\frac{35}{18}+\frac{39}{4\pi^{2}}+\frac{\pi^{2}}{6}+\frac{4}{3}\ln
2-\frac{1}{\pi^{2}}\zeta_{3}-\frac{2}{3}\ln
v\right)+C_{A}\left(\frac{179}{72}-\frac{151}{36\pi^{2}}-\left(\frac{8}{3}+\frac{22}{3\pi^{2}}\right)\ln
2-\frac{13}{2\pi^{2}}\zeta_{3}-\ln v\right)$
$\displaystyle\qquad+T_{F}n_{l}\left(\frac{11}{9\pi^{2}}+\frac{8}{3\pi^{2}}\ln
2\right)+T_{F}n_{h}\left(-\frac{4}{9}+\frac{44}{9\pi^{2}}\right)$
$\displaystyle\qquad+\left[\frac{839}{216}+\frac{173}{18}\ln
2+\frac{101}{36}\ln v+T_{F}n_{l}\left(-\frac{13}{18}+\frac{1}{3}\ln
v\right)+T_{F}n_{h}\left(\frac{1}{6}+\frac{1}{3}\ln
2\right)\right]v+\mathcal{O}(v^{2})\bigg{\\}},$ (39)

where ${\rm Li}_{n}(z)=\sum_{k=1}^{\infty}z^{k}/k^{n}$ is polylogarithm
function. As for the N3LO coefficient $r_{3}$, the analytic expression of its
coefficient $r_{3,+}$ with non-negative power terms of $v$ is very tedious and
hard to give its analytic form, so we directly adopt its numerical form given
by Ref.[19].

As for the N3LO-level coefficients $r_{3,v}$ and $r_{3,v^{2}}$ with negative
power terms $1/v$ and $1/v^{2}$ that are helpful for fixing the PMC scale of
the Coulomb part and for getting the SGS factor, we adopt the PSLQ algorithm
with the help of Mathematica package PolyLogTools [60] so as to derive their
analytic form, e.g.

$\displaystyle r_{3,v^{2}}=$
$\displaystyle\,C_{F}\left[-\frac{\pi}{3}C_{F}^{2}+C_{F}C_{A}\left(\frac{31\pi}{216}-\frac{11}{6\pi}\zeta_{3}-\frac{11\pi}{36}\ln
v\right)+C_{F}T_{F}n_{l}\left(-\frac{5\pi}{54}+\frac{2}{3\pi}\zeta_{3}+\frac{\pi}{9}\ln
v\right)\right],$ (40) $\displaystyle r_{3,v}=$
$\displaystyle\,C_{F}\bigg{\\{}C_{F}^{2}\left(\frac{39}{8\pi}-\frac{35\pi}{36}+\pi\ln
2-\frac{\zeta_{3}}{2\pi}-\frac{\pi}{3}\ln(2v)\right)+C_{A}^{2}\left(\frac{4343}{5184\pi}+\frac{175\pi}{432}-\frac{\pi^{3}}{128}+\frac{11}{48\pi}\zeta_{3}-\frac{247}{108\pi}\ln
v+\frac{121}{72\pi}\ln^{2}v\right)$
$\displaystyle\qquad+C_{A}C_{F}\left[-\frac{275}{72\pi}+\frac{179\pi}{144}-\left(\frac{4\pi}{3}+\frac{11}{3\pi}\right)\ln
2-\frac{13}{4\pi}\zeta_{3}+\left(\frac{11}{3\pi}-\frac{\pi}{2}\right)\ln
v\right]+C_{F}T_{F}n_{h}\left(\frac{22}{9\pi}-\frac{2\pi}{9}\right)$
$\displaystyle\qquad+T_{F}n_{l}\left[C_{F}\left(\frac{331}{288\pi}+\frac{\zeta_{3}}{2\pi}+\frac{4}{3\pi}\ln
2-\frac{13}{12\pi}\ln
v\right)+C_{A}\left(-\frac{899}{1296\pi}-\frac{11\pi}{54}-\frac{7\zeta_{3}}{12\pi}+\frac{217}{108\pi}\ln
v-\frac{11}{9\pi}\ln^{2}v\right)\right]$
$\displaystyle\qquad+T_{F}^{2}n_{l}^{2}\left[\frac{25}{162\pi}+\frac{\pi}{27}-\frac{10}{27\pi}\ln
v+\frac{2}{9\pi}\ln^{2}v\right]\bigg{\\}}.$ (41)

It should be pointed out that the correctness of the PSLQ algorithm can be
further tested by comparing the analytical results of $r_{1}$ and $r_{2}$
given by Ref.[11, 12] with the numerical results given in Ref.[19]; and by
comparing with the above $r_{3}$ with the semi-numerical expression of $r_{3}$
given in Refs.[13, 14].

## References

  * [1] J. Jersak, E. Laermann and P. M. Zerwas, “Electroweak Production of Heavy Quarks in e+ e- Annihilation,” Phys. Rev. D 25, 1218 (1982) [erratum: Phys. Rev. D 36, 310 (1987)].
  * [2] W. Beenakker, S. C. van der Marck and W. Hollik, “e+ e- annihilation into heavy fermion pairs at high-energy colliders,” Nucl. Phys. B 365, 24 (1991).
  * [3] K. G. Chetyrkin, J. H. Kuhn and M. Steinhauser, “Three loop polarization function and O (alpha-s**2) corrections to the production of heavy quarks,” Nucl. Phys. B 482, 213 (1996).
  * [4] J. Gao and H. X. Zhu, “Electroweak prodution of top-quark pairs in e+e- annihilation at NNLO in QCD: the vector contributions,” Phys. Rev. D 90, 114022 (2014).
  * [5] J. Gao and H. X. Zhu, “Top Quark Forward-Backward Asymmetry in $e^{+}e^{-}$ Annihilation at Next-to-Next-to-Leading Order in QCD,” Phys. Rev. Lett. 113, 262001 (2014).
  * [6] L. Chen, O. Dekkers, D. Heisler, W. Bernreuther and Z. G. Si, “Top-quark pair production at next-to-next-to-leading order QCD in electron positron collisions,” JHEP 12, 098 (2016).
  * [7] Z. Capatti, V. Hirschi and B. Ruijl, “Local unitarity: cutting raised propagators and localising renormalisation,” JHEP 10, 120 (2022).
  * [8] Y. Kiyo, A. Maier, P. Maierhofer and P. Marquard, “Reconstruction of heavy quark current correlators at O(alpha(s)**3),” Nucl. Phys. B 823, 269 (2009).
  * [9] M. Beneke, Y. Kiyo, P. Marquard, A. Penin, J. Piclum and M. Steinhauser, “Next-to-Next-to-Next-to-Leading Order QCD Prediction for the Top Antitop $S$-Wave Pair Production Cross Section Near Threshold in $e^{+}e^{-}$ Annihilation,” Phys. Rev. Lett. 115, 192001 (2015).
  * [10] M. Beneke, Y. Kiyo, A. Maier and J. Piclum, “Near-threshold production of heavy quarks with $\tt{QQbar\\_threshold}$,” Comput. Phys. Commun. 209, 96 (2016).
  * [11] A. Czarnecki and K. Melnikov, “Two loop QCD corrections to the heavy quark pair production cross-section in e+ e- annihilation near the threshold,” Phys. Rev. Lett. 80, 2531 (1998).
  * [12] R. N. Lee, A. V. Smirnov, V. A. Smirnov and M. Steinhauser, “Three-loop massive form factors: complete light-fermion corrections for the vector current,” JHEP 03, 136 (2018).
  * [13] M. Fael, F. Lange, K. Schönwald and M. Steinhauser, “Massive Vector Form Factors to Three Loops,” Phys. Rev. Lett. 128, 172003 (2022).
  * [14] M. Fael, F. Lange, K. Schönwald and M. Steinhauser, “Singlet and nonsinglet three-loop massive form factors,” Phys. Rev. D 106, 034029 (2022).
  * [15] R. Schwienhorst, D. Wackeroth, K. Agashe, S. Airen, S. Alioli, J. Aparisi, G. Bevilacqua, H. Y. Bi, R. Brock and A. G. Camacho, et al. “Report of the Topical Group on Top quark physics and heavy flavor production for Snowmass 2021,” [arXiv:2209.11267 [hep-ph]].
  * [16] H. Baer et al. [ILC], “The International Linear Collider Technical Design Report - Volume 2: Physics,” [arXiv:1306.6352 [hep-ph]].
  * [17] J. B. Guimarães da Costa et al. [CEPC Study Group], “CEPC Conceptual Design Report: Volume 2 - Physics & Detector,” [arXiv:1811.10545 [hep-ex]].
  * [18] A. Abada et al. [FCC], “FCC Physics Opportunities: Future Circular Collider Conceptual Design Report Volume 1,” Eur. Phys. J. C 79, 474 (2019).
  * [19] X. Chen, X. Guan, C. Q. He, X. Liu and Y. Q. Ma, “Heavy-quark-pair production at lepton colliders at NNNLO in QCD,” [arXiv:2209.14259 [hep-ph]].
  * [20] D. J. Gross and F. Wilczek, “Asymptotically Free Gauge Theories - I,” Phys. Rev. D 8, 3633 (1973).
  * [21] H. D. Politzer, “Asymptotic Freedom: An Approach to Strong Interactions,” Phys. Rept. 14, 129 (1974).
  * [22] W. E. Caswell, “Asymptotic Behavior of Nonabelian Gauge Theories to Two Loop Order,” Phys. Rev. Lett. 33, 244 (1974).
  * [23] O. V. Tarasov, A. A. Vladimirov and A. Y. Zharkov, “The Gell-Mann-Low Function of QCD in the Three Loop Approximation,” Phys. Lett. B 93, 429 (1980).
  * [24] S. A. Larin and J. A. M. Vermaseren, “The Three loop QCD Beta function and anomalous dimensions,” Phys. Lett. B 303, 334 (1993).
  * [25] T. van Ritbergen, J. A. M. Vermaseren and S. A. Larin, “The Four loop beta function in quantum chromodynamics,” Phys. Lett. B 400, 379 (1997).
  * [26] K. G. Chetyrkin, “Four-loop renormalization of QCD: Full set of renormalization constants and anomalous dimensions,” Nucl. Phys. B 710, 499 (2005).
  * [27] M. Czakon, “The Four-loop QCD beta-function and anomalous dimensions,” Nucl. Phys. B 710, 485 (2005).
  * [28] P. A. Baikov, K. G. Chetyrkin and J. H. Kühn, “Five-Loop Running of the QCD coupling constant,” Phys. Rev. Lett. 118, 082002 (2017).
  * [29] F. Herzog, B. Ruijl, T. Ueda, J. A. M. Vermaseren and A. Vogt, “The five-loop beta function of Yang-Mills theory with fermions,” JHEP 02, 090 (2017).
  * [30] S. J. Brodsky and X. G. Wu, “Scale Setting Using the Extended Renormalization Group and the Principle of Maximum Conformality: the QCD Coupling Constant at Four Loops,” Phys. Rev. D 85, 034038 (2012).
  * [31] S. J. Brodsky and X. G. Wu, “Eliminating the Renormalization Scale Ambiguity for Top-Pair Production Using the Principle of Maximum Conformality,” Phys. Rev. Lett. 109, 042002 (2012).
  * [32] S. J. Brodsky and L. Di Giustino, “Setting the Renormalization Scale in QCD: The Principle of Maximum Conformality,” Phys. Rev. D 86, 085026 (2012).
  * [33] M. Mojaza, S. J. Brodsky and X. G. Wu, “Systematic All-Orders Method to Eliminate Renormalization-Scale and Scheme Ambiguities in Perturbative QCD,” Phys. Rev. Lett. 110, 192001 (2013).
  * [34] S. J. Brodsky, M. Mojaza and X. G. Wu, “Systematic Scale-Setting to All Orders: The Principle of Maximum Conformality and Commensurate Scale Relations,” Phys. Rev. D 89, 014027 (2014).
  * [35] X. G. Wu, S. J. Brodsky and M. Mojaza, “The Renormalization Scale-Setting Problem in QCD,” Prog. Part. Nucl. Phys. 72, 44 (2013).
  * [36] X. G. Wu, Y. Ma, S. Q. Wang, H. B. Fu, H. H. Ma, S. J. Brodsky and M. Mojaza, “Renormalization Group Invariance and Optimal QCD Renormalization Scale-Setting,” Rept. Prog. Phys. 78, 126201 (2015).
  * [37] X. G. Wu, J. M. Shen, B. L. Du, X. D. Huang, S. Q. Wang and S. J. Brodsky, “The QCD Renormalization Group Equation and the Elimination of Fixed-Order Scheme-and-Scale Ambiguities Using the Principle of Maximum Conformality,” Prog. Part. Nucl. Phys. 108, 103706 (2019).
  * [38] L. Di Giustino, S. J. Brodsky, P. G. Ratcliffe, X. G. Wu and S. Q. Wang, arXiv:2307.03951 [hep-ph], accepted for publication in Prog. Part. Nucl. Phys.
  * [39] J. M. Shen, X. G. Wu, B. L. Du and S. J. Brodsky, “Novel All-Orders Single-Scale Approach to QCD Renormalization Scale-Setting,” Phys. Rev. D 95, 094006 (2017).
  * [40] J. Yan, Z. F. Wu, J. M. Shen and X. G. Wu, “Precise perturbative predictions from fixed-order calculations,” J. Phys. G 50, 045001 (2023).
  * [41] S. J. Brodsky, G. P. Lepage and P. B. Mackenzie, “On the Elimination of Scale Ambiguities in Perturbative Quantum Chromodynamics,” Phys. Rev. D 28, 228 (1983).
  * [42] S. J. Brodsky, L. Di Giustino, P. G. Ratcliffe, S. Q. Wang and X. G. Wu, [arXiv:2311.17360 [hep-ph]].
  * [43] X. G. Wu, J. M. Shen, B. L. Du and S. J. Brodsky, “Novel demonstration of the renormalization group invariance of the fixed-order predictions using the principle of maximum conformality and the $C$-scheme coupling,” Phys. Rev. D 97, 094030 (2018).
  * [44] E. C. G. Stueckelberg and A. Petermann, “Normalization of constants in the quanta theory,” Helv. Phys. Acta 26, 499 (1953).
  * [45] A. Peterman, “Renormalization Group and the Deep Structure of the Proton,” Phys. Rept. 53, 157 (1979).
  * [46] M. Gell-Mann and F. E. Low, “Quantum electrodynamics at small distances,” Phys. Rev. 95, 1300 (1954).
  * [47] N.N. Bogoliubov and D.V. Shirkov, “Application of the renormalization group to improve the formulae of perturbation theory,” Dok. Akad. Nauk SSSR 103, 391 (1955).
  * [48] S. J. Brodsky and X. G. Wu, “Self-Consistency Requirements of the Renormalization Group for Setting the Renormalization Scale,” Phys. Rev. D 86, 054018 (2012).
  * [49] K. Hagiwara, Y. Sumino and H. Yokoya, Phys. Lett. B 666, 71 (2008).
  * [50] Y. Kiyo, J. H. Kuhn, S. Moch, M. Steinhauser and P. Uwer, Eur. Phys. J. C 60, 375 (2009).
  * [51] S. J. Brodsky and X. G. Wu, Phys. Rev. D 86, 014021 (2012).
  * [52] S. Q. Wang, S. J. Brodsky, X. G. Wu, L. Di Giustino and J. M. Shen, Phys. Rev. D 102, 014005 (2020).
  * [53] A. B. Arbuzov and T. V. Kopylova, “On relativization of the Sommerfeld-Gamow-Sakharov factor,” JHEP 04, 009 (2012).
  * [54] M. Cacciari and N. Houdeau, “Meaningful characterisation of perturbative theoretical uncertainties,” JHEP 09, 039 (2011).
  * [55] E. Bagnaschi, M. Cacciari, A. Guffanti and L. Jenniches, “An extensive survey of the estimation of uncertainties from missing higher orders in perturbative calculations,” JHEP 02, 133 (2015).
  * [56] M. Bonvini, “Probabilistic definition of the perturbative theoretical uncertainty from missing higher orders,” Eur. Phys. J. C 80, 989 (2020).
  * [57] C. Duhr, A. Huss, A. Mazeliauskas and R. Szafron, “An analysis of Bayesian estimates for missing higher orders in perturbative calculations,” JHEP 09, 122 (2021).
  * [58] H. R. P. Ferguson, and D. H. Bailey, “A Polynomial Time, Numerically Stable Integer Relation Algorithm,” RNR Technical Report RNR-91-032 (1992).
  * [59] H. R. P. Ferguson, D. H. Bailey, and S. Arno, “Analysis of PSLQ, an Integer Relation Finding Algorithm,” Math. Comp. 68 (1999), 351.
  * [60] C. Duhr and F. Dulat, “PolyLogTools — polylogs for the masses,” JHEP 08, 135 (2019).
  * [61] T. Appelquist, M. Dine and I. J. Muzinich, The Static Potential in Quantum Chromodynamics, Phys. Lett. 69B, 231 (1977).
  * [62] W. Fischler, Quark - anti-Quark Potential in QCD, Nucl. Phys. B 129, 157 (1977).
  * [63] M. Peter, The Static quark - anti-quark potential in QCD to three loops, Phys. Rev. Lett. 78, 602 (1997).
  * [64] Y. Schroder, The Static potential in QCD to two loops, Phys. Lett. B 447, 321 (1999).
  * [65] A. V. Smirnov, V. A. Smirnov and M. Steinhauser, Fermionic contributions to the three-loop static potential, Phys. Lett. B 668, 293 (2008).
  * [66] A. V. Smirnov, V. A. Smirnov and M. Steinhauser, Three-loop static potential Phys. Rev. Lett. 104, 112002 (2010).
  * [67] C. Anzai, Y. Kiyo and Y. Sumino, Static QCD potential at three-loop order, Phys. Rev. Lett. 104, 112003 (2010).
  * [68] A. L. Kataev and V. S. Molokoedov, “The generalized Crewther relation and V-scheme: analytic $O(\alpha^{4}_{s})$ results in QCD and QED,” Teor. Mat. Fiz. 217, 44 (2023)
  * [69] H. Y. Bi, X. G. Wu, Y. Ma, H. H. Ma, S. J. Brodsky and M. Mojaza, “Degeneracy Relations in QCD and the Equivalence of Two Systematic All-Orders Methods for Setting the Renormalization Scale,” Phys. Lett. B 748, 13 (2015).
  * [70] J. M. Shen, Z. J. Zhou, S. Q. Wang, J. Yan, Z. F. Wu, X. G. Wu and S. J. Brodsky, “Extending the predictive power of perturbative QCD using the principle of maximum conformality and the Bayesian analysis,” Eur. Phys. J. C 83, 326 (2023)
  * [71] S. J. Brodsky and H. J. Lu, “Commensurate scale relations in quantum chromodynamics,” Phys. Rev. D 51, 3652 (1995).
  * [72] X. D. Huang, X. G. Wu, Q. Yu, X. C. Zheng, J. Zeng and J. M. Shen, “Generalized Crewther relation and a novel demonstration of the scheme independence of commensurate scale relations up to all orders,” Chin. Phys. C 45, 103104 (2021).
  * [73] A. Petermann, “La normalisation des constantes dans la therie des quantaNormalization of constants in the quanta theory,” Helv. Phys. Acta 26, 499 (1953).
  * [74] M. Gell-Mann and F. E. Low, “Quantum electrodynamics at small distances,” Phys. Rev. 95, 1300 (1954).
  * [75] N.N. Bogolyubov and D.V. Shirkov, “Group of Charge Renormalization in Quantum Field Theory”, Dokl. Akad. Nauk SSSR 103, 391 (1955).
  * [76] C. G. Callan, Jr., “Broken scale invariance in scalar field theory,” Phys. Rev. D 2, 1541 (1970).
  * [77] K. Symanzik, “Small distance behavior in field theory and power counting,” Commun. Math. Phys. 18, 227 (1970).
  * [78] A. Peterman, “Renormalization Group and the Deep Structure of the Proton,” Phys. Rept. 53, 157 (1979).