ML-based tunnel detection and tunneled application classification

Johan Mazel

Matthieu Saudrais

Antoine Hervieu

Authors     Refs.     1lYear
8cMM methods

1c90Decision Tree
1c90K-Nearest Neighbor
1c90Naive Bayes
1c90Support Vector
1c90Gaussian Mixture Model
1c90Markov Chains
1c90Neural Networks

Bernaille et al.       [1]     2006
-     -     -     -     -                 x

Bernaille et al.       [2]       2006
-     -     -     -     -           x                 x           x

Kim et al.             [3]          2008
-     -     -     -     -     x     x     x     x                       x

Jaber et al.          [4]             2011
-     -     -     -     -                                   x

classification inside tunnel

Wright et al.          [5]      2006
x           x                       x                             x

Bernaille et al.       [6]       2007
x                                                     x

Maiolini et al.        [7]         2009
x                                               x

Sun et al.             [8]             2010
x     x                                   x

Okada et al.           [9]     2011
x     x                 x

Okada et al.           [10]     2011
x     x     x           x     x

Korczyński et al.     [11]     2014
x                                                                 x

Kumano et al.          [12]        2014
x     x     x                 x

Tunnelling procotocols and machine learning techniques used in
previous work.

Type     Authors     Ref.     Year

4c1.4cmN first

(lr)5-8 (lr)11-14

90# packet

90Packet size

90TCP-based (eg. flag)
90TLS type

90Packet size
90Elapsed time
90Packet/byte burst

Bernaille et al.                             [1]

Bernaille et al.                             [2]

Kim et al.                                   [3]
x         {T,m,M,$\mu$,$\sigma$}/    x
{T,m,M,$\mu$,$\sigma$}/    x               x

Jabber et al.                                [4]

classification inside tunnel

Wrigh et al.                                 [5]

x/         x

Bernaille et al.                             [6]

Maiolini et al.                              [7]

x     x     x

Sun et al.                                   [8]
x         $\mu$,m,M

Okada et al.                                 [9]
x/    {T,m,M,,$\mu$,$\sigma^2$}/    x

Okada et al.                                 [10]
x/    {T,m,M,,$\mu$,$\sigma^2$}/    x

Korczyński et al.                            [11]

Kumano et al.                                [12]
x/    {T,m,M,,$\mu$,$\sigma^2$}/    x

Netflow v5                                   [13]

x         T,$\mu$                                x

Netflow v9                                   [14]

x/    T/,{m,M}/𝕀    x

Argus                                        [15]

x/    T/    x

Zeek (conn.log)                               [16]

x/    T/    x

Our work

x/    {T,m,M,$\mu$,$\sigma$}/    x
{m,M,$\mu$,$\sigma$}/                x     x/    x     x     x/

Network traffic features used for protocol classification and tunnel
protocol identification.
Abbreviations: x: used, : sign encodes packet direction
; $\mu$: mean, $\sigma$: standard deviation, $\sigma^2$: variance ; :
; 𝕀: incoming only ; : all direction
; T: total ; m: minimum ; M: maximum ; : $25^{th}$ centile, $50^{th}$
centile and $75^{th}$ centile.

Encrypted tunneling protocols are widely used.
Beyond business and personal uses, malicious actors also deploy
tunneling to hinder the detection of Command and Control and data exfiltration.
A common approach to maintain visibility on tunneling is to rely on network
traffic metadata and machine learning to analyze tunnel occurrence without
actually decrypting data.
Existing work that address tunneling protocols however exhibit several weaknesses:
their goal is to detect application inside tunnels and not tunnel identification,
they exhibit limited protocol coverage (e.g. OpenVPN and Wireguard are not
addressed), and both inconsistent features and diverse machine learning techniques
which makes performance comparison difficult.

Our work makes four contributions that address these limitations and provide
further analysis.
First, we address OpenVPN and Wireguard.
Second, we propose a complete pipeline to detect and classify tunneling
protocols and tunneled applications.
Third, we present a thorough analysis of the performance of both network traffic
metadata features and machine learning techniques.
Fourth, we provide a novel analysis of domain generalization regarding background
untunneled traffic, and, both domain generalization and adversarial learning regarding
Maximum Transmission Unit (MTU).

§ INTRODUCTION

Off-the-shelf encrypted tunneling protocols such as SSH or OpenVPN provide
a reliable way for users to protect their network traffic from passive analysis.
Companies often deploy such tools between locations or to access remote
General public also use tunneling to avoid passive monitoring or tampering
These same tools are also used by attackers to hide C2 (Command and Control)
and data exfiltration.
Although many tools already exists, such as IPsec, OpenVPN or SSH, new ones
also recently appeared (e.g. Wireguard).
The encrypted tunneling protocol landscape thus keep getting more complex
and makes malicious actors monitoring more difficult.

Port-based heuristic relying on tunnel default configuration are easy to
circumvent by attackers.
Supervised machine learning (ML) applied to network traffic metadata (usually
byte- and time-based) is thus a natural port-agnostic approach to identify some
protocols without actually inspecting encrypted data.
Several work use this approach to identify applications inside
tunneling tools such as SSL/TLS [5, 6, 8, 11],
both IPsec and PPTP [9, 10, 12],
and SSH [5, 7].
Their tunneling coverage is however limited: they do not address
common protocols such as OpenVPN or Wireguard.
These work also use non-overlapping network traffic metadata feature sets.
Furthermore, they use different machine learning techniques.
These last two aspects makes it difficult to compare their results.
Although machine learning is now easier to use than ever thanks to libraries
such as Scikit-learn [18], clear guidelines on how to use
machine learning to analyze tunneling protocols are missing.
Beyond application classification inside tunnel, up to our knowledge,
the preliminary task of identifying tunnels remains unaddressed.

Our goal is to address these limitations and design an approach to detect and classify
existing tunneling protocols.
We also want to provide a detailed analysis of network traffic metadata features
and machine learning techniques that provides clear and actionable guidelines for
real-world deployment.
Finally, we want to address domain generalization as it is paramount to assess
the whether a previously learned model will be usable beyond its initial training
We here address two aspects: background untunneled traffic and Maximum
Transmission Unit (MTU).

Our work provide four main contributions.
First, we include the two tunneling protocols that are not addressed in
existing work: OpenVPN and Wireguard.
Second, we design a complete pipeline to detect and classify tunneling
protocols and tunneled applications.
Third, we provide a thorough comparisons of network traffic metadata features
and machine learning techniques in the context of tunneling traffic detection
and classification.
Fourth, we address domain generalization regarding background untunneled
traffic, and then, both domain generalization and adversarial learning regarding MTU.
Up to our knowledge, both of these aspects have never been addressed before.

Our paper is structured as follows.
<Ref> present existing work on tunneling protocols analysis.
<Ref> details our methodology.
<Ref> presents our detection and classification results.

§ RELATED WORK

The idea of using network traffic metadata can be traced back to the original
work of Wright et al. [5] who applied it to application
classification inside SSL/TLS and SSH tunnels.
Several works extend this initial contribution regarding other
tunneling tools such as Tor [8] or IPsec and PPTP
[9, 10, 12].
Some works revisit previously addressed protocols but are relevant
considering the quick evolution of the tunneling and encryption landscapes
[6, 7, 11].

Another use case of network traffic metadata that generated many contributions
is protocol classification.
Bernaille et al. [1] is the seminal work for this use
Several other works then improve on this initial contribution
[2, 3, 4].

<Ref> describes the targeted tunneling tools and ML techniques used.
ML techniques and targeted tunneling protocols widely
vary across studies which make their results extremely difficult to compare.
These work also do not address tunneling tools that are common such OpenVPN or
<Ref> describes features used in existing work and tools.
They are very diverse.
The first papers on protocol classification inside tunnel used N first packet size
or IAT, while later work usually leverage five tuple flow-related feature such as
total number of packet or byte exchanged.
Features that were more recently proposed such as bursts [19]
have also not been evaluated in the context of tunneling protocol analysis.

Up to our knowledge, there is no existing work on the preliminary steps to
application classification inside tunnels which are tunnelling protocol
detection and/or classification.

Our goal is thus to 1) provide a complete pipeline for the analysis of tunelling
protocol (detection, classification, and application classification within tunnels),
and 2) address the diversity of existing work by providing a unified
point of view on performance of both network traffic features, and machine
learning techniques.

TODO: add velan paper

Processing pipeline.
Functional boxes have a purple background and rounded corners.
Data boxes have a white background and square corners.

Traffic generation topology.

1rTarget nb
1rMax application run
1rMax target nb

Web         Firefox/browsh     Alexa Top                     1M     20      1

wget        wget               Debian packages subset     42669      4      1

FTP get     lftp/vsftp         Debian packages subset      1000     20     10

FTP put     lftp/vsftp         Debian packages subset      1000     20     10

Deployed applications.
Max application run is the maximum number of times an application is run inside a tunnel.
Max target nb is the maximum number of target (website or debian package) retrieved by each application run.
In both cases, the actual number is randomly picked between 1 and this maximum value.
Timewise, browsing a website lasts between 30 and 45s and file downloads are consecutive.

Example of used N first feature.

UN-B [20]      66,585         27       -       -       -       -      66,612

UN-S [20]           0     27,978       -       -       -       -      27,978

UPC [21, 22]
525,537      4,093       -       -       -       -     529,630

This work             241,941       2,400     2,400     2,400     2,400     2,400     253,941

Data used in our work.
UN-B (resp. UN-S) is the UNIBS dataset with untunneled (resp. SSH) traffic.

4emScikit learn function
4emShort name
Parameter name
5emParameter default value
Parameter values for grid search

AdaBoost                 [23]     AB
learning rate        1.0     $10^{-3,-2,-1,0,1}$

Gradient boosting        [24]     GB
max depth              3     1, 2, 3, 4, 5, 6, 7, 10, 20, 50

2*7emGaussian naive Bayes
[25]     GNB
-                       -     -

K nearest neighbors      [26]     KNN     n neighbors             5     1, 5, 10, 50

2*7emSupport vector with liblinear
[27]     SV LL
C                     1.0     $10^{-1,0,1,2,3,4,5,6,7}$

Logistic regression      [28]     LR
solver               saga     saga

C                     1.0     $10^{-1,0,1,2,3,4,5,6,7}$

2*7emLogistic regression with SGD
[29]     LR SGD
loss                  log     log

alpha                $10^{-4}$     $10^{-8,-7,-6,-5,-4,-3,-2,-1,0,1}$

2*7emSupport vector with SGD
[29]     SV SGD
loss                 hinge     hinge

alpha                $10^{-4}$     $10^{-8,-7,-6,-5,-4,-3,-2,-1,0,1}$

Random forest and        [30]     RF
criterion            gini     gini, entropy

decision tree            [31]     DT
min samples split       2     2,3,4,5,10,50,100

Machine learning algorithms used and parameters explored.
Random forest and decision tree both use the criterion and min samples split parameters.
Short names are used to identify algorithms in the rest of the paper.

§ METHODOLOGY

In this section, we present our approach.
First, we address the general machine learning-based pipeline in <Ref>,
and then, we detail the main steps: traffic generation in <Ref>,
used datasets in <Ref>,
feature extraction in <Ref>, and, machine learning in

§.§ General pipeline

<Ref> picture the whole pipeline that we use to analyze tunneling
First, we generate network traffic (see <Ref>) to
complete the external datasets that we obtained (see <Ref>).
Then, using generated PCAP files, we extract network traffic metadata feature
(see <Ref>).
The next step is the use of ML to determine if a given flow is a tunnel.
This task is named tunnel detection in the remainder of the paper.
If it is the case, we then identify which tunnel is present.
This step is called tunnel classification.
If a tunnel is present, we classify the application used inside the tunnel.
We call this phase application classification inside tunnel.

§.§ Traffic generation

Although some existing dataset provide some tunneling protocol network traffic
(eg. SSH in UNIBS [20]), we could not find any data for IPsec, OpenVPN
or Wireguard.
We thus design a generic topology to collect tunnel traffic.
Here, a host uses a tunnel to access some resources on Internet or on an
FTP file server.
This topology is presented on <Ref>.
The tunnel is setup between two hosts (tunnel start and tunnel end) and is
routed through a machine called gateway.
This topology is setup using Vagrant with the libvirt provider.
All hosts use Debian Testing with the version 20210228.1
(this is equivalent to Debian Bullseye).
We deactivate NIC offloading on all machines in order to avoid packet size
alteration during capture.

We deploy each tunnel between the tunnel start and tunnel end machines.
We setup the following tunnels: IPsec, OpenVPN TCP/UDP, SSH and Wireguard.

We use several application inside each tunnel: web browsing, FTP get/put and wget.
We actually use a single application inside each tunnel instance.
<Ref> provides additional details regarding applications used
inside tunnels.

§.§ Datasets

<Ref> present the data we use in this work.
The UNIBS datasets [20] were collected in 2009 and contain both SSH
tunneled and untunneled traffic.
The UPC dataset [22] was gathered in 2014 and contains mostly untunneled
traffic and SSH and Tor tunnels.
In this work, we only use the untunneled and SSH tunneled traffic.
The network traffic generation in <Ref> create 100 flows
for six MTU values (see <Ref> for additional details on MTU
choice) and four applications used with 5 tunnels and without tunnel.
We thus obtain 2400 flows for each tunneling protocol.
Untunneled traffic does not yield the same number of flows for each MTU due
to varying interactions with third-party during website browsing.

§.§ Feature extraction

Once network traffic has been generated, we extract five tuple flow (source and
destination IP, transport protocol, and source and destination port) and their
associated features.
We extract the following features about each five tuple flow: transport protocol
(e.g. TCP, UDP), packet number, duration, and,
mean, minimum, maximum, standard deviation in both direction and the total of
both packet size and Inter-Arrival Times (IAT).
We also extract the following N first feature (see example in <Ref>):
elapsed time since the start of the flow, IAT, packet direction, packet size,
packet size with direction in sign (called packet size with direction in the
remainder of the paper), packet burst and byte burst.
Burst are build by grouping consecutive packet in the same direction.
Packet bursts are the grouping sizes and byte bursts are the grouping total byte
We also build packet size from source and destination only, and byte burst from
source and destination only.
The <Ref> presents our features and those of previous work.

TODO: byte number is total of payload $\rightarrow$ replace by unique bytes, ie drop retransmission

§.§ Machine learning

add XGBoost

We use algorithms available in scikit-learn [18] for multiclass
classification: random forest, decision tree, AdaBoost, logistic regression
and support vector with
SGD, logistic regression with the saga solver, Gaussian naive Bayes, and K
nearest neighbors.
Some algorithms are omitted in the following experiments due to excessive
running time (e.g. K nearest neighbors for tunnel detection).

In terms of performance metric, we use the F1 score for tunnel detection and
the average F1 score across classes (macro in scikit learn) for tunnel
classification and application classification inside tunnels.

Performance comparison between algorithms is performed using a nested
cross-validation [32].
The inner loop uses a grid search to find the best parameters for each
algorithm (see <Ref>).
The best parameters are determined using F1 score as specified above.
The outer loop then reuse these optimal parameters to provide a performance lower
bound for comparison with other algorithms.

50 first packet sizes     50 first packet directions     50 first packet sizes with direction

50 first packet bursts     50 first byte bursts     50 first IAT

Tunnel detection using 50 first feature.
Shaded areas below and above the curve represent a confidence interval with 99% confidence level.

50 first packet sizes     50 first packet directions     50 first packet sizes with direction

50 first packet bursts     50 first byte bursts     50 first IAT

Tunnel classification using 50 first feature.
Shaded areas below and above the curve represent a confidence interval with 99% confidence level.

50 first packet sizes     50 first packet directions     50 first packet sizes with direction

50 first packet bursts     50 first byte bursts     50 first IAT

Application classification inside SSH using 150 first features.
Shaded areas below and above the curve represent a confidence interval with 99% confidence level.

§ RESULTS

In this section, we present the performance evaluation for all three steps
defined in <Ref>.
We first present a study of the impact of N where N is the number N first
feature used, e.g. packet size (see <Ref> and
We then compare N first features inside feature families (here packet direction/size-related
and burst-related) in <Ref>.
Next, we reuse the best N first features inside each family and compare their
performance in <Ref>.
We then provide a performance evaluation of these best features used
together along with feature importance analysis (see <Ref>).
Finally, we address domain generalization regarding change of background network
traffic in <Ref>, and both domain generalization and adversarial
learning regarding MTU in <Ref>.

We only present the results of application classification inside tunnel for the
SSH tunnel due to the lack of space.
This tunnel however exhibits the worst performance across all experiments.
Presented results thus always provide a performance lower bound for application
classification inside tunnels.

In the next sections, we often use a subset of all extracted features (see
<Ref>) to compare them.
We however always use the one-hot encoded transport protocol (TCP or UDP) as a
additional feature to any of these feature subsets.

Packet direction     Packet size     9emPacket size from src
9emPacket size from dst
10emPacket size with direction

90 Tunnel detection

90 Tunnel classification

909emApplication classification in SSH

Comparison of packet size, packet size from source, packet size from destination, packet direction and packet size with direction for all pipeline steps.
Error bars represent a confidence interval with 99% confidence level.

Packet burst     Byte burst     9emByte burst from src     9emByte burst from dst

90 Tunnel detection

90 Tunnel classification

909emApplication classification in SSH

Comparison of packet burst, byte burst, byte burst from source, and byte burst from destination for all pipeline steps.
Error bars represent a confidence interval with 99% confidence level.

§.§ N first feature

TODO: compare feature importance with F1 score jump on nfirst

In this section we analyze the impact of N on performance, where is N the number of N first
feature used.
N first features are: packet direction, packet size, packet size with direction
encoded in sign, packet burst, byte burst, and inter-arrival time.
We do not present elapsed time since flow start due to the lack of space and
because its performances are always lower than IAT.
check this again
Due to the computational cost of testing mutliple values of N, we do not
perform a grid search to find optimal parameters inside the inner loop, and
simply use default parameters from <Ref>.
We address all pipeline steps: tunnel detection, tunnel classification, and
application classification.

TODO: add text of pso/psr bbo/bbr when data is ready

TODO: split netflow v5/v9 into base and mod (with bonus)

§.§.§ Tunnel detection and classification

<Ref> (resp. <Ref>) presents the F1 score obtained with N values
between 1 and 50 for tunnel detection (resp. classification).
For both tasks, we use N between 1 and 50 because N values above 20 do
not provide significant performance improvement for the best algorithm (random forest).

Algorithms exhibits widely varying performance.
Random forest is the best one.
When N increases, decision tree and random forest exhibit increasing performance
for all features.
Other algorithms' performances sometimes decrease when N increases (e.g. N
between 20 and 25 with packet direction and Gaussian naive Bayes for tunnel
detection on <Ref>).

We then compare feature performance using the F1 score for N equal to 50 with
the best performing algorithms (decision tree and random forest).
Packet size, packet size with direction, byte burst and IAT provide very good
overall results.
Packet direction and burst exhibit slightly lower performance.

In terms of smallest N to reach maximum performance, ie smallest possible
memory usage to reach optimal performance, byte burst is the best
feature, followed by packet size with direction and IAT, then packet size,
and finally packet burst and packet direction.
Byte burst is the most efficient feature in terms of memory usage.
It however does not mean that byte burst is quick time-wise to reach its best
performance level because burst position in a flow is not directly linked
to packet position.
Indeed, the fourth byte burst may actually be located at the 30th packet which
would make byte burst slower than packet size with direction.

When one compares tunnel detection and tunnel classification, the smallest N
values to reach maximum performance are similar.
Algorithm performances are however better and vary less across algorithms for
tunnel classification than for tunnel detection.

§.§.§ Application classification

<Ref> presents the F1 score obtained with N values
between 1 and 150 for application classification inside SSH tunnels.
We use N between 1 and 150 because N values above 120 do not provides
significant performance improvement for the best algorithm (random forest).
We do not present result for other tunnels due to the lack of space.
TODO: add test about TR ;TR????

Results are similar to the tunnel detection and classification use cases.
Algorithm performance are diverse.
Random forest exhibits the best F1 scores, and its performances increase
when N increases for all features.
This is not the case for other algorithms (see N between 100 and 150 with
packet direction and Gaussian naive Bayes).
Overall, byte burst is the best feature.
Packet size, packet size with direction, packet burst are close second.
Packet direction exhibits the worst performance.
In terms of smallest N to reach maximum performance, byte burst is the best
feature, closely followed by packet burst, then packet direction, packet size
with direction, packet direction and IAT, in that order.

§.§.§ Summary

In this preliminary comparison, we determine that using 50 (resp. 150) features
is enough for tunnel detection and tunnel classification (resp. application
Overall, random forest is the best algorithm.
Byte burst is the best feature in terms of smallest N to reach the best

10emPacket size with direction     Byte burst     IAT     Netflow v5     Netflow v9

90 Tunnel detection

90 Tunnel classification

909emApplication classification in SSH

Comparison of packet size with direction, byte burst and IAT for all pipeline steps.
Error bars represent a confidence interval with 99% confidence level.
TODO: fix figure for tunnel detection and Netflow v9

§.§ Features

TODO: add comparison with different payload size?

In this section, we compare features: first, among feature families (e.g. packet
size vs packet size with direction or packet burst vs byte burst), then across
families using the best features from each families (e.g. packet size with
direction vs byte burst).

§.§.§ Preliminary feature comparison

We here analyze performance inside each feature families: packet-related
(with packet size, packet direction, packet size from source, packet size
from destination, packet size with direction in sign), and burst-related
(packet burst, byte burst, byte burst from source, byte burst from destination).

Packet direction, packet size from source, from destination, from both source and destination, and packet size with direction

<Ref> pictures the F1 scores obtained
with ML algorithms for packet direction and feature related to packet size.
F1 scores obtained with packet size from source are always better than packet
size from destination.
Packet size from source is also always better than packet size, except for
tunnel detection with algorithms such as logistic regression, support vector
and AdaBoost.
Packet direction is almost always worse than packet size from source, except for
application classification in SSH with AdaBoost.
Packet size with direction and packet size from source exhibits very similar
In the remainder of the paper, we arbitrarily choose to only consider packet size
with direction.

fix missing value for app classif/adaboost/packet direction

Packet burst and byte burst from source/destination/both

<Ref> pictures the F1 scores obtained
with ML algorithms for features related to bursts.
F1 scores obtained with packet burst are always worse than byte burst
(except for SSH application classification with Gaussian naive Bayes).
Byte burst from source and destination does not yield better results than byte burst.

Byte burst bigram and trigram

We now explore the use of bigram and trigram on byte burst values.
We do not use a figure due to the lack space but one is provided in appendix
(see <Ref>).

Byte burst bigram used alone and in combination with byte burst does not
improve results, except for IPsec-ESP with Gaussian Naive Bayes and OpenVPN TCP
with AdaBoost and Decision Tree.
Byte burst trigram used alone and in combination) does not improve results
(except for IPsec-ESP with all methods and OpenVPN TCP with AdaBoost and Decision
Considering the limited performance improvements and costly feature generation,
we do not use byte burst bigrams or trigrams in the remainder of the paper.

remove figure and keep text?

§.§.§ Feature comparison

Next we compare features selected from the first two preliminary comparisons,
packet size with direction and byte burst, with others features such as N first
IAT and Netflow features (v5 and v9) on <Ref>.
Netflow features always yield worse results than any other features.
IAT are worse than packet size with direction and byte burst, except for tunnel
classification and some ML algorithms such as both logistic regressions.
We do not present results for elapsed time since the start of the flow because its
results are worse than IAT.
check this again
Previous work [33] show that tampering with IAT and elapsed time
since the flow start
using a proxy yields a performance decrease.
We thus only consider packet size with direction and byte burst as relevant
features for the remainder of the paper.

§.§.§ Summary

Overall, packet size with direction and byte burst offer the best performance.
We emphasize that Netflow v5 and v9 feature sets do not provide good performances.

§.§ General results

[Tunnel detection]

[Tunnel classification]

Tunnel detection and tunnel classification using selected features (open world).
Error bars represent a confidence interval with 99% confidence level.

20 most important feature regarding Mean Decreased Impurity (MDI)
using Random Forest with 100 trees for tunnel detection.
Error bars represent a confidence interval with 99% confidence level.

20 most important feature regarding Mean Decreased Impurity (MDI)
using Random Forest with 100 trees for tunnel classification.
Error bars represent a confidence interval with 99% confidence level.

[OpenVPN TCP]

[OpenVPN UDP]

[in SSH]

[in Wireguard]

Application classification detection using selected features.
Error bars represent a confidence interval with 99% confidence level.

TODO: add feature importance for fs2

We now provide a performance evaluation of packet size with direction and byte
burst used together.
<Ref> pictures the F1 scores obtained for tunnel detection
and classification.
TODO: complete after update for ow
We observe very good performance for both tunnel detection and tunnel
classification using the best algorithms (e.g. random forest or decision tree).
TODO: check 15 and 20 features
<Ref> displays the 10 most feature importance based on Mean Decreased
Impurity (MDI) for tunnel detection.
Overall, IP protocol field-related one-hot-encoded features are not relevant.
This is consistent with the fact that tunnels here use both TCP and UDP.
16 out of the 20 most important feature are related to byte burst.
Byte bursts are thus here more important than packet sizes with direction.
<Ref> displays the 10 most feature importance based on Mean Decreased
Impurity (MDI) for tunnel classification.
IP protocol field-related one-hot-encoded features are here the first and fourth
most important features.
Beyond these two attributes, packet size with direction and byte burst are similarly important.

<Ref> pictures the F1 scores obtained for tunnel detection
and classification.
The most difficult tunneling protocol to classify application inside is SSH,
followed by both OpenVPN tunnels, and finally, IPsec and Wireguard.
Similarly to <Ref>, using the best algorithms, random forest
or decision tree, provides very good performances.

add text for <Ref>

§.§ Domain generalization regarding untunneled traffic

tableDomain generalization performance between the two
background traffic datasets regarding F1 score using random forest for
tunnel detection.

Feature     UNIBS     UPC

Netflow v5                  0.28±0.00

Netflow v9                  0.37±0.02

Byte burst                  0.61±0.03

Packet size w/ dir.         0.75±0.02

Netflow v5              0.54±0.01

Netflow v9              0.74±0.01

Byte burst              0.98±0.00

Packet size w/ dir.     0.99±0.00

[Tunnel detection]

[Tunnel classification]

[Application classification in SSH]

Random forest performance for training using MTU equals to 1500 and several MTU values for testing data.
Error bars represent a confidence interval with 99% confidence level.

[Tunnel detection]

[Tunnel classification]

[Application classification in SSH]

Packet sizes with direction in byte of the 10 most important feature regarding
MDI using Random Forest with 100 trees for all pipeline steps.
Boxplots encode median as the line inside the box, first and third quartile
as upper and lower limit of the box, whiskers as 1st and 99th quantiles, and
flier as points outside of these quantiles.

[Tunnel detection]

[Tunnel classification]

[Application classification in SSH]

Random forest performance for training using MTU equals to 1500 and several MTU values for testing data.
Error bars represent a confidence interval with 99% confidence level.

We here address the impact of the untunneled traffic on tunnel detection
We here only consider tunnel detection because untunneled traffic is
only used in this step of the ML pipeline (see <Ref>).
To this end, we use the UNIBS dataset for training and the UPC dataset for
testing, and then, the other way around.
We here do not use tunneling protocol-related flows from the UNIBS and UPC
<Ref> present F1 scores for packet sizes with direction,
byte bursts, and Netflow v5 and Netflow v9 features.
Although we previously show in  <Ref> that Netflow features
exhibit poor performances, we still include them here as classic flow-based features.
When we train a random forest model on the UPC (resp. UNIBS) dataset and test
it on the UNIBS (resp. UPC) dataset using Netflow v5, we obtain a 0.54 (0.28) F1 score.
These results are very poor, especially compared to what we previously obtained in
Netflow v9 exhibit the second worst performance decrease, followed by byte burst
and then, packet size with direction.
We hypothesize that training with UPC yield better results than training with
UNIBS because of two reasons.
First, the number of instances in the UPC dataset is much bigger than in the
UNIBS one (see <Ref>) which helps generalization of models
trained with UPC data.
Second, UPC dataset is more recent than UNIBS (2014 vs 2009).
This means that traffic categories present in UNIBS are probably also in UPC while
new protocols, applications, or usages that appeared between 2009 and 2014 may be
present in UPC but not in UNIBS.
Although we do not actually try to detect these new applications,
we hypothesize that their presence inside training data may be enough to modify
learned model.

The observed performance decrease is consistent with existing work targeting
traffic classification [34].
This section shows that training data must be as diverse as possible to ensure
a limited performance decrease when the trained model is deployed in a new

§.§ Domain generalization and adversarial learning regarding MTU

The Maximum Transmission Unit (MTU) is the maximum amount of data that can be
sent on a physical medium.
Its value is 1500 for Ethernet and much of Internet [35].
This value is configurable by OSes on network interfaces.
Reducing MTU decreases the size of sent and received packets.
This in turns impacts both the number sent and received packets, and the total
amount of data sent (because of the increasing overhead).
We here address the impact of MTU modification on our three pipeline steps.
We first consider common MTU used on Internet beyond the default value of 1500
which are, according to [35], 1472 and 1420.
Then, we also consider an attacker that leverages a purposely lowered MTU
to generate adversarial samples.
We here use MTU values of 1400, 1300, and 1200.
We thus obtain six MTU values (1500, 1472, 1420, 1400, 1300 and 1200) that we use to generate
traffic (see <Ref>) for the following experiments.
We here do not use the SSH flows from UNIBS because we do not know which MTU was used.

explain specific k-fold cross validation

<Ref> presents our results for all three pipeline steps,
packet size with direction, byte burst and Netflow v5 and v9 features when a
1500 bytes MTU is used for training.
Results for MTU testing equal to 1500 are different than
<Ref> because UNIBS SSH is not used.
Similarly to <Ref>, we include Netflow feature as classic flow-based features.
Performances modifications regarding testing MTU change across pipeline steps
and features: tunnel classification performances only decrease for Netflow v5,
tunnel detection is impacted with both Netflow feature sets,
while application classification's F1 scores deteriorate for all features
except packet sizes with direction.
MTU lowering reduces the amount of data sent in a single packet.
This, in turn, increases the packet number in a flow.
Netflow v5 and v9 feature sets are impacted because they contains packet
number-related features, either for all packet of flow or a subset (e.g. only
from the source in Netflow v9).
We observe that byte burst performances decrease as the testing MTU is smaller
for application classification, while this does not happen with packet sizes
with direction.
We hypothesize that byte burst is actually able to look further in time into
the flow, and is thus more sensitive to traffic changes caused by MTU reduction.
TODO: dig into actual position of byte burst in terms of packets using packet bursts
Contrary to byte burst, we observe that performances for packet sizes with
direction do not deteriorate as testing MTU decreases.
<Ref> displays boxplots and violinplots of packet
sizes with direction which exhibit the ten biggest MDI feature importance when
random forest is used.
It shows that most feature values are far from the testing MTU.
We thus hypothesize that packet size with direction performances do not
decrease because feature values are not modified when testing MTU changes.

We do not display the results obtained with a training MTU of 1420 and 1470
due to the lack of space.
The results are similar to what is obtained with a training MTU of 1500:
best results are obtained when testing MTU is equal to the training one, and
performances deteriorate otherwise.

§ DISCUSSIONS AND FUTURE WORK

We identify two main limitations to our work.
The first limitation is related to tunneling protocols and application coverage.
Although our tunneling protocol set is more complete than existing work, we
did not address tools such as PPTP.
We also did not deploy DNS-based tunelling which is commonly used for obfuscation.
Regarding applications, we use web browsing and file transfer as applications
inside tunnels.
Common usage such as command line is not included in our work.
We want to improve both deployed tunneling and application in future work.
The other limitation is that we only use a single application inside each tunnel.
Users however do not use tunelling this way, and often use several applications
inside a tunnel.
Our results regarding tunnel detection and classification using N first features
such as packet size or byte burst are however not impacted by this limitation
because they need a very small of features N to perform well (see

Many recent on network traffic metadata analysis work use deep learning-related
We however achieve very good performance with classical ML algorithms such as
random forest or decision tree and thus decided to not use deep learning at all.
TODO: add motivation on runtime

We also plan to extend this work along two other axes.
First, we use each tunelling protocol with a single cipher suite.
We want to extend this work by exploring other cipher suites (e.g. a different
block size) which may alter network traffic metadata.
Second, we plan to investigate prediction intervals for both aleatoric and epistemic
uncertainties using [36].

We think that our contribution on domain generalization and adversarial learning
regarding MTU is very important for all use cases that leverage network traffic
Regarding security and privacy, we hypothesize that use cases such as malware C2
detection, website fingerprinting, DoH fingerprinting, Tor detection are impacted by MTU

§ CONCLUSION

In this work, we design a methodology to detect tunneling protocols, classify
them once detected, and identify application used inside tunnels.
We provide a thorough analysis of both features and machine learning algorithms.
We also address domain generalization evaluation regarding the untunneled traffic
for tunnel detection, and domain generalization and adversarial
learning analysis regarding MTU for all three pipeline steps.

§ APPENDIX

§.§ Learning curve

[Tunnel detection]

[Tunnel classification]

[Application classification in SSH]

Random forest performance for several training dataset sizes.
Error bars represent a confidence interval with 99% confidence level.

<Ref> pictures the F1 scores obtained with random forest for
several trainign dataset sizes.
For all scenario, performance increase more important for smaller increase of
training data size than
This performance increase is negligible for the biggest increase in training
data size (cf point to the right of each figure and log scale).

§.§ Netflow

Netflow v5 base     Netflow v5 ext.     Netflow v9 base     Netflow v9 ext.

90 Tunnel detection

90 Tunnel classification

909emApplication classification in SSH

Comparison of packet size with direction, byte burst and IAT for all pipeline steps.
Error bars represent a confidence interval with 99% confidence level.
TODO: fix figure for tunnel detection and Netflow v9

§.§ Byte burst bigram and trigram

Byte burst     Byte burst bigram only     Byte burst with bigram     Byte burst trigram only     Byte burst with trigram

90 Tunnel detection

90 Tunnel classification

909emApplication classification in SSH

Comparison byte burst, byte burst bigram only, byte burst with
bigram, byte burst trigram only, and byte burst with trigram for all pipeline steps.
Error bars represent a confidence interval with 99% confidence level.

<Ref> pictures the F1 scores obtained
with ML algorithms for features from the packet burst family.

[1]
L. Bernaille, R. Teixeira, I. Akodkenou, A. Soule, and K. Salamatian, “Traffic
classification on the fly,” ACM SIGCOMM Computer Communication
Review, vol. 36, no. 2, pp. 23–26, 2006.

[2]
L. Bernaille, R. Teixeira, and K. Salamatian, “Early application
identification,” in Proceedings of the 2006 ACM CoNEXT
conference.1em plus 0.5em minus 0.4emACM, 2006, p. 6.

[3]
H. Kim, K. C. Claffy, M. Fomenkov, D. Barman, M. Faloutsos, and K. Lee,
“Internet traffic classification demystified: myths, caveats, and the best
practices,” in Proceedings of the 2008 ACM CoNEXT conference.1em plus 0.5em minus 0.4emACM, 2008, p. 11.

[4]
M. Jaber, R. G. Cascella, and C. Barakat, “Can we trust the inter-packet time
for traffic classification?” in Communications (ICC), 2011 IEEE
International Conference on.1em plus 0.5em minus 0.4emIEEE,
2011, pp. 1–5.

[5]
C. V. Wright, F. Monrose, and G. M. Masson, “On inferring application protocol
behaviors in encrypted network traffic,” Journal of Machine Learning
Research, vol. 7, no. Dec, pp. 2745–2769, 2006.

[6]
L. Bernaille and R. Teixeira, “Early recognition of encrypted applications,”
Passive and Active Network Measurement, pp. 165–175, 2007.

[7]
G. Maiolini, A. Baiocchi, A. Iacovazzi, and A. Rizzi, “Real time
identification of ssh encrypted application flows by using cluster analysis
techniques,” NETWORKING 2009, pp. 182–194, 2009.

[8]
G.-L. Sun, Y. Xue, Y. Dong, D. Wang, and C. Li, “An novel hybrid method for
effectively classifying encrypted traffic,” in Global
Telecommunications Conference (GLOBECOM 2010), 2010 IEEE.1em plus
0.5em minus 0.4emIEEE, 2010, pp. 1–5.

[9]
Y. Okada, S. Ata, N. Nakamura, Y. Nakahira, and I. Oka, “Application
identification from encrypted traffic based on characteristic changes by
encryption,” in Communications Quality and Reliability (CQR), 2011
IEEE International Workshop Technical Committee on.1em plus 0.5em
minus 0.4emIEEE, 2011, pp. 1–6.

[10]
——, “Comparisons of machine learning algorithms for application
identification of encrypted traffic,” in Machine Learning and
Applications and Workshops (ICMLA), 2011 10th International Conference on,
vol. 2.1em plus 0.5em minus 0.4emIEEE, 2011, pp. 358–361.

[11]
M. Korczyński and A. Duda, “Markov chain fingerprinting to classify
encrypted traffic,” in Infocom, 2014 Proceedings IEEE.1em plus
0.5em minus 0.4emIEEE, 2014, pp. 781–789.

[12]
Y. Kumano, S. Ata, N. Nakamura, Y. Nakahira, and I. Oka, “Towards real-time
processing for application identification of encrypted traffic,” in
Computing, Networking and Communications (ICNC), 2014 International
Conference on.1em plus 0.5em minus 0.4emIEEE, 2014, pp.

[13]
“Cisco NetFlow v1, v5, v7 and v8,”
accessed: 2017-01-31.

[14]
“Cisco NetFlow v9,”
accessed: 2017-01-31.

[15]
“Audit Record Generation and Utilization System,”
<https://www.qosient.com/argus/index.shtml>, accessed: 2017-01-31.

[16]
“Zeek IDS,” <https://www.zeek.org/>, accessed: 2017-01-31.

[17]
“Verizon Injecting Perma-Cookies to Track Mobile Customers, Bypassing Privacy
Controls,” <https://www.eff.org/deeplinks/2014/11/verizon-x-uidh>,
accessed: 2021-09-29.

[18]
F. Pedregosa, G. Varoquaux, A. Gramfort, V. Michel, B. Thirion, O. Grisel,
M. Blondel, P. Prettenhofer, R. Weiss, V. Dubourg, J. Vanderplas, A. Passos,
D. Cournapeau, M. Brucher, M. Perrot, and E. Duchesnay, “Scikit-learn:
Machine learning in Python,” Journal of Machine Learning Research,
vol. 12, pp. 2825–2830, 2011.

[19]
T. Wang, X. Cai, R. Nithyanand, R. Johnson, and I. Goldberg, “Effective
attacks and provable defenses for website fingerprinting,” in 23rd
$\{$USENIX$\}$ Security Symposium ($\{$USENIX$\}$ Security 14), 2014, pp.

[20]
“UNIBS: Data sharing,” <http://netweb.ing.unibs.it/ ntw/tools/traces/>,
accessed: 2021-09-29.

[21]
T. Bujlow, V. Carela-Español, and P. Barlet-Ros, “Independent comparison
of popular dpi tools for traffic classification,” Computer Networks,
vol. 76, pp. 75–89, 2015.

[22]
“UPC Traffic Classification ,”
<https://cba.upc.edu/monitoring/traffic-classification>, accessed:

[23]
“Scikit-learn - AdaBoost classifier,”
accessed: 2021-09-29.

[24]
“Scikit-learn - Gradient boosting classifier,”
accessed: 2021-09-29.

[25]
“Scikit-learn - Gausian naive Bayes,”
accessed: 2021-09-29.

[26]
“Scikit-learn - K nearest neighbors classifier,”
accessed: 2021-09-29.

[27]
“Scikit-learn - Linear SVC,”
accessed: 2021-09-29.

[28]
“Scikit-learn - Logistic regression classifier,”
accessed: 2021-09-29.

[29]
“Scikit-learn - SGD classifier,”
accessed: 2021-09-29.

[30]
“Scikit-learn - Random forest classifier,”
accessed: 2021-09-29.

[31]
“Scikit-learn - Decision tree classifier,”
accessed: 2021-09-29.

[32]
S. Raschka, “Model evaluation, model selection, and algorithm selection in
machine learning,” arXiv preprint arXiv:1811.12808, 2018.

[33]
C. Novo and R. Morla, “Flow-based detection and proxy-based evasion of
encrypted malware c2 traffic,” in Proceedings of the 13th ACM Workshop
on Artificial Intelligence and Security, 2020, pp. 83–91.

[34]
M. Pietrzyk, J.-L. Costeux, G. Urvoy-Keller, and T. En-Najjary, “Challenging
statistical classification for operational usage: the adsl case,” in
Proceedings of the 9th ACM SIGCOMM conference on Internet
measurement.1em plus 0.5em minus 0.4emACM, 2009, pp. 122–135.

[35]
A. Custura, G. Fairhurst, and I. Learmonth, “Exploring usable path mtu in the
internet,” in Network Traffic Measurement and Analysis
Conference.1em plus 0.5em minus 0.4emIFIP Open Digital
Library, 2018.

[36]
“MAPIE - Model Agnostic Prediction Interval Estimator,”
<https://mapie.readthedocs.io/en/latest/index.html>, accessed:


# Effectively Counting $s-t$ Simple Paths in Directed Graphs

Mostafa Haghir Chehreghani
Department of Computer Engineering
Amirkabir University of Technology (Tehran Polytechnic)  Iran
<EMAIL_ADDRESS>
###### Abstract

An important tool in analyzing complex social and information networks is
$s-t$ simple path counting, which is known to be $\\#$P-complete. In this
paper, we study efficient $s-t$ simple path counting in directed graphs. For a
given pair of vertices $s$ and $t$ in a directed graph, first we propose a
pruning technique that can efficiently and considerably reduce the search
space. Then, we discuss how this technique can be adjusted with exact and
approximate algorithms, to improve their efficiency. In the end, by performing
extensive experiments over several networks from different domains, we show
high empirical efficiency of our proposed technique. Our algorithm is not a
competitor of existing methods, rather, it is a friend that can be used as a
fast pre-processing step, before applying any existing algorithm.

###### keywords:

Complex network analysis, directed graphs, simple paths, $s-t$ (simple) path
counting, exact algorithm, approximate algorithm

Effectively Counting $s-t$ Simple Paths in Directed Graphs

##  1 Introduction

Graphs (networks) are powerful tools that are used to model data in different
domains, including social networks, information networks, road networks and
the world wide web. A property seen in most of these real-world networks is
that the links between the vertices do not always represent reciprocal
relations [18]. Consequently, the generated networks are directed graphs,
where any edge has a direction and the edges are not necessarily symmetric.

Path counting is an important problem in analyzing large networks. For
example, some of network indices, such as $k$-path centrality are based on
counting paths [7]. A computationally demanding path counting problem is $s-t$
simple path counting. In this problem, given two vertices $s$ and $t$, the
goal is to count the number of all simple paths from $s$ to $t$. A simple path
is path that contains each node of the graph at most once. For both directed
and undirected graphs, this problem is known to be $\\#$P-complete [25]111It
is important to not confuse our studied problem with the problem wherein the
constraint ”simple” is removed, i.e., we count the number of all paths (simple
or non-simple) between two vertices $s$ and $t$. In this problem, a node may
appear for several times in a path. Unlike our studied problem (which is
$\\#$P-complete), this problem can be efficiently solved in polynomial time,
using e.g., dynamic programming or matrix multiplication. Another relevant
problem, that unlike our studied problem can be solved in polynomial time, is
counting the number of shortest paths between two vertices. . Hence, it is
critical to develop algorithms that work efficiently in practice.

In this paper, we focus on efficient $s-t$ simple path counting in directed
graphs. Our key contributions are as follows.

  * •

First, for a given pair of vertices $s$ and $t$, we propose a pruning
technique that can reduce the search space significantly. Moreover, it can be
computed very efficiently.

  * •

Then, we discuss how this technique can be adjusted with exact and approximate
algorithms, to improve their performance. We show, for example, that under
some conditions, applying this pruning technique yields a polynomial time
exact algorithm for enumerating $s-t$ simple paths.

  * •

Finally, by performing extensive experiments over several networks from
different domains (social, peer-to-peer, communication, citation, stack
exchange, product co-purchasing, …), we show high empirical efficiency of our
pruning technique.

Our algorithm is not a competitor of existing $s-t$ simple path counting
algorithms. Rather, it is a friend as it can be used as a fast pre-processing
step, before applying any of existing algorithms.

The rest of this paper is organized as follows. In Section 2, we introduce
preliminaries and necessary definitions used in the paper. In Section 3, we
give an overview on related work. In Section 4, we present our pruning
technique and discuss how it improves exact and approximate algorithms. In
Section 5, we empirically evaluate our pruning technique. Finally, the paper
is concluded in Section 6.

##  2 Preliminaries

We assume that the reader is familiar with basic concepts in graph theory.
Throughout the paper, $G$ refers to a directed graph. For simplicity, we
assume that $G$ is a connected and loop-free graph without multi-edges. By
default, we assume that $G$ is an unweighted graph, unless it is explicitly
mentioned that $G$ is weighted. $V(G)$ and $E(G)$ refer to the set of vertices
and the set of edges of $G$, respectively. For a vertex $v\in V(G)$, the
number of head ends adjacent to $v$ is called its in degree and the number of
tail ends adjacent to $v$ is called its out degree. For a vertex $v$, by
$\mathcal{ON}_{G}(v)$ we denote the set of outgoing neighbors of $v$ in the
graph $G$. A (directed) walk in a directed graph is a sequence of edges
directed in the same direction and joins a sequence of vertices. A directed
trail is a (directed) walk wherein all the edges are distinct. A (simple and
directed) path is a (directed) trail wherein all the vertices are distinct.
Let $A\subset V(G)$. The subgraph of $G$ induced by $A$ is the graph whose
vertices are $A$ and edges are all the edges in $E(G)$ that have both
endpoints in $A$.

##  3 Related work

There exist a number of algorithms in the literature for
counting/enumerating/listing $s-t$ simple paths. In one of the earlier
studies, Valiant [25] showed that it is $\\#$P-complete to enumerate simple
$s-t$ paths in both directed and undirected graphs. In [10], Hura utilized
Petri nets to enumerate all $s-t$ simple paths in a directed graph. Bax [2]
exploited Hamiltonian path algorithms to derive an
$O\left(2^{|V(G)|}poly(|V(G)|)\right)$ time and $O\left(poly(|V(G)|)\right)$
space algorithm for $s-t$ simple path counting. Knuth [11] presented the
simpath algorithm that generates a (not-always-reduced) binary decision
diagram for all simple paths from $s$ to $t$. The constructed diagram is
called zero-suppressed decision diagram (ZDD). In recent years, more efficient
algorithms have been proposed to construct ZDD [26]. Roberts and Kroese [22]
presented an approximate algorithm for $s-t$ simple path counting, which is
based on sequential importance sampling (we discuss their algorithm in more
details in Section 4.3, and explain how our pruning technique can improve it).
In [17], Mihalák et.al. studied approximately counting approximately shortest
paths in directed graphs. However, their algorithm is restricted to directed
acyclic graphs (DAGs).

In the literature, there are problems that are close to our studied problem,
or are a restricted form of it. In the following, we briefly review them:

  * •

In a restricted form of our studied problem, the number of $s-t$ simple paths
of a (fixed) length $l$ are counted. Flum and Grohe [8] proved that counting
simple paths of length $l$ on both directed and undirected graphs,
parameterized by $l$, is #W[1]-complete. Recently, Giscard et.al. [9] proposed
an algorithm for counting $s-t$ simple paths of length up to $l$. Time
complexity of their algorithm is
$O\left(|V(G)|+|E(G)|+(l^{w}+l\Delta)|S_{l}|\right)$, where $\Delta$ is the
maximum degree of the graph, $|S_{l}|$ is the number of weakly connected
induced subgraphs of $G$ on at most $l$ vertices, and $w$ is the exponent of
matrix multiplication.

  * •

The other problem which is close to our studied problem is counting the number
of all paths from a vertex $s$ to another vertex $t$. However, unlike our
studied problem (which is $\\#$P-complete), this problem can be solved in
polynomial time, using e.g., dynamic programming or fast matrix
multiplication. For example and as described in e.g., [3], the number of all
paths of size $k$ (simple or non-simple) from $s$ to $t$ is equal to the
$st$-th entry of $A^{k}$, where $A$ is the adjacency matrix of the graph and
$A^{k}$ is the $k$-th power of $A$. Arenas et.al. [1] study the problem of
counting the number of all paths between two vertices $s$ and $t$ of a length
at most $k$, and show that it admits a fully polynomial-time randomized
approximation scheme (FPRAS).

  * •

Another relevant problem is counting the number of shortest paths between a
given pair of vertices. Unlike our studied problem, this problem can be solved
efficiently in polynomial time, too. Using breadth first search (BFS), this
problem can be solved in $O\left(|V(G)|+|E(G)|\right)$ time [3]. Recently, a
number of faster approximate algorithms have been proposed to address this
problem in static [16] and dynamic graphs [24]. The notion of shortest paths
is used to develop several tools of network analysis, including betweenness
centrality [4, 6, 5].

In the current paper, we propose a pruning technique to improve
exact/approximate algorithms of the general form of $s-t$ simple path counting
in directed graphs. Our algorithm is not a competitor of these existing
algorithms. Rather, it is a friend that can be used as a fast pre-processing
step before applying any of them.

##  4 Counting $s-t$ paths in directed graphs

In this section, given a source vertex $s$ and a target vertex $t$ in a
directed graph $G$, we propose algorithms to count the number of all simple
paths from $s$ to $t$. First in Section 4.1, we present a pruning technique
used to reduce search space for $s-t$ path counting. Then in Section 4.2, we
discuss how this pruning technique improves exact $s-t$ path counting. Finally
in Section 4.3, we study how our pruning technique can improve approximate
algorithms.

###  4.1 A pruning technique

In this section, we describe our pruning technique. We start with introducing
the sets $\mathcal{S}_{st}$ (Definition 4.1) and $\mathcal{IS}_{ts}$
(Definition 4.3), that are used to reduce the search space.

###### Definition 4.1

Let $G$ be a directed graph and $s,t,v\in V(G)$ such that $s\neq t$. We say
$v$ is in the scope of $s$ with respect to $t$ iff either $v=s$, or there is
at least one directed path from $s$ to $v$ in $G$ that does not pass over $t$.
222The path may end to $t$, however, it can not pass over $t$. The set of
vertices $v$ that are in the scope of $s$ with respect to $t$ is denoted with
$\mathcal{S}_{st}$.

###### Definition 4.2

Let $G$ be a directed graph. Inverse graph of $G$, denoted with
$\mathcal{I}(G)$, is a directed graph such that: (i) $V(\mathcal{I}(G))=V(G)$,
and (ii) $(u,v)\in E(\mathcal{I}(G))$ if and only if $(v,u)\in E(G)$ [6].

###### Definition 4.3

Let $G$ be a directed graph and $s,t,v\in V(G)$ such that $s\neq t$. We say
$v$ is inversely in the scope of $t$ with respect to $s$ iff either $v=t$, or
there is at least one directed path from $t$ to $v$ in $\mathcal{I}(G)$ that
does not pass over $s$ 333The path may end to $s$, however, it can not pass
over $s$.. The set of vertices $v$ that are inversely in the scope of $t$ with
respect to $s$ is denoted with $\mathcal{IS}_{ts}$.

For example, in the graph of Figure 1, $\mathcal{S}_{st}$ consists of vertices
$s$, $v_{1}$, $v_{2}$, $v_{3}$, $v_{4}$, $v_{5}$ and $t$. It does not contain
vertex $v_{6}$, as the path from $s$ to $v_{6}$ passes over $t$. Moreover,
$\mathcal{IS}_{ts}$ consists of vertices $t$, $v_{5}$ and $s$. For a given
pair of vertices $s,t\in V(G)$, we can compute $\mathcal{S}_{st}$ and
$\mathcal{IS}_{ts}$ efficiently in $O(|E(G)|)$ time. To compute
$\mathcal{S}_{st}$, we act as follows:

  1. 1.

First, if $G$ is weighted, weights of the edges of $G$ are discarded,

  2. 2.

Then, a (revised) breadth first search (BFS) or a depth-first search (DFS) on
$G$ starting from $s$ is conducted, with an small change: when $t$ is met, the
traversal is not expanded to the children (and descendants) of $t$. All the
vertices that are met during the traversal are appended to $\mathcal{S}_{st}$.

To compute $\mathcal{IS}_{ts}$, we act as follows:

  1. 1.

First, by flipping the direction of the edges of $G$, we construct
$\mathcal{I}(G)$,

  2. 2.

Then, if $G$ is weighted, weights of the edges are discarded,

  3. 3.

Finally, a (revised) breadth-first search or depth-first search on
$\mathcal{I}(G)$ starting from $t$ is conducted, with an small change: when
$s$ is met, the traversal is not expanded to its children (and descendants).
All the vertices that are met during the traversal are appended to
$\mathcal{IS}_{ts}$.

It is easy to see that both $\mathcal{S}_{st}$ and $\mathcal{IS}_{ts}$ can be
computed in $O(|E(G)|)$ time, for both unweighted and weighted graphs.
Furthermore, we have the following lemma.

###### Lemma 4.4

Given a directed graph $G$ and vertices $s,t\in V(G)$, the number of simple
paths from $s$ to $t$ is equal to the number of simple paths from $s$ to $t$
whose vertices belong to $\mathcal{S}_{st}\cap\mathcal{IS}_{ts}$.

###### Proof 4.5

Each path should start from $s$ and end with $t$, therefore its vertices
should belong to both $\mathcal{S}_{st}$ and $\mathcal{IS}_{ts}$.

Lemma 4.4 says that in order to compute $s-t$ paths, we require to consider
only those paths that their vertices belong to both $\mathcal{S}_{st}$ and
$\mathcal{IS}_{ts}$. We can use this to prune many vertices of the graph that
do not belong to either $\mathcal{S}_{st}$ or $\mathcal{IS}_{ts}$ or both. As
a result, if in a graph, $\mathcal{S}_{st}\cap\mathcal{IS}_{ts}$ is
(considerably) smaller than the number of vertices of the graph, this pruning
technique can enormously improve the efficiency of the $s-t$ path counting
algorithm. Note that time complexity of $s-t$ simple path counting algorithms
is exponential in terms of the number of vertices of the graph. Hence,
discarding a considerable part of the vertices can hugely improve the
efficiency of the algorithms.

In Algorithm 1, we use this optimization technique to make $s-t$ path counting
algorithms more efficient. If $out\\_degree$ of $s$ is $0$ or $in\\_degree$ of
$t$ is $0$, there is no path from $s$ to $t$, hence, Algorithm 1 returns $0$.
Then, it computes $\mathcal{S}_{st}$ and $\mathcal{IS}_{ts}$ and stores them
respectively in the sets $S$ and $IS$. Then, for each vertex $v$ in the graph,
we check if it belongs to both $S$ and $IS$. We form the subgraph induced by
all such vertices and call it $\rho_{st}(G)$. Finally, we apply the $s-t$ path
counting algorithm on $\rho_{st}(G)$. the graph $\rho_{st}(G)$ is usually much
smaller than $G$, therefore, the algorithm can be run much faster.

Figure 1:  The graph on left is a directed graph $G$ and the graph on right
shows $\rho_{st}(G)$. Over $G$ and $s$ and $t$, while the variance of our
algorithm is $0$, the variance of the algorithm of [22] is $4$.

(a) A graph $G$

(b) $\rho_{st}(G)$

Algorithm 1 High level pseudo code of computing the number of $s-t$ paths in a
directed graph.

1: Input. A directed network $G$, a pair of vertices $s,t\in V(G)$.

2: Output. The number of simple paths from $s$ to $t$.

3: if out degree of $s$ is $0$ or in degree of $t$ is $0$ then

4: return $0$.

5: end if

6: $S\leftarrow$ compute $\mathcal{S}_{st}$.

7: $IS\leftarrow$ compute $\mathcal{IS}_{ts}$.

8: for all vertices $v\in V(G)$ do

9: if $v\in S$ and $v\in IS$ then

10: Mark $v$.

11: end if

12: end for

13: $\rho_{st}(G)\leftarrow$ the induced subgraph of $G$ consisting of the
marked vertices.

14: $path\\_no\leftarrow$ call an (exact/approximate) $s-t$ path counting
algorithm on $\rho_{st}(G)$.

15: return $path\\_no$.

###  4.2 Exact algorithm

In an exact algorithm for counting $s-t$ simple paths that exhaustively
enumerates all simple paths from $s$ to $t$, we can use backtracking: we start
from $s$, take a path and walk it (so that the path does not contain repeated
vertices). If the path ends to $t$, we count it and backtrack and take some
other path. If the path does not reach to $t$, we discard it and take some
other path. This algorithm enumerates all the paths from $s$ to $t$. We refer
to this algorithm as Exhaustive Simple Path Enumerator, or ESPE in short. As
discussed in [25], enumerating $s-t$ paths for both directed and undirected
graphs is $\\#$P-complete. In the following, we discuss that using our pruning
technique and for a certain type of vertices in a directed graph, the simple
path enumeration problem can be solved in polynomial time.

In the $s-t$ simple path counting problem, the input size consists of the
number of vertices and the number of edges of $G$. In Theorem 4.6 we show that
if for $s$ and $t$, $\left|\mathcal{S}_{st}\cap\mathcal{IS}_{ts}\right|$ is a
constant, then the $s-t$ simple counting problem can be solved efficiently.
Later in Section 5 we empirically show that in most of real-world networks,
vertices $s$ and $t$ have usually a small quantity for
$\left|\mathcal{S}_{st}\cap\mathcal{IS}_{ts}\right|$. It should be highlighted
that even for such pairs of vertices, applying the above mentioned exhaustive
simple path enumeration algorithm without using our pruning technique, does
not yield a polynomial time algorithm.

###### Theorem 4.6

Let $G$ be a directed graph, where for a given pair of vertices $s$ and $t$,
$|V\left(\rho_{st}(G)\right)|)$ is a constant. The list of simple paths from
$s$ to $t$ can be enumerated in polynomial time (in terms of $|V(G)|$ and
$|E(G)|$).

###### Proof 4.7

We can compute the sets $\mathcal{S}_{st}$, $\mathcal{IS}_{ts}$, and their
intersection in $O(|E(G)|)$ time. So, $\rho_{st}(G)$ can be computed in
$O(|E(G)|)$ time. After that, we work with $\rho_{st}(G)$ which has only a
constant number of vertices. Hence, all the paths can be enumerated in a
constant time.

In fact, the proof of Theorem 4.6 presents a much better result than
polynomial time: it says that under the mentioned condition and after a linear
time spent to construct $\rho_{st}(G)$, the whole list can be enumerated in a
constant time (in terms of $|V(G)|$ and $|E(G)|$).

###  4.3 Approximate algorithm

Roberts and Kroese [22] proposed randomized algorithms for estimating the
number of $s-t$ simple paths in a graph. In this section, we investigate how
the first randomized algorithm of [22] can be improved using our pruning
technique 444The second randomized algorithm of [22] and most of the other
randomized algorithms in the literature can be improved in a similar way..
First, we briefly describe how this algorithm changes if we apply our pruning
technique on it. Then, we analyze this revised algorithm and compare it with
the first algorithm of [22].

In the approximate algorithm, we sample $N$ independent (random) paths and
compute the probability of sampling each path. The following procedure is used
to sample each path:

  1. 1.

Start with vertex $s$. Initialize the following variables: $c\leftarrow s$
(current vertex), $h\leftarrow 1$ (probability), $k\leftarrow 1$ (counter),
and $p\leftarrow\\{s\\}$ (path).

  2. 2.

Mark $s$ as $visited$ (to be sure that $p$ will not visit $s$ again).

  3. 3.

Let $V^{\prime}=\\{v\in\mathcal{ON}_{\rho_{st}(G)}(c)|v\text{ is not
}visited\\}$ be the set of possible vertices for the next step of the path.

  4. 4.

Choose the next vertex $v\in V^{\prime}$ of the path uniformly at random and
set $x_{k+1}$ to $v$.

  5. 5.

Set $c\leftarrow v$, $h\leftarrow\frac{h}{|V^{\prime}|}$, $k\leftarrow k+1$,
and $p\leftarrow p\cup\\{v\\}$. Mark $v$ as $visited$.

  6. 6.

If $c=t$, then stop. Otherwise go to step $3$.

If we do not use our pruning technique, instead of working with $\rho(G)$, we
should work (e.g., in step 3) with $G$. This has consequences. For example,
while a path sampled from $G$ may never reach to $t$ (so after some
iterations, $V^{\prime}$ becomes empty), this never happens when we work with
$\rho(G)$. Hence, as stated in [22], when working with $G$ we need to check
the following condition in step 3: if $V^{\prime}=\emptyset$, we do not
generate a valid path, so stop.

Let $x_{1},x_{2},\ldots,x_{m}$ be the set of vertices of each path $p$, with
$x_{1}=s$. In the above mentioned procedure, path $p$ is chosen with
probability:

$\displaystyle h(p)=$ $\displaystyle h(x_{2}|x_{1})\cdot
h(x_{3}|x_{1},x_{2})\ldots h(x_{m}|x_{1},\ldots,x_{m-1})$ $\displaystyle=$
$\displaystyle\left|\mathcal{ON}_{\rho_{st}(G)}(x_{1})\setminus\\{x_{1}\\}\right|\cdot\left|\mathcal{ON}_{\rho_{st}(G)}(x_{2})\setminus\\{x_{1},x_{2}\\}\right|\ldots$
$\displaystyle\left|\mathcal{ON}_{\rho_{st}(G)}(x_{m-1})\setminus\\{x_{1},\ldots,x_{m-1}\\}\right|$
(1)

Let $\mathcal{X}$ be the set of simple paths from $s$ to $t$. Let also
$p^{1},\ldots,p^{N}$ be the set of sampled paths. For each $p^{k}$, we
estimate the number of $s-t$ simple paths as $X^{k}=\frac{1}{h(p^{k})}$. The
final estimation is the average of estimations of different samples:

$X=\frac{1}{N}\sum_{k=1}^{N}\frac{1}{h(p^{k})}.$ (2)

It is easy to see that $X$ is an unbiased estimator for the number of simple
paths from $s$ to $t$:

$\mathbb{E}[X^{k}]=\sum_{p\in\mathcal{X}}\frac{1}{h(p)}\cdot
h(p)=\sum_{p\in\mathcal{X}}1=|\mathcal{X}|.$

Since $X$ is the average of $N$ random variables whose expected values is
$|\mathcal{X}|$, the expected value of $X$ is $|\mathcal{X}|$, too.

For the variance of $X$, we have:

$\displaystyle\mathbb{V}ar\left[X^{k}\right]$
$\displaystyle=\mathbb{E}\left[(X^{k})^{2}\right]-\mathbb{E}\left[X^{k}\right]^{2}$
$\displaystyle=\sum_{p\in\mathcal{X}}\frac{1}{h(p)^{2}}\cdot
h(p)-\mathcal{X}^{2}=\sum_{p\in\mathcal{X}}\frac{1}{h(p)}-\mathcal{X}^{2}.$

Since $X$ is the average of $N$ independent random variable
$X^{1},\ldots,X^{N}$, we have:

$\mathbb{V}ar\left[X\right]=\frac{1}{N}\sum_{p\in\mathcal{X}}\frac{1}{h(p)}-\frac{\mathcal{X}^{2}}{N}.$

Now let compare this variance with the variance of the estimator of the
algorithm of [22]. Unlike our algorithm, in the algorithm of [22] it is
possible that a sampled path does not belong to $\mathcal{X}$. Let
$\mathbf{X}$ be the sample space of the algorithm of [22] that consists of
those paths in $G$ that start from $s$ and end to $t$ or to a vertex whose
$V^{\prime}$ is empty. In a way similar to our analysis, its variance is:

$\displaystyle\mathbb{V}ar\left[X\right]$
$\displaystyle=\frac{1}{N}\sum_{p\in\mathbf{X}}\frac{I(p\in\mathbf{X})}{h^{\prime}(p)}-\frac{\mathcal{X}^{2}}{N}$
$\displaystyle=\frac{1}{N}\sum_{p:p\in\mathcal{X}}\frac{1}{h^{\prime}(p)}-\frac{\mathcal{X}^{2}}{N},$
(3)

where $h^{\prime}(p)$ is defined as follows ($x_{1},\ldots,x_{m}$ are vertices
of $p$):

$\displaystyle h^{\prime}(p)=$
$\displaystyle\left|\mathcal{ON}_{G}(x_{1})\setminus\\{x_{1}\\}\right|\cdot\left|\mathcal{ON}_{G}(x_{2})\setminus\\{x_{1},x_{2}\\}\right|\ldots$
$\displaystyle\left|\mathcal{ON}_{G}(x_{m-1})\setminus\\{x_{1},\ldots,x_{m-1}\\}\right|.$
(4)

Moreover in Equation 4.3, $I(p\in\mathcal{X})$ is an indicator function, which
is $1$ if $p$ is a path from $s$ to $t$, and $0$ otherwise.

Comparing Equation 4.3 with Equation 4.3 reveals that our algorithm has a
better (lower) variable. In order to have a lower variance, an algorithm
should assign as low as possible probabilities to the paths that it may sample
but they do not belong to $\mathcal{X}$. Because in this case their inverse
will be a very large number and the contribution of these large numbers to the
variance will be $0$. Our algorithm assigns probability $0$ to such paths, so
it manages them very efficiently. However, the algorithm of [22] may assign a
(large) non-zero value to many of them. For example, in Figure 1 assume that
we want to estimate the number of simple paths from vertex $s$ to vertex $t$.
The sample space of our algorithm consists of only only one path which is the
path $s\rightarrow v_{5}\rightarrow t$. Hence, the variance of its estimation,
even when $N$ is $1$, is $0$. In contrast, the sample space of the algorithm
of [22] consists of $5$ paths: $s\rightarrow v_{1}$, $s\rightarrow v_{2}$,
$s\rightarrow v_{3}$, $s\rightarrow v_{4}$, and $s\rightarrow v_{5}\rightarrow
t$, where it assigns an equal probability $1/5$ to each path. Among them, only
the last one belongs to $\mathcal{X}$. While this algorithm assigns the same
probability to each one, our algorithm gives probability $0$ to the paths that
do not belong to $\mathcal{X}$. With $N=1$, the variance of the estimation of
the algorithm of [22] is $4$. Since our algorithm assigns probability $0$ to
the paths that are not in $\mathcal{X}$, it never finds a variance worse than
the algorithm of [22].

Table 1:  Specifications of our large networks.

Dataset | #vertices | #edges | domain
---|---|---|---
soc-sign-Slashdot090221 [13] | 82140 | 549202 | social
soc-sign-epinions [13] | 131828 | 841372 | social
p2p-Gnutella08 [15] | 6301 | 20777 | peer-to-peer file sharing
p2p-Gnutella30 [15] | 36682 | 88328 | peer-to-peer file sharing
p2p-Gnutella31 [15] | 62586 | 147892 | peer-to-peer file sharing
wikiVote [13] | 7115 | 103689 | Wikipedia vote
collegeMsg [19] | 1899 | 59835 | message exchanging
soc-Epinions1 [21] | 75879 | 508837 | who-trust-whom
email-EuAll [15] | 265214 | 420045 | communication
cit-HepPh [14] | 34546 | 421578 | citation
cit-HepTh [14] | 27770 | 352807 | citation
sx-mathoverflow [20] | 24818 | 506550 | stack exchange (Math Overflow)
sx-askubuntu [20] | 159316 | 964437 | stack exchange (Ask Ubuntu)
sx-superuser [20] | 194085 | 1443339 | stack exchange (Super User)
amazon0302 [12] | 262111 | 1234877 | product co-purchasing
amazon0312 [12] | 400727 | 3200440 | product co-purchasing

##  5 Experimental results

In order to evaluate the empirical behaviour of our pruning technique, we
perform extensive experiments over several real-world networks from different
domains. The objective is to empirically show that in real-world directed
networks, for a given pair of vertices $s$ and $t$,
$\mathcal{S}_{st}\cap\mathcal{IS}_{ts}$ is small (compared to $V(G)$), and it
can be computed very efficiently. These two automatically improve efficiency
of any (exact/approximate) simple path counting algorithm.

Table 2: Empirical evaluation of our proposed pruning technique. Column $max_{sz}$ presents the maximum size of $\mathcal{S}_{st}\cap\mathcal{IS}_{ts}$ divided by the number of vertices (over different pairs of $s,t$ vertices). Column $max_{tm}$ presents the maximum time to compute $\mathcal{S}_{st}\cap\mathcal{IS}_{ts}$, where the reported times are in second(s).  Dataset |  $max_{sz}$ |  $max_{tm}$
---|---|---
soc-sign-Slashdot090221 |  0.3333 |  0.0889
soc-sign-epinions |  0.3144 |  0.1172
p2p-Gnutella08 |  0.3296 |  0.0014
p2p-Gnutella31 |  0.2261 |  0.0150
p2p-Gnutella30 |  0.2315 |  0.0081
wikiVote |  0.1570 |  0.0068
collegeMsg |  0.6821 |  0.0004
email-EuAll |  0.1289 |  0.0618
soc-Epinions1 |  0.4246 |  0.0761
cit-HepPh |  0.0000002 |  0.2833
cit-HepTh |  0.0000002 |  0.2765
sx-mathoverflow |  0.1478 |  0.0495
sx-askubuntu |  0.1160 |  0.1399
sx-superuser |  0.1477 |  0.2418
amazon0302 |  0.9229 |  0.2255
amazon0312 |  0.9488 |  0.5525
Table 3:  Specifications of our smaller networks. Dataset | #vertices | #edges | domain
---|---|---|---
bio-CE-GN [23] | 2220 | 53683 | Gene functional associations
bio-SC-LC [23] | 2004 | 20452 | Gene functional associations
bn-mouse_retina_1 [23] | 1123 | 90811 | Brain network
C1000-9 [23] | 1001 | 450081 | DIMACS
C2000-9 [23] | 2001 | 1799534 | DIMACS
SmaGri [23] | 1060 | 4921 | Citation network
Table 4: Empirical comparison of our proposed pruning technique. All the reported times are in seconds and ”not terminated” means the algorithm does not terminate within 12 hours!  Dataset | max $\left|\mathcal{S}\cap\mathcal{IS}\right|$ | max total running time | max $\\#$simple paths
---|---|---|---
|  | pruned-ESPE | ESPE |
bio-CE-GN | 135 | 5.253 | not terminated | 1,180,840
bio-SC-LC | 144 | 414.605 | not terminated | 119,738,817
bn-mouse_retina_1 | 29 | 183.987 | not terminated | 17,012,976
C1000-9 | 35 | 195.333 | 448.225 | 919,343,373
C2000-9 | 34 | 171.182 | 253.073 | 887,306,870
SmaGri | 52 | 0.001 | 5.973 | 1,030

The experiments are done on an Intel processor clocked at $2.4$ GHz with $3$
GB main memory, running Ubuntu Linux $18.04.2$ LTS. We perform our tests over
16 real-world datasets from the SNAP
repository555https://snap.stanford.edu/data/. The specifications of the
datasets are summarized in Table 1. Over each dataset, we choose $500$ pairs
of vertices $s$ and $t$, uniformly at random. Then, we run our proposed
pruning technique for each pair. In the end, we report maximum size of
$\mathcal{S}_{st}\cap\mathcal{IS}_{ts}$ divided by $|V(G)|$ (in the $max_{sz}$
column) and the maximum time to compute
$\mathcal{S}_{st}\cap\mathcal{IS}_{ts}$ (in the $max_{tm}$ column), over all
the sampled pairs.

The results are reported in Table 2. On the one hand, since $max_{sz}$ is
usually considerably smaller than the number of vertices of the graph, our
pruning technique usually makes the real-world graphs considerably smaller.
Hence, any $s-t$ path counting problem (exact, approximate, …) can be solved
much faster on the smaller graph. Compared to the number of vertices, we may
consider $\mathcal{S}_{st}\cap\mathcal{IS}_{ts}$ as a constant. Therefore
after pruning, the $s-t$ path counting algorithm finds a constant time. Note
that the reported values for $|\mathcal{S}_{st}\cap\mathcal{IS}_{ts}|/|V(G)|$
are the maximum ones and most of the other values are much smaller.

On the other hand, the values reported in the $max_{tm}$ column indicate that
computing $\mathcal{I}_{st}$, $\mathcal{IS}_{ts}$, and their intersection can
be done very quickly. Note that $max_{tm}$ contains all the time required to
compute $\mathcal{S}_{st}$, $\mathcal{IS}_{ts}$, and their intersection. As
reported in the table, over all the datasets, $max_{tm}$ is less than $0.6$
seconds! This time is quite ignorable, compared to the time of counting all
$s-t$ simple paths. Note that while worst case time complexity of our pruning
technique is $O(|E(G)|)$, in practice it is usually performed much faster. The
reason is that for most of vertices, the direction of edges makes the sets
$\mathcal{S}$ and $\mathcal{IS}$ small. As a result, BFSs (or DFSs) conducted
on these small sets, can be conducted very efficiently.

In our experiments, we tried to compare the ESPE algorithm (discussed in
Section 4.2), against the method wherein first we apply our pruning technique
and then we use ESPE (we refer to this method as pruned-ESPE). However, over
any of the above mentioned networks, none of the methods ESPE and pruned-ESPE
finish during a reasonable time (24 hours!). Therefore and to compare the
efficiency of pruned-ESPE against ESPE, we switch to smaller graphs. Our
smaller real-world graphs have between 1K-2.3K vertices. Their specifications
are summarized in Table 3 (they can be downloaded from
http://networkrepository.com/networks.php). We treat all these networks as
directed graphs.

When comparing pruned-ESPE against ESPE, we perform the following scenario:
over each dataset, we sample uniformly at random 500 pairs of nodes $s$ and
$t$, such that there exists at least one path from $s$ to $t$. Then for each
pair, we run pruned-ESPE and ESPE and compute the following statistics:
$\left|\mathcal{S}\cap\mathcal{IS}\right|$, total running time of pruned-ESPE
(which includes the pre-processing/pruning phase), running time of ESPE, and
the number of simple paths. In the end and among all the sampled pairs, we
choose the pair that yields the maximum number of simple paths and report the
results of this chosen pair. For both algorithms, this chosen pair usually
takes the most amount of time to process. Thus, ”max
$\left|\mathcal{S}\cap\mathcal{IS}\right|$” depicts the value of
$\left|\mathcal{S}\cap\mathcal{IS}\right|$ for the chosen pair, ”max total
running time” depicts the total time to process the chosen pair (by each of
the algorithms), and ”max #simple paths” presents the number of simple paths
between the two nodes of the chosen pair.

Table 4 reports the results. As seen in the table, using our proposed pruning
techniques significantly improves the running time. All the times reported in
this table are in seconds, and ”not terminated” means the algorithm does not
terminate within 12 hours! So, for the chosen pair and over bio-CE-GN, bio-SC-
LC and bn-mouseretina1, while ESPE doe not produce the results during a
reasonable time (12 hours!), pruned-ESPE processes the chosen pairs quickly.
Over the other graphs that both algorithms terminate within a reasonable time,
pruned-ESPE acts quite faster. As evidenced by the values of the ”max
$\left|\mathcal{S}\cap\mathcal{IS}\right|$” column, this is due to the
efficiency of our pruning technique in considerably reducing the search space.

##  6 Conclusion

In this paper, we studied efficient $s-t$ simple path counting in directed
graphs. For a given pair of vertices $s$ and $t$ in a directed graph, first we
presented a pruning technique that can efficiently reduce the search space.
Then, we investigated how this pruning technique can be adjusted with exact
and approximate algorithms, to improve their efficiency. In the end, we
performed extensive experiments over several real-world networks to show the
empirical efficiency of our proposed pruning technique.

## References

  * [1] Marcelo Arenas, Luis Alberto Croquevielle, Rajesh Jayaram, and Cristian Riveros.  Efficient logspace classes for enumeration, counting, and uniform generation.  In Dan Suciu, Sebastian Skritek, and Christoph Koch, editors, Proceedings of the 38th ACM SIGMOD-SIGACT-SIGAI Symposium on Principles of Database Systems, PODS 2019, Amsterdam, The Netherlands, June 30 - July 5, 2019, pages 59–73. ACM, 2019.
  * [2] Eric T. Bax.  Algorithms to count paths and cycles.  Inf. Process. Lett., 52(5):249–252, 1994.
  * [3] U. Brandes.  A faster algorithm for betweenness centrality.  Journal of Mathematical Sociology, 25(2):163–177, 2001.
  * [4] Mostafa Haghir Chehreghani.  An efficient algorithm for approximate betweenness centrality computation.  Comput. J., 57(9):1371–1382, 2014.
  * [5] Mostafa Haghir Chehreghani, Talel Abdessalem, and Albert Bifet.  Metropolis-hastings algorithms for estimating betweenness centrality.  In Melanie Herschel, Helena Galhardas, Berthold Reinwald, Irini Fundulaki, Carsten Binnig, and Zoi Kaoudi, editors, Advances in Database Technology - 22nd International Conference on Extending Database Technology, EDBT 2019, Lisbon, Portugal, March 26-29, 2019, pages 686–689. OpenProceedings.org, 2019.
  * [6] Mostafa Haghir Chehreghani, Albert Bifet, and Talel Abdessalem.  Efficient exact and approximate algorithms for computing betweenness centrality in directed graphs.  In 22nd Pacific-Asia Conference on Knowledge Discovery and Data Mining, PAKDD, 2018.
  * [7] Mostafa Haghir Chehreghani, Albert Bifet, and Talel Abdessalem.  Adaptive algorithms for estimating betweenness and _k_ -path centralities.  In Wenwu Zhu, Dacheng Tao, Xueqi Cheng, Peng Cui, Elke A. Rundensteiner, David Carmel, Qi He, and Jeffrey Xu Yu, editors, Proceedings of the 28th ACM International Conference on Information and Knowledge Management, CIKM 2019, Beijing, China, November 3-7, 2019, pages 1231–1240. ACM, 2019.
  * [8] Jörg Flum and Martin Grohe.  The parameterized complexity of counting problems.  SIAM J. Comput., 33(4):892–922, 2004.
  * [9] Pierre-Louis Giscard, Nils M. Kriege, and Richard C. Wilson.  A general purpose algorithm for counting simple cycles and simple paths of any length.  Algorithmica, 81(7):2716–2737, 2019.
  * [10] G.S. Hura.  Enumeration of all simple paths in a directed graph using petri net: A systematic approach.  Microelectronics Reliability, 23(1):157 – 159, 1983.
  * [11] Donald Knuth.  The Art of Computer Programming, volume 4A.  Addison-Wesley Professional, Boston, MA, USA, 2011.
  * [12] Jure Leskovec, Lada A. Adamic, and Bernardo A. Huberman.  The dynamics of viral marketing.  ACM Transactions on the Web (TWEB), 1(1), 2007.
  * [13] Jure Leskovec, Daniel P. Huttenlocher, and Jon M. Kleinberg.  Signed networks in social media.  In Proceedings of the 28th International Conference on Human Factors in Computing Systems, Atlanta, Georgia, USA, April 10-15, pages 1361–1370, 2010.
  * [14] Jure Leskovec, Jon M. Kleinberg, and Christos Faloutsos.  Graphs over time: densification laws, shrinking diameters and possible explanations.  In Robert Grossman, Roberto J. Bayardo, and Kristin P. Bennett, editors, Proceedings of the Eleventh ACM SIGKDD International Conference on Knowledge Discovery and Data Mining, Chicago, Illinois, USA, August 21-24, 2005, pages 177–187. ACM, 2005.
  * [15] Jure Leskovec, Jon M. Kleinberg, and Christos Faloutsos.  Graph evolution: Densification and shrinking diameters.  ACM Transactions on Knowledge Discovery from Data (TKDD), 1(1), 2007.
  * [16] Dennis Nii Ayeh Mensah, Hui Gao, and Liang Wei Yang.  Approximation algorithm for shortest path in large social networks.  Algorithms, 13(2):36, 2020.
  * [17] Matús Mihalák, Rastislav Srámek, and Peter Widmayer.  Approximately counting approximately-shortest paths in directed acyclic graphs.  Theory Comput. Syst., 58(1):45–59, 2016.
  * [18] M. E. J. Newman.  The structure and function of complex networks.  SIAM REVIEW, 45:167–256, 2003.
  * [19] Pietro Panzarasa, Tore Opsahl, and Kathleen M. Carley.  Patterns and dynamics of users’ behavior and interaction: Network analysis of an online community.  JASIST, 60(5):911–932, 2009.
  * [20] Ashwin Paranjape, Austin R. Benson, and Jure Leskovec.  Motifs in temporal networks.  In Maarten de Rijke, Milad Shokouhi, Andrew Tomkins, and Min Zhang, editors, Proceedings of the Tenth ACM International Conference on Web Search and Data Mining, WSDM 2017, Cambridge, United Kingdom, February 6-10, 2017, pages 601–610. ACM, 2017.
  * [21] Matthew Richardson, Rakesh Agrawal, and Pedro M. Domingos.  Trust management for the semantic web.  In Dieter Fensel, Katia P. Sycara, and John Mylopoulos, editors, The Semantic Web - ISWC 2003, Second International Semantic Web Conference, Sanibel Island, FL, USA, October 20-23, 2003, Proceedings, volume 2870 of Lecture Notes in Computer Science, pages 351–368. Springer, 2003.
  * [22] Ben Roberts and Dirk P. Kroese.  Estimating the number of s-t paths in a graph.  J. Graph Algorithms Appl., 11(1):195–214, 2007.
  * [23] Ryan A. Rossi and Nesreen K. Ahmed.  The network data repository with interactive graph analytics and visualization.  In Blai Bonet and Sven Koenig, editors, Proceedings of the Twenty-Ninth AAAI Conference on Artificial Intelligence, January 25-30, 2015, Austin, Texas, USA, pages 4292–4293. AAAI Press, 2015.
  * [24] Konstantin Tretyakov, Abel Armas-Cervantes, Luciano García-Bañuelos, Jaak Vilo, and Marlon Dumas.  Fast fully dynamic landmark-based estimation of shortest path distances in very large graphs.  In Craig Macdonald, Iadh Ounis, and Ian Ruthven, editors, Proceedings of the 20th ACM Conference on Information and Knowledge Management, CIKM 2011, Glasgow, United Kingdom, October 24-28, 2011, pages 1785–1794. ACM, 2011.
  * [25] Leslie G. Valiant.  The complexity of enumeration and reliability problems.  SIAM J. Comput., 8(3):410–421, 1979.
  * [26] Norihito Yasuda, Teruji Sugaya, and Shin-ichi Minato.  Fast compilation of s-t paths on a graph for counting and enumeration.  In Antti Hyttinen, Joe Suzuki, and Brandon M. Malone, editors, Proceedings of the 3rd Workshop on Advanced Methodologies for Bayesian Networks, AMBN 2017, Kyoto, Japan, September 20-22, 2017, volume 73 of Proceedings of Machine Learning Research, pages 129–140. PMLR, 2017.