# Distributed Online Anomaly Detection for Virtualized Network Slicing
Environment

Weili Wang, Chengchao Liang, Qianbin Chen, , Lun Tang, Halim Yanikomeroglu W.
Wang, C. Liang, Q. Chen, and L. Tang are with the School of Communication and
Information Engineering and the Key Laboratory of Mobile Communication,
Chongqing University of Posts and Telecommunications, Chongqing 400065, China
(e-mail<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>tangl@cqupt.edu.cn).H. Yanikomeroglu is with the Department of Systems and
Computer Engineering, Carleton University, Ottawa, ON, Canada (e-mail:
halim@sce.carleton.ca).

###### Abstract

As the network slicing is one of the critical enablers in communication
networks, one anomalous physical node (PN) or physical link (PL) in substrate
networks that carries multiple virtual network elements can cause significant
performance degradation of multiple network slices. To recover the substrate
networks from anomaly within a short time, rapid and accurate identification
of whether or not the anomaly exists in PNs and PLs is vital. Online anomaly
detection methods that can analyze system data in real-time are preferred.
Besides, as virtual nodes and links mapped to PNs and PLs are scattered in
multiple slices, the distributed detection modes are required to adapt to the
virtualized environment. According to those requirements, in this paper, we
first propose a distributed online PN anomaly detection algorithm based on a
decentralized one-class support vector machine (OCSVM), which is realized
through analyzing real-time measurements of virtual nodes mapped to PNs in a
distributed manner. Specifically, to decouple the OCSVM objective function, we
transform the original problem to a group of decentralized quadratic
programming problems by introducing the consensus constraints. The alternating
direction method of multipliers is adopted to achieve the solution for the
distributed online PN anomaly detection. Next, by utilizing the correlation of
measurements between neighbor virtual nodes, another distributed online PL
anomaly detection algorithm based on the canonical correlation analysis is
proposed. The network only needs to store covariance matrices and mean vectors
of current data to calculate the canonical correlation vectors for real-time
PL anomaly analysis. The simulation results on both synthetic and real-world
network datasets show the effectiveness and robustness of the proposed
distributed online anomaly detection algorithms.

###### Index Terms:

Virtualized network slicing, online anomaly detection, distributed learning,
one-class support vector machine (OCSVM), canonical correlation analysis
(CCA).

##  I Introduction

Network slicing has been regarded as an important technology to satisfy
diverse requirements for emerging service types [1], [2]. Network function
virtualization (NFV) facilitates the deployment of network slices through
unleashing network functions from the dedicated hardware and implementing them
on general-purpose servers [3]. By NFV technique, network slicing can
customize virtual network functions based on demands to generate service
paths, i.e., customized service function chains (SFCs). Generally, multiple
network slices can be instantiated in a shared substrate network to support
diversified applications [4], [5].

A substrate network consists of multiple physical nodes (PNs) and physical
links (PLs), and a virtualized network slice consists of multiple virtual
nodes (VNs) and virtual links (VLs). Due to the many-to-one embedding
relationship between VNs (VLs) and PNs (PLs), network slices are prone to
unexpected and silent forms of failures caused by anomalies arising in the
shared substrate network [6]. As the performance of multiple network slices
hinges on the normal running of the shared substrate network, the accurate and
rapid anomaly detection for PNs and PLs in the substrate network is the
prerequisite for ensuring the service quality of network slices.

Anomaly detection refers to the problem of finding patterns in data that do
not conform to the expected behaviors [7]. Existing efforts [8, 9, 10, 11, 12,
13, 14, 15] on anomaly detection usually model all data of a time slot as a
long vector and consider it as a single sample in implementing data analysis.
Based on this idea, in virtualized network slicing environment, we can obtain
the traning dataset for PN anomaly detection through modeling all measurements
of VNs mapped to a PN within a period as a series of vectors. Then, normal
profiles can be trained to differentiate anomalies arising in PNs.
Furthermore, the VN measurements can also be utilized to implement PL anomaly
detection. Since the flow passes through each VN of an SFC in sequence, the
measurements between neighbor VNs are naturally related [16]. Therefore, we
can implement PL anomaly detection by analyzing the correlation of
measurements between neighbor VNs, which are mapped to both ends of the PL.
Through modeling all measurements of VNs mapped to both ends of a PL as two
sets of vectors, we can obtain the training dataset for the PL anomaly
detection. Note that the key difference between PL and PN anomaly detection is
that the anomaly detection for PLs is required to find anomalous correlation
patterns between two vectors.

However, since VNs are distributed in multiple slices, if we collect all
relevant VN data of a time slot to form vectors and analyze them in a
centralized mode, the data privacy of different slices will be compromised.
Besides, additional communication and storage cost can be introduced by
collecting all data to a central manager. Therefore, distributed detection
modes[17], [18], which can detect anomalies in PNs and PLs by analyzing
distributed VN datasets in respective managers, are more effective for
virtualized network slicing scenarios.

To recover the substrate networks from anomaly within a short time, the
identification of whether or not the anomaly exists in PNs and PLs should be
rapid and accurate. However, classical anomaly detection methods are usually
based on the batch formulation, which requires storing all historical data
within a period and trains the normal profiles through offline learning. The
process will not only introduce high storage and computation cost but prevent
timely detection of anomalies [19], [20]. Besides, in anomaly detection, the
training datasets are assumed to contain only one-class samples, i.e., normal
samples. If the number of anomalous data contained in training samples
increases, the performance of the detection model trained by the offline
learning can be greatly degraded [21]. Therefore, online anomaly detection
methods that can analyze system data in real-time are preferred.

Based on the above analysis, we propose a distributed online PN anomaly
detection algorithm based on a decentralized one-class support vector machine
(OCSVM) and a distributed online PL anomaly detection based on the canonical
correlation analysis (CCA). By the distributed online mode, the real-time
measurements of VNs can be analyzed in respective managers distributedly to
detect anomalies in PNs and PLs. To the best of our knowledge, this is the
first work which proposes the distributed online anomaly detection algorithms
for the virtualized network slicing environment. Our main contributions are
summarized as follows:

Firstly, based on a decentralized OCSVM, we propose a distributed online PN
anomaly detection algorithm for the virtualized network slicing environment.
Specifically, to decouple the OCSVM objective function, we transform the
original OCSVM objective function to decentralized quadratic programming
problems by introducing the consensus constraints. To realize a distributed
online implementation for PN anomaly detection, we establish a distributed
online augmented Lagrange function and solve it by the alternating direction
method of multipliers (ADMM).

Secondly, we propose a CCA-based distributed online PL anomaly detection
algorithm, which is realized through analyzing the real-time correlation of
measurements between neighbor VNs in respective managers distributedly. To
facilitate the real-time analysis, we derive a new method to calculate the
canonical correlation vectors, with which the managers only need to store
covariance matrices and mean vectors of current VN data instead of all
historical data.

Finally, the effectiveness and robustness of the proposed distributed online
anomaly detection algorithms are verified on both synthetic and real-world
network datasets.

The rest of the paper is organized as follows. Section II presents the related
works and Section III describes the system model. The distributed online PN
and PL anomaly detection algorithms are elucidated in Section IV and Section
V, respectively. Simulation results on both synthetic and real-world network
datasets are shown in Section VI. Section VII concludes the paper.

##  II Related Works

The key feature of the anomaly detection problem is that, a labeled dataset,
where each training sample is labeled as anomalous or normal, is usually
prohibitive to obtain [22]. Therefore, we usually deal with the anomaly
detection as a special classification problem by assuming that the entire
dataset contains only one-class samples, i.e., normal samples. OCSVM [23] is a
popular unsupervised method to detect anomalies in data. In literature [12,
13, 14, 24], it has been proved that the OCSVM algorithm has good
applicability to a wide range of anomaly detection problems.

CCA algorithm is an important method for analyzing the correlation between two
sets of data [25]. In [15, 26, 27, 28], the CCA algorithm has identified
whether or not the anomaly exits in the whole system successfully by detecting
the variation of correlation between two subsystems, which is also an
unsupervised learning method.

However, the classical OCSVM and CCA algorithms are both based on the batch
formulation, which requires storing all historical data within a period and
trains the normal profiles through offline learning. The process will not only
introduce high storage and computation cost but prevent timely detection of
anomalies. To conquer this difficulty, online anomaly detection algorithms
have attracted much attention in many domains since they can analyze the
system data in real-time and enable the timely detection of anomalies. For
instance, [20] proposed a bilateral principle component analysis algorithm to
realize the online and accurate traffic anomaly detection for Internet
management. In [29], an adaptive OCSVM method was designed to realize online
novelty detection for time series scenarios. [30] proposed an online subspace
tracking algorithm to supervise the flight data anomalies of the unmanned
aerial vehicle system. An online and unsupervised anomaly detection algorithm
for streaming data was designed in [31] using an array of sliding windows and
the probability density-based descriptors. These online algorithms showed good
anomaly detection performance with relatively short running time and low CPU
resources.

Centralized detection, where the training dataset is available in its entirety
to one centralized detector, is a well-studied area. However, if the dataset
is distributed over more than one location, different approaches need to be
taken. Therefore, some distributed algorithms have been designed to achieve
the global anomaly detection with datasets distributed in the networks. In
[17], a distributed version of principal component analysis algorithm was
developed to identify the global anomalies in distributed local datasets.
Based on the hierarchical temporal memory, a distributed anomaly detection
system was designed in [32] for the security of the in-vehicle network. [33]
proposed a collaborative intrusion detection framework for Vehicular Ad hoc
Networks, which enabled multiple controllers jointly train a global intrusion
detection model without direct sub-network flow exchange. In wireless sensor
networks, because data used for anomaly detection were distributedly
collected, [34] proposed a distributed online one-class support vector machine
algorithm for anomaly detection. This algorithm showed a good anomaly
detection performance with requiring relatively short running time and low CPU
resources. As the virtualized network slices utilize the network resources
spanning the whole substrate networks, the training dataset used for PN and PL
anomaly detection has a similar distributed property as wireless sensor
networks. Therefore, inspired by [34], we develop a distributed online anomaly
detection method that can be applied to the virtualized network slicing
environment.

##  III System Model

Figure 1: Network slice management and orchestration (MANO) overview.

As shown in Fig. 1, when the networks receive a network slice service, its
specific demand and required virtual network functions are analyzed by the
Service Manager to generate customized SFCs. The Network Slice Manager is
responsible for the configuration, life-cycle management and performance
monitoring of all VNs in SFCs. During the instantiation process of network
slices, multiple SFCs are embedded to a shared substrate network to provide
customized services for different slice requests [35]. The Infrastructure
Manager takes charge of the management and maintenance of PNs and PLs in
substrate networks.

Many unexpected behaviors can cause the anomaly of PNs and PLs in substrate
networks, including device malfunction, malicious attackers, fault
propagation, external environment changes, etc. Since multiple VNs and VLs
from different slices can be mapped to the same PN and PL, one anomalous PN or
PL in substrate networks can degrade the performance of multiple network
slices. Therefore, the accurate and rapid anomaly detection for PNs and PLs in
substrate networks is the prerequisite for ensuring the performance of network
slices.

We present a substrate network by an undirected graph $G=(Q,R)$, where $Q$ and
$R$ represent the PN and PL sets contained in the substrate network. Besides,
the set of VNs and VLs mapped to PN $q\;(q\in Q)$ and PL $r(r\in R)$ are
assumed to be $J_{q}$ and $J_{r}$. In the virtualized network slicing
environment, the Network Slice Manager will monitor the performance and
resource consumption and store the measurements of each customized VN
contained in an SFC. VN measurements include processing rate, data flow,
queuing delay, processing delay, and the usage of CPU, memory, etc. For each
VN $j_{q}\in J_{q}$ mapped to PN $q$, its measurements in a time slot can be
represented by
${\bm{x}}_{j_{q}}=\\{x_{j_{q}}^{1},x_{j_{q}}^{2},...,x_{j_{q}}^{p}\\}$, where
$p$ is the number of features. Therefore, the measurements related to the
anomaly detection of PN $q$ can be expressed as

${\bm{x}}_{q}=\left[{\begin{array}[]{*{20}c}{{\bm{x}}_{1}}\\\
{{\bm{x}}_{\text{2}}}\\\ \vdots\\\ {{\bm{x}}_{|J_{q}|}}\\\
\end{array}}\right]=\left({\begin{array}[]{*{20}c}{x_{1}^{1}}&{x_{1}^{2}}&\ldots&{x_{1}^{p}}\\\
{x_{\text{2}}^{1}}&{x_{\text{2}}^{2}}&\ldots&{x_{\text{2}}^{p}}\\\
\vdots&\vdots&\ddots&\vdots\\\
{x_{|J_{q}|}^{1}}&{x_{|J_{q}|}^{2}}&\ldots&{x_{|J_{q}|}^{p}}\\\
\end{array}}\right),$ (1)

where $|J_{q}|$ represents the number of VNs mapped to PN $q$.

The VN measurements can also be utilized to implement PL anomaly detection.
Since the flow passes through each VN of an SFC in sequence, the measurements
between neighbor VNs are naturally related. Therefore, we can implement PL
anomaly detection by analyzing the correlation of measurements between
neighbor VNs, which are mapped to both ends of the PL. The measurements
related to the anomaly detection of PL $r$ can be expressed as

$\begin{gathered}{\bm{x}}_{m}=\left[{\begin{array}[]{*{20}c}{{\bm{x}}_{1}}\\\
{{\bm{x}}_{\text{2}}}\\\ \vdots\\\ {{\bm{x}}_{|J_{r}|}}\\\
\end{array}}\right]=\left({\begin{array}[]{*{20}c}{x_{1}^{1}}&{x_{1}^{2}}&\ldots&{x_{1}^{p}}\\\
{x_{\text{2}}^{1}}&{x_{\text{2}}^{2}}&\ldots&{x_{\text{2}}^{p}}\\\
\vdots&\vdots&\ddots&\vdots\\\
{x_{|J_{r}|}^{1}}&{x_{|J_{r}|}^{2}}&\ldots&{x_{|J_{r}|}^{p}}\\\
\end{array}}\right),\hfill\\\
{\bm{x}}_{l}=\left[{\begin{array}[]{*{20}c}{{\bm{x}}_{1}}\\\
{{\bm{x}}_{\text{2}}}\\\ \vdots\\\ {{\bm{x}}_{|J_{r}|}}\\\
\end{array}}\right]=\left({\begin{array}[]{*{20}c}{x_{1}^{1}}&{x_{1}^{2}}&\ldots&{x_{1}^{d}}\\\
{x_{\text{2}}^{1}}&{x_{\text{2}}^{2}}&\ldots&{x_{\text{2}}^{d}}\\\
\vdots&\vdots&\ddots&\vdots\\\
{x_{|J_{r}|}^{1}}&{x_{|J_{r}|}^{2}}&\ldots&{x_{|J_{r}|}^{d}}\\\
\end{array}}\right),\hfill\\\ \end{gathered}$ (2)

where $m,l\in Q$ represent the both ends of PL $r$, $|J_{r}|$ is the number of
VLs mapped to PL $r$, and $p$ and $d$ denote the number of features for VNs
mapped to PN $m$ and $l$, respectively.

To realize the real-time detection for substrate networks in virtualized
network slicing scenarios, the Infrastructure Manager is required to collect
VN measurements from the Network Slice Manager sequentially if the centralized
detection modes are used. As a central manager, the management module of PN
$q$ can utilize measurements of all VNs mapped to it to realize anomalous PN
detection. Existing anomaly detection algorithms [8, 9, 10, 11, 12, 13, 14,
15] usually model all data of a time slot as a long vector and consider it as
a single sample in implementing data analysis, so the collected data in a time
slot used for the PN anomaly detection can be represented by

$\begin{gathered}{\bm{x}}_{q}=\\{{\bm{x}}_{1},...,{\bm{x}}_{|J_{q}|}\\}\hfill\\\
\;\;\;\;=\\{\underbrace{x_{1}^{1},x_{1}^{2},...,x_{1}^{p}}_{{\text{VN}}1\;{\text{mapped}}\;{\text{to}}\;{\text{PN}}\;q},\;...,\underbrace{x_{|J_{q}|}^{1},x_{|J_{q}|}^{2},...,x_{|J_{q}|}^{p}}_{{\text{VN
}}|J_{q}|\;{\text{mapped}}\;{\text{to}}\;{\text{PN}}\;q}\\}.\hfill\\\
\end{gathered}$ (3)

Through modeling all measurements of VNs in a time slot mapped to both ends of
a PL as two sets of vectors, we can obtain the data used for the PL anomaly
detection. The set of VLs mapped to PL $r$ is assumed to be $J_{r}$, so the
anomaly detection for PL $r$ can be realized through analyzing the correlation
between the following vectors:

$\begin{gathered}{\bm{x}}_{m}=\\{{\bm{x}}_{1},...,{\bm{x}}_{|J_{r}|}\\}\hfill\\\
\;\;\;\;=\\{\underbrace{x_{1}^{1},x_{1}^{2},...,x_{1}^{p}}_{{\text{VN}}1\;{\text{mapped}}\;{\text{to}}\;{\text{PN}}\;m},...,\underbrace{x_{{}_{|J_{r}|}}^{1},x_{{}_{|J_{r}|}}^{2},...,x_{{}_{|J_{r}|}}^{p}}_{{\text{VN
}}_{|J_{r}|}\;{\text{mapped}}\;{\text{to}}\;{\text{PN}}\;m}\\}\hfill\\\
{\bm{x}}_{l}=\\{{\bm{x}}_{1},...,{\bm{x}}_{|J_{r}|}\\}\hfill\\\
\;\;\;\;=\\{\underbrace{x_{1}^{1},x_{1}^{2},...,x_{1}^{d}}_{{\text{VN}}1\;{\text{mapped}}\;{\text{to}}\;{\text{PN}}\;l},...,\underbrace{x_{{}_{|J_{r}|}}^{1},x_{{}_{|J_{r}|}}^{2},...,x_{{}_{|J_{r}|}}^{d}}_{{\text{VN
}}_{|J_{r}|}\;{\text{mapped}}\;{\text{to}}\;{\text{PN}}\;l}\\}.\hfill\\\
\end{gathered}$ (4)

However, the centralized detection modes will introduce additional
communication and storage cost by collecting all VN measurements to the
Infrastructure Manager. Besides, the data privacy of different slices can be
compromised. Therefore, to determine whether or not the anomaly exists in PNs
and PLs rapidly and accurately with low communication and storage cost, the
distributed online PN and PL anomaly detection methods are required to analyze
the real-time measurements of VNs in respective Network Slice Manager
distributedly.

Notation: In this paper, boldface letters are used to denote vectors and
$(\cdot)^{\rm{T}}$ is the operator of transposition. $||\cdot||$ represents
the Euclidean norm. $\bm{0}$ and $\bm{1}$ denote the vectors of all zeros and
ones, respectively. Other notations will be defined if necessary.

##  IV Distributed Online PN Anomaly Detection Based on Decentralized OCSVM

From the perspective of the centralized mode, to implement the anomaly
detection for a PN, its management module will model measurements of all VNs
mapped to it within a period as a series of vectors to form the training set.
For each VN $j_{q}\in J_{q}$, its training set can be expressed as
$S_{j_{q}}=\\{({\bm{x}}_{j_{q}n},y_{j_{q}n}),n=1,2,....,N_{j_{q}}\\}$, where
${\bm{x}}_{j_{q}n}\in\mathbb{R}^{p}$, $N_{j_{q}}$ is the number of training
data, and $y_{j_{q}n}=+1$ indicates that all training data are normal samples.
As the training set of VN is $S_{j_{q}}$, the collected data of PN $q$ within
a period can be represented by $S_{q}{\rm{=}}\bigcup\nolimits_{j_{q}\in
J_{q}}{S_{j_{q}}}=\\{({\bm{x}}_{qn},y_{qn}),n=1,2,...,N_{q}\\}$, where
${\bm{x}}_{qn}=\bigcup\nolimits_{j_{q}\in J_{q}}{{\bm{x}}_{j_{q}n}}$ and
$y_{qn}=y_{1n}=,...,=y_{j_{q}n}=+1$. Generally, the training data are complex
and linearly inseparable. Therefore, OCSVM algorithm uses the feature mapping
function to map the training data from input space to high-dimensional feature
space. The purpose of the mapping is to make the data linearly separable in
the new space. Through finding a hyperplane in the feature space, OCSVM
algorithm can isolate the mapped samples from the origin with maximum margin
[23].

We use the feature mapping function $\bm{\phi}(\cdot)$ to map the input data
$\bm{x}_{q}$ into a reproducing kernel Hilbert space, and then calculate the
inner product of the feature mapping function by the kernel function
$k(\bm{x}_{q1},\bm{x}_{q2})$. The commonly used kernel function is the
Gaussian kernel function [23]

$\begin{gathered}k({\bm{x}}_{q1},{\bm{x}}_{q2})=\langle{\bm{\phi}}({\bm{x}}_{q1}),{\bm{\phi}}({\bm{x}}_{q2})\rangle\hfill\\\
\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;=\exp\left[{-||{\bm{x}}_{q1}-{\bm{x}}_{q2}||^{2}/\sigma^{2}}\right],\end{gathered}$
(5)

where $\sigma$ denotes the kernel width. According to formula (1), we can get
that
$k({\bm{x}}_{q},{\bm{x}}_{q})=<{\bm{\phi}}({\bm{x}}_{q}),{\bm{\phi}}({\bm{x}}_{q})>=1$,
which means that the mapping function $\bm{\phi}(\cdot)$ will map the original
samples into a hypersphere with center $\bm{0}$ and radius 1. As

$\begin{gathered}k({\bm{x}}_{q1},{\bm{x}}_{q2})=\exp\left[{-||{\bm{x}}_{q1}-{\bm{x}}_{q2}||^{2}/\sigma^{2}}\right]\hfill\\\
\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;=||{\bm{\phi}}({\bm{x}}_{q1})||\cdot||{\bm{\phi}}({\bm{x}}_{q2})||\cdot\cos\langle{\bm{\phi}}({\bm{x}}_{q1}),{\bm{\phi}}({\bm{x}}_{q2})\rangle,\end{gathered}$
(6)

if sample ${\bm{x}}_{q1}$ is normal while ${\bm{x}}_{q2}$ is distant from
${\bm{x}}_{q1}$ in input space, the value of $k({\bm{x}}_{q1},{\bm{x}}_{q2})$
will be small, which means that the value of
$\cos\langle{\bm{\phi}}({\bm{x}}_{q1}),{\bm{\phi}}({\bm{x}}_{q2})\rangle$ is
small and the angle between ${\bm{\phi}}({\bm{x}}_{q1})$ and
${\bm{\phi}}({\bm{x}}_{q2})$ is large. In this case, we usually consider
${\bm{x}}_{q2}$ as anomalous data, which can be isolated from normal data
through finding ($\bm{w}_{q}$, $\rho_{q}$), which are a weight vector and an
offset parameterizing a hyperplane
${\bm{w}}_{q}^{\rm{T}}{\bm{\phi}}({\bm{x}})-\rho_{q}=0$ in the feature space
associated with the kernel. The hyperplane can isolate the training samples
from the origin with maximum margin $(\rho_{q})/(||{\bm{w}}_{q}||)$.

Based on the above reasoning, the problem of the anomaly detection for PN $q$
can be formulated as the following optimization objective:

$\begin{gathered}\mathop{\min}\limits_{{\bm{w}}_{q},\rho_{q}}\;\frac{1}{2}||{\bm{w}}_{q}^{\rm{T}}||^{2}+C\sum\limits_{n=1}^{N_{q}}{\xi_{qn}}-\rho_{q}\hfill\\\
{\rm{s.t.}}\;{\bm{\phi}}^{\rm{T}}({\bm{x}}_{qn}){\bm{w}_{q}}\geq\rho_{q}-\xi_{qn},\;\;\;\;n=1,2,....,N_{q},\hfill\\\
\;\;\;\;\;\;\xi_{qn}\geq
0,\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;n=1,2,....,N_{q},\;\;\;\;\;\;\end{gathered}$
(7)

where $C$ is the penalty parameter, and $\xi_{qn}$ denotes the slack variable
that allows some samples to appear in the margin of the hyperplane. After
obtaining the optimal ${\bm{w}}_{q}^{*}$ and $\rho_{q}^{\rm{*}}$ for anomaly
detection, we can define the discriminant function as

$g({\bm{x}}_{qn})={\mathop{\rm
sgn}}({\bm{\phi}}^{\rm{T}}({\bm{x}}_{qn}){{\bm{w}_{q}}^{\rm{*}}}-\rho_{q}^{*}),\;\forall
q\in Q,\;n\in N_{q}.$ (8)

For any data ${\bm{x}}_{qn}$, it can be determined as a normal data if
$g({\bm{x}}_{qn})\geq 0$, and an anomalous data otherwise.

###  IV-A Decentralization of OCSVM

The above-mentioned centralized detection mode requires the Infrastructure
Manager collecting measurements of all VNs sequentially from the Network Slice
Manager and analyze them as a whole, which will not only bring high
communication and storage cost but compromise the data privacy of different
slices. Therefore, we distribute the OCSVM-based PN anomaly detection problem
to the Network Slice Manager instead of transmitting all VN measurements to
the Infrastructure Manager. The Network Slice Manager is responsible to
analyze the measurements of each VN.

To present the centralized problem (7) in a distributed form, a group of
decentralized quadratic programming problems that satisfy the consensus
constraints are established [36].

The global variables $\\{{\bm{w}}_{q},\rho_{q}\\}$ are replaced by auxiliary
local variables $\\{{\bm{w}}_{j_{q}},\rho_{j_{q}}\\}_{j_{q}=1}^{J_{q}}$ of
each VN. To ensure the consistency among local variables, extra consensus
constraints are added to the new objective. Therefore, the distributed form of
(7) can be expressed as

$\begin{gathered}\mathop{\min}\limits_{{\\{\bm{w}}_{j_{q}},\rho_{j_{q}}\\}}\;\frac{1}{2}\sum\limits_{j_{q}=1}^{|J_{q}|}{||{\bm{w}}_{j_{q}}||^{2}}+|J_{q}|C\sum\limits_{j_{q}=1}^{|J_{q}|}{\sum\limits_{n=1}^{N_{j_{q}}}{\xi_{j_{q}n}}}-\sum\limits_{j_{q}=1}^{|J_{q}|}{\rho_{j_{q}}}\hfill\\\
\;{\rm{s.t.}}\;\bm{\phi}^{\rm{T}}({\bm{x}}_{j_{q}n}){\bm{w}}_{j_{q}}\geq\rho_{j_{q}}-\xi_{j_{q}n},\;\;\forall
j_{q}\in J_{q},\;n=1,2,....,N_{j_{q}},\hfill\\\ \;\;\;\;\;\;\;\xi_{j_{q}n}\geq
0,\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\forall j_{q}\in
J_{q},\;n=1,2,....,N_{j_{q}},\hfill\\\
\;\;\;\;\;\;\;{\bm{w}}_{j_{q}}=\;{\bm{w}}_{i_{q}},\rho_{j_{q}}=\rho_{i_{q}},\;\;\;\;\;\;\;\;\forall
j_{q},i_{q}\in
J_{q}.\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\end{gathered}$ (9)

As the local variables $\\{{\bm{w}}_{j_{q}},\rho_{j_{q}}\\}_{j_{q}=1}^{J_{q}}$
of VNs mapped to the same PN satisfy
${\bm{w}}_{1}={\bm{w}}_{2}=,...,={\bm{w}}_{J_{q}}={\bm{w}}_{q}$ and
$\rho_{1}=\rho_{2}=,...,=\rho_{J_{q}}=\rho_{q}$, we can rewrite (9) as

$\begin{gathered}\mathop{\min}\limits_{{\\{\bm{w}}_{q},\rho_{q}\\}}\;|J_{q}|\;\left({\frac{1}{2}||{\bm{w}}_{q}||^{2}+C\sum\limits_{j_{q}=1}^{|J_{q}|}{{\bm{1}}_{j_{q}}^{\rm{T}}{\bm{\xi}}_{j_{q}}}-\rho_{q}}\right)\hfill\\\
\;{\rm{s.t.}}\;{\bm{\Phi}}^{\rm{T}}({\bm{X}}_{j_{q}}){\bm{w}}_{q}\geq\rho_{q}{\bm{1}}_{j_{q}}-{\bm{\xi}}_{j_{q}},\;\;\;\;\forall
j_{q}\in J_{q},\hfill\\\
\;\;\;\;\;\;\;\;{\bm{\xi}}_{j_{q}}\geq\bm{0}_{j_{q}},\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\forall
j_{q}\in J_{q},\;\;\end{gathered}$ (10)

where
$\bm{\Phi}{\rm{(}}{\bm{X}}_{j_{q}}{\rm{)=[}}\bm{\phi}({\bm{x}}_{j_{q}1}),...,\bm{\phi}({\bm{x}}_{j_{q}N_{j_{q}}}){\rm{]}}$
and $\\{{\bm{w}}_{q},\rho_{q}\\}$ is the feasible solution for (7). Since the
factor $|J_{q}|$ is a constant, the objective (10) is equivalent to (7).
Therefore, the distributed objective (9) is also equivalent to (7).

Similar to the objective (7),
$k({\bm{x}}_{j_{q}n},{\bm{x}}_{i_{q}m})=\langle\bm{\phi}({\bm{x}}_{j_{q}n}),\bm{\phi}({\bm{x}}_{i_{q}m})\rangle,\;\forall
j_{q},i_{q}\in J_{q},\;n\in N_{j_{q}},\;m\in N_{i_{q}},$ should be calculated
instead of $\bm{\phi}({\bm{x}}_{j_{q}n})$ or $\bm{\phi}({\bm{x}}_{i_{q}m})$ in
solving the distributed problem (9), where $\bm{\phi}({\bm{x}}_{j_{q}n})$ and
$\bm{\phi}({\bm{x}}_{i_{q}m})$ are unknown [34]. Since the measurements are
distributed in many VNs, it is difficult to calculate
$k({\bm{x}}_{j_{q}n},{\bm{x}}_{i_{q}m})$ without information exchange among
VNs. Therefore, the random approximation method [37] is adopted to approximate
$\bm{\phi}({\bm{x}}_{j_{q}n})$ with $\bm{z}({\bm{x}}_{j_{q}n})$.

Through introducing a random approximation function
${\bm{z}}:\mathbb{R}^{p}\to\mathbb{R}^{D}$, we can map the input data to a
random feature space, where $D$ is the dimension of the random feature space
and satisfies $D>p$. Using this technique, $\bm{\phi}({\bm{x}}_{j_{q}n})$ can
be approximated by $\bm{z}({\bm{x}}_{j_{q}n})$. The inner product calculation
between $\bm{\phi}({\bm{x}}_{j_{q}n})$ and $\bm{\phi}({\bm{x}}_{i_{q}m})$ can
be expressed as

$k({\bm{x}}_{j_{q}n},{\bm{x}}_{i_{q}m})=\langle\bm{\phi}({\bm{x}}_{j_{q}n}),{\bm{\phi}}({\bm{x}}_{i_{q}m})\rangle\approx{\bm{z}}({\bm{x}}_{j_{q}n})^{\rm{T}}{\bm{z}}({\bm{x}}_{i_{q}m}),$
(11)

where
${\bm{z}}({\bm{x}}_{j_{q}n})={\rm{[}}z_{{\bm{\omega}}_{1}}{\rm{(}}{\bm{x}}_{j_{q}n}{\rm{),}}...{\rm{,}}z_{{\bm{\omega}}_{D}}{\rm{(}}{\bm{x}}_{j_{q}n}{\rm{)]}}^{\rm{T}}$,
and ${\bm{z}}_{{\bm{\omega}}_{i}}({\bm{x}}_{j_{q}n})$ is a mapping function
and given by

$\begin{gathered}\bm{z}_{{\bm{\omega}}_{\bm{i}}}{\rm{(}}{\bm{x}}_{j_{q}n}{\rm{)}}=\sqrt{\frac{2}{D}}\cos({\bm{\omega}}_{i}^{\rm{T}}{\bm{x}}_{j_{q}n}+\vartheta_{i}),\;\;\hfill\\\
\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;i=1,2,...,D,\;\forall j_{q}\in
J_{q},\;n\in N_{j_{q}},\hfill\\\ \end{gathered}$ (12)

where $\vartheta_{i}$ is uniformly drawn from $[0,2\pi]$ and
${\bm{\omega}}_{i}$ is drawn from
$p({\bm{\omega}})=(2\pi)^{-(D/2)}e^{-[||{\bm{\omega}}||^{2}/2]}$, which
denotes the Fourier transform of Gaussian kernel function [37].

Note that the random approximate function $\bm{z}({\bm{x}}_{j_{q}n})$ and the
mapping feature function ${\bm{\phi}}({\bm{x}}_{j_{q}n})$ share the same
properties in reproducing kernel Hilbert space. $\bm{z}({\bm{x}}_{j_{q}n})$
also can map the original data into a hypersphere with center $\bm{0}$ and
radius 1 similar to ${\bm{\phi}}({\bm{x}}_{j_{q}n})$. Besides, as proved in
[37], if the dimension of $\bm{z}({\bm{x}}_{j_{q}n})$ is proper,
$k({\bm{x}}_{j_{q}n},{\bm{x}}_{i_{q}m})$ can be well approximated by
${\bm{z}}({\bm{x}}_{j_{q}n})^{\rm{T}}{\bm{z}}({\bm{x}}_{i_{q}m})$. Then, we
have

$\begin{gathered}k({\bm{x}}_{j_{q}n},{\bm{x}}_{i_{q}m})\hfill\\\
=||{\bm{\phi}}({\bm{x}}_{j_{q}n})||\cdot||{\bm{\phi}}({\bm{x}}_{i_{q}m})||\cdot\cos\langle{\bm{\phi}}({\bm{x}}_{j_{q}n}),{\bm{\phi}}({\bm{x}}_{i_{q}m})\rangle\hfill\\\
\approx||{\bm{z}}({\bm{x}}_{j_{q}n})||\cdot||{\bm{z}}({\bm{x}}_{i_{q}m})||\cdot\cos\langle{\bm{z}}({\bm{x}}_{j_{q}n}),{\bm{z}}({\bm{x}}_{i_{q}m})\rangle,\end{gathered}$
(13)

which means that
$\cos\langle{\bm{\phi}}({\bm{x}}_{j_{q}n}),{\bm{\phi}}({\bm{x}}_{i_{q}m})\rangle\approx\cos\langle{\bm{z}}({\bm{x}}_{j_{q}n}),{\bm{z}}({\bm{x}}_{i_{q}m})\rangle$.
Therefore, if data ${\bm{x}}_{j_{q}n}$ is close to ${\bm{x}}_{i_{q}m}$, their
corresponding mapping vectors ${\bm{z}}_{j_{q}n}$ and ${\bm{z}}_{i_{q}m}$ will
be close in the new random feature space.

Based on the above reasoning, the problem (9) can be transformed as follows:

$\begin{gathered}\mathop{\min}\limits_{{\\{\bm{w}}_{j_{q}},\rho_{j_{q}}\\}}\;\sum\limits_{j_{q}=1}^{J_{q}}{\left({\;\frac{1}{2}||{\bm{w}}_{j_{q}}||^{2}{\rm{+|}}J_{q}|C\sum\limits_{n=1}^{N_{j_{q}}}{\xi_{j_{q}n}}-\rho_{j_{q}}}\right)}\hfill\\\
\;{\rm{s.t.}}\;{\bm{z}}^{\rm{T}}({\bm{x}}_{j_{q}n}){\bm{w}}_{j_{q}}\geq\rho_{j_{q}}-\xi_{j_{q}n},\;\forall
j_{q}\in J_{q},\;n=1,2,....,N_{j_{q}},\hfill\\\ \;\;\;\;\;\;\xi_{j_{q}n}\geq
0,\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\forall j_{q}\in
J_{q},\;n=1,2,....,N_{j_{q}},\hfill\\\
\;\;\;\;\;\;{\bm{w}}_{j_{q}}=\;{\bm{w}}_{i_{q}},\rho_{j_{q}}=\rho_{i_{q}},\;\;\;\;\;\;\;\forall
j_{q},i_{q}\in J_{q},\hfill\\\ \end{gathered}$ (14)

where ${\bm{x}}_{j_{q}n}\in S_{j_{q}}$ denotes the training samples of VN
$j_{q}$. It should be noted that the symbol ${\bm{w}}_{j_{q}}$ in (14) has
different dimension with that in (9), but for simplicity, we choose to use the
same symbol to represent them.

###  IV-B Distributed Online PN Anomaly Detection

For each VN, new unlabeled measurements will be generated at each time. If we
keep all historical data of each VN within a period for offline training, it
will not only introduce high storage and computation cost but prevent timely
detection of anomalies. Therefore, unlabeled training data are expected to be
processed in an online mode. Based on the above analysis, we propose a
distributed online PN anomaly detection algorithm based on the decentralized
OCSVM.

Through introducing the Lagrange multipliers $\\{\kappa_{j}\\}$,
$\\{\lambda_{j}\\}$, $\\{\bm{\alpha}_{ji}\\}$ and $\\{\beta_{ji}\\}$, the
augmented Lagrange function for the problem (14) is expressed as

$\begin{gathered}L(\\{{\bm{w}}_{j_{q}}\\},\\{\rho_{j_{q}}\\},\\{\xi_{j_{q}n}\\},\\{\kappa_{j}\\},\\{\lambda_{j}\\},\\{\bm{\alpha}_{ji}\\},\\{\beta_{ji}\\})=\hfill\\\
\sum\limits_{j_{q}=1}^{|J_{q}|}{\left(\begin{gathered}\;\frac{1}{2}||{\bm{w}}_{j_{q}}||^{2}+|J_{q}|C\sum\limits_{n=1}^{N_{j_{q}}}{\xi_{j_{q}n}}-\rho_{j_{q}}-\frac{{\kappa_{j}}}{{N_{j_{q}}}}\sum\limits_{n=1}^{N_{j_{q}}}{\xi_{j_{q}n}}\hfill\\\
-\frac{{\lambda_{j}}}{{N_{j_{q}}}}\sum\limits_{n=1}^{N_{j_{q}}}{\left[{{\bm{z}}^{\rm{T}}({\bm{x}}_{j_{q}n}){\bm{w}}_{j_{q}}-\rho_{j_{q}}+\xi_{j_{q}n}}\right]}\hfill\\\
+\sum\limits_{i_{q}=1}^{|J_{q}|}{\bm{\alpha}_{ji}^{\rm{T}}({\bm{w}}_{j_{q}}-{\bm{w}}_{i_{q}})+}\sum\limits_{i_{q}=1}^{|J_{q}|}{\beta_{ji}(\rho_{j_{q}}-\rho_{i_{q}})}\hfill\\\
+\frac{\eta}{2}\sum\limits_{i_{q}=1}^{|J_{q}|}{[||{\bm{w}}_{j_{q}}-{\bm{w}}_{i_{q}}||^{2}+||\rho_{j_{q}}-\rho_{i_{q}}||^{2}]}\hfill\\\
\end{gathered}\right)}.\hfill\\\ \end{gathered}$ (15)

Here, the last term is the regularization term, which plays two roles: (1) It
eliminates the condition that $L(\bm{\theta}_{j})$ must be differentiable,
where $\bm{\theta}_{j}$ represents the set of variables; (2) The convergence
speed of ADMM could be controlled by adjusting the augmented Lagrange
parameter $\eta$. Then, $L(\bm{\theta}_{j})$ could be minimized in a cycle
fashion: at each iteration, we minimize $L(\bm{\theta}_{j})$ with respect to
one variable while keeping all other variables fixed [36]. The ADMM solution
at each iteration $t+1$ takes the form

$\begin{gathered}\\{{\bm{w}}_{j_{q}}(t+1),\rho_{j_{q}}(t+1),\xi_{j_{q}n}(t+1)\\}=\hfill\\\
\mathop{\arg\min}\limits_{\\{{\bm{w}}_{j_{q}}\\},\\{\rho_{j_{q}}\\},\\{\xi_{j_{q}n}\\}}L\left(\begin{gathered}\\{{\bm{w}}_{j_{q}}\\},\\{\rho_{j_{q}}\\},\\{\xi_{j_{q}n}\\},\\{\kappa_{j}\\}\hfill\\\
\\{\lambda_{j}\\},\bm{\alpha}_{ji}(t)\\},\\{\beta_{ji}(t)\\}\end{gathered}\right),\hfill\\\
\end{gathered}$ (16a)
$\bm{\alpha}_{j}(t+1)=\bm{\alpha}_{j}(t)+\frac{\eta}{2}\sum\limits_{i_{q}=1}^{|J_{q}|}{({\bm{w}}_{j_{q}}(t+1)-{\bm{w}}_{i_{q}}(t+1))},$
(16b)
$\beta_{j}(t+1)=\beta_{j}(t)+\frac{\eta}{2}\sum\limits_{i_{q}=1}^{|J_{q}|}{(\rho_{j_{q}}(t+1)-\rho_{i_{q}}(t+1))},$
(16c)

where $\bm{\alpha}_{j}(t)=\sum\nolimits_{i=1}^{|J_{q}|}{\bm{\alpha}_{ji}(t)}$
and $\beta_{j}(t)=\sum\nolimits_{i=1}^{|J_{q}|}{\beta_{ji}(t)}$.

Model parameters $\\{{\bm{w}}_{j_{q}}\\}$, $\\{\rho_{j_{q}}\\}$ and
$\\{\xi_{j_{q}n}\\}$ of decentralized OCSVM could be obtained by solving the
problem (16a). It is obvious that the problem (16a) is a batch formulation of
distributed algorithm, where all data of each VN are required to be available.
This will bring a big challenge to the storage resources of the networks. To
overcome this problem, a distributed online augmented Lagrange function is
defined in equation (17) by replacing $\rho_{i_{q}}$ and ${\bm{w}}_{i_{q}}$ at
time $t$ with $\frac{1}{2}(\rho_{j_{q}}(t)+\rho_{i_{q}}(t))$ and
$\frac{1}{2}({\bm{w}}_{j_{q}}(t)+{\bm{w}}_{i_{q}}(t))$ as in [36].

$\begin{gathered}L^{\prime}\left({\\{{\bm{w}}_{j_{q}}\\},\\{\rho_{j_{q}}\\},\xi_{j_{q}}(t),\\{\kappa_{j}\\},\\{\lambda_{j}\\},\\{{\bm{w}}_{j_{q}}(t)\\},\\{\rho_{j_{q}}(t)\\},\\{\bm{\alpha}_{j}(t)\\},\\{\beta_{j}(t)\\}}\right)=\hfill\\\
\sum\limits_{j_{q}=1}^{|J_{q}|}{\left(\begin{gathered}\;\frac{1}{2}||{\bm{w}}_{j_{q}}||^{2}+|J_{q}|C\xi_{j_{q}}(t)-\rho_{j_{q}}-\kappa_{j}\xi_{j_{q}}(t)-\lambda_{j}[{\bm{z}}^{\rm{T}}_{j_{q}}(t){\bm{w}}_{j_{q}}-\rho_{j_{q}}+\xi_{j_{q}}(t)]+2\bm{\alpha}_{j}^{\rm{T}}(t){\bm{w}}_{j_{q}}\hfill\\\
+2\beta_{j}(t)\rho_{j_{q}}+\frac{\eta}{2}\sum\limits_{i_{q}=1}^{|J_{q}|}{[||{\bm{w}}_{j_{q}}-\frac{1}{2}({\bm{w}}_{j_{q}}(t)+{\bm{w}}_{i_{q}}(t))||^{2}+||\rho_{j_{q}}-\frac{1}{2}(\rho_{j_{q}}(t)+\rho_{i_{q}}(t))||^{2}]}\hfill\\\
\end{gathered}\right)}.\hfill\\\ \end{gathered}$ (17)

In (17), $t$ is the time instant of online learning,
${\bm{z}}_{j_{q}}(t)={\rm{[}}z_{{\bm{\omega}}_{1}}({\bm{x}}_{j_{q}}(t)){\rm{,}}...{\rm{,}}z_{{\bm{\omega}}_{D}}({\bm{x}}_{j_{q}}(t)){\rm{]}}^{\rm{T}}$,
where ${\bm{x}}_{j_{q}}(t)$ is the training data of VN $j_{q}$ at time $t$.
$\xi_{j_{q}}(t)$ is the slack variable for ${\bm{x}}_{j_{q}}(t)$.

From the KKT conditions for (17) it follows that

${\bm{w}}_{j_{q}}(t+1)=\frac{1}{A}\left(\begin{gathered}\bm{z}_{j_{q}}(t+1)\lambda_{j}(t+1)-2\bm{\alpha}_{j}(t)\hfill\\\
+\frac{\eta}{2}\sum\limits_{i_{q}=1}^{|J_{q}|}{({\bm{w}}_{j_{q}}(t)+{\bm{w}}_{i_{q}}(t))}\hfill\\\
\end{gathered}\right),$ (18a)
$\rho_{j_{q}}(t+1)=\frac{1}{A-1}\left(\begin{gathered}1-\lambda_{j}(t+1)-2\beta_{j}(t)\hfill\\\
+\frac{\eta}{2}\sum\limits_{i_{q}=1}^{|J_{q}|}{(\rho_{j_{q}}(t)+\rho_{i_{q}}(t))}\hfill\\\
\end{gathered}\right),$ (18b) $0=|J_{q}|C-\lambda_{j}-\kappa_{j},$ (18c)

where $A=\eta|J_{q}|+1$. The KKT conditions require $\lambda_{j}\geqslant
0,\;\kappa_{j}\geqslant 0$, so (18c) is allowed to be replaced by
$0\leqslant\lambda_{j}\leqslant|J_{q}|C$. To implement the update (18a) and
(18b) for every VN, the optimal Lagrange multipliers $\\{\lambda_{j}(t+1)\\}$
are obtained by solving the dual problem of (17). The corresponding dual
function is given by

$\begin{gathered}L^{\prime}_{\lambda}\\{\lambda_{j}(t+1)\\}=\hfill\\\
-[\bm{z}_{j_{q}}^{\rm{T}}(t+1)A^{-1}\bm{z}_{j_{q}}(t+1)+(A-1)^{-1}]\lambda^{2}_{j}\hfill\\\
+[A^{-1}l_{j}(t)+(A-1)^{-1}(1-h_{j}(t))]\lambda_{j},\hfill\\\ \end{gathered}$
(19)

where
$l_{j}(t)=2\bm{\alpha}_{j}(t)-\frac{\eta}{2}\sum\limits_{i_{q}=1}^{|J_{q}|}{({\bm{w}}_{j_{q}}(t)+{\bm{w}}_{i_{q}}(t))}$
and
$h_{j}(t)=2\beta_{j}(t)-\frac{\eta}{2}\sum\limits_{i_{q}=1}^{|J_{q}|}{(\rho_{j_{q}}(t)+\rho_{i_{q}}(t))}$.
$\\{\lambda_{j}(t+1)\\}$ is given by

$\begin{gathered}\lambda_{j}{\text{(}}t+1{\text{)}}=\mathop{\arg\max}\limits_{0\leqslant\lambda_{j}\leqslant|J_{q}|C}-[\bm{z}_{j_{q}}^{\rm{T}}(t+1)A^{-1}\bm{z}_{j_{q}}(t+1)+\hfill\\\
(A-1)^{-1}]\lambda^{2}_{j}+[A^{-1}l_{j}(t)+(A-1)^{-1}(1-h_{j}(t))]\lambda_{j}.\hfill\\\
\end{gathered}$ (20)

Algorithm 1 presents the detailed steps of the distributed online PN anomaly
detection algorithm. In Algorithm 1, steps (2)-(8) are implemented in Network
Slice Manager and steps (9)-(13) are implemented in Infrastructure Manager.
Since the distributed online algorithm can eliminate the impact of anomalous
data on estimates in each iteration (step 12), it can hold high detection
accuracy without any labeled data.

Algorithm 1 OCSVM-based distributed online PN anomaly detection algorithm

0: Initialize the iteration number $T$, dimension $D$, multipliers
$\bm{\alpha}_{j}(0)$ and $\beta_{j}(0)$, and estimates ${\bm{w}}_{j_{q}}(0)$
and $\rho_{j_{q}}(0)$ of each VN

1: for $t=1,2,...,T$ do

2: for $j_{q}\in J_{q}$ do

3: VN $j_{q}$ obtains a new sample ${\bm{x}}_{j_{q}}(t)$, and calculate the
approximation ${\bm{z}}_{j_{q}}(t)$ of $\bm{\phi}(\bm{x}_{j_{q}}(t))$ by the
random approximate function

4: Compute $\lambda_{j}(t)$, ${\bm{w}}_{j_{q}}(t)$ and $\rho_{j_{q}}(t)$
according to (20), (18a) and (18b)

5: Compute $\bm{\alpha}_{j}(t)$ and $\beta_{j_{q}}(t)$ according to (16b) and
(16c)

6: Compute $g({\bm{x}}_{j_{q}}(t))={\mathop{\rm
sgn}}({\bm{z}}_{j_{q}}^{\rm{T}}(t){\bm{w}}_{j_{q}}(t)-\rho_{j_{q}}(t))$

7: Send $g({\bm{x}}_{j_{q}}(t))$, ${\bm{w}}_{j_{q}}(t)$ and $\rho_{j_{q}}(t)$
to the management module of PN $q$

8: end for

9: if $\prod\nolimits_{j_{q}\in J_{q}}{g({\bm{x}}_{j_{q}}(t))}==1$ then

10: PN $q$ is detected as normal at time $t$, and update estimates
${\bm{w}}_{j_{q}}(t)$ and $\rho_{j_{q}}(t)$. Then, broadcast them to related
VNs

11: else

12: PN $q$ is detected as anomaly at time $t$, and reserve estimates at time
$t-1$ and discard current ones

13: end if

14: end for

From formulas (18a) and (18b), we can know that VN $j_{q}$ requires estimates
${\bm{w}}_{i_{q}}(t)$ and $\rho_{i_{q}}(t)(i_{q}\in J_{q})$ to be available
when updating estimates ${\bm{w}}_{j_{q}}(t+1)$ and $\rho_{j_{q}}(t+1)$.
According to [23], the time complexity of classical OCSVM algorithm is
$O(N_{q}^{3})$. Since all training samples and their labels need to be
available in OCSVM, its storage complexity is $O(N_{q}|J_{q}|(p+1))$. The
calculation of the online OCSVM algorithm depends on the number of iterations
$T$ and the number of VNs $|J_{q}|$ mapped to PN $q$, so its time complexity
is $O(N_{q}|J_{q}|)$ when the number of iterations $T=N_{q}$. Since the online
algorithm only needs current measurements to be available, its storage
complexity is $O(|J_{q}|p)$. Due to the distributed online mode, the time and
storage complexity of the proposed PN anomaly detection algorithm are
$O(N_{q})$ and $O(p)$, respectively.

##  V CCA-based distributed online PL Anomaly Detection

The basic principle of PL anomaly detection: Since the flow passes through
each VN of an SFC in sequence, the measurements between neighbor VNs are
naturally related [16]. Assume that the virtual path
${\text{VN}}_{m}\mathop{\to}\limits^{{\text{VL}}_{m,l}}{\text{VN}}_{l}$ is
instantiated into the path
${\text{PN}}_{m}\mathop{\to}\limits^{{\text{PL}}_{m,l}}{\text{PN}}_{l}$, where
virtual link ${\text{VL}}_{m,l}$ is instantiated into physical link
${\text{PL}}_{m,l}$, and ${\rm{VN}}_{m}$ and ${\rm{VN}}_{l}$ are instantiated
into ${\rm{PN}}_{m}$ and ${\rm{PN}}_{l}$, respectively. Then, the correlation
of measurements between ${\rm{VN}}_{m}$ and ${\rm{VN}}_{l}$ will be stable
within a certain range if ${\rm{PL}}_{m,l}$ is in a normal state, and will
change if an anomaly occurs in ${\rm{PL}}_{m,l}$. As shown in Fig. 2, if the
working state of PL 34 becomes anomalous, the correlation of measurements
between VN 2 and 3 in Slice 2 will change. Therefore, we can implement PL
anomaly detection based on the correlation of measurements between neighbor
VNs, which are mapped to both ends of this PL.

Figure 2: The PL mapping schematic.

CCA is a widely used multivariate analysis algorithm [38]. Given two sets of
random variables $\bm{U}$ and $\bm{Y}$ where

$\begin{array}[]{l}\bm{U}=\left[{\begin{array}[]{*{20}c}{{\bm{u}}_{1}}\\\
{{\bm{u}}_{2}}\\\ \vdots\\\ {{\bm{u}}_{t}}\\\
\end{array}}\right]=\left[{\begin{array}[]{*{20}c}{u_{11}}&{u_{12}}&\ldots&{u_{1p}}\\\
{u_{21}}&{u_{22}}&\ldots&{u_{2p}}\\\ \vdots&\vdots&\ddots&\vdots\\\
{u_{t1}}&{u_{t2}}&\ldots&{u_{tp}}\\\ \end{array}}\right]\\\
\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;(t\;{\rm{samples}}\;p\;{\rm{variables}}),\\\
\bm{Y}=\left[{\begin{array}[]{*{20}c}{{\bm{y}}_{1}}\\\ {{\bm{y}}_{2}}\\\
\vdots\\\ {{\bm{y}}_{t}}\\\
\end{array}}\right]=\left[{\begin{array}[]{*{20}c}{y_{11}}&{y_{12}}&\ldots&{y_{1d}}\\\
{y_{21}}&{y_{22}}&\ldots&{y_{2d}}\\\ \vdots&\vdots&\ddots&{\rm{\ldots}}\\\
{y_{t1}}&{y_{t2}}&\ldots&{y_{td}}\\\ \end{array}}\right]\\\
\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;(t\;{\rm{samples}}\;d\;{\rm{variables}}),\\\
\end{array}$ (21)

CCA algorithm attempts to find canonical correlation vectors ${\bm{J}}$ and
${\bm{L}}$, which maximize the correlation $\rho$ between
${\bm{J}}^{\rm{T}}\bm{U}$ and ${\bm{L}}^{\rm{T}}\bm{Y}$, i.e.,

$\begin{gathered}({\bm{J}},{\bm{L}})=\mathop{\arg\;\max}\limits_{({\bm{J}},{\bm{L}})}\;\rho_{({\bm{J}}^{\rm{T}}\bm{U})({\bm{L}}^{\rm{T}}\bm{Y})}\hfill\\\
\;\;\;\;\;\;\;\;\;=\mathop{\arg\;\max}\limits_{({\bm{J}},{\bm{L}})}\;\frac{{{\bm{J}}^{\rm{T}}\bm{\Sigma}_{\bm{UY}}{\bm{L}}}}{{({\bm{J}}^{\rm{T}}\bm{\Sigma}_{\bm{U}}{\bm{J}})^{1/2}({\bm{L}}^{\rm{T}}\bm{\Sigma}_{\bm{Y}}{\bm{L}})^{1/2}}},\hfill\\\
\end{gathered}$ (22)

where $\bm{\Sigma}_{(\cdot)}$ denotes the covariance matrix. The solution to
problem (22) can be obtained by implementing singular value decomposition on
the matrix $\bm{K}$, which is shown as

${\bm{K}}=\bm{\Sigma}_{\bm{U}}^{-1/2}\bm{\Sigma}_{\bm{UY}}\bm{\Sigma}_{\bm{Y}}^{-1/2}={\bm{R}}\bm{\Sigma}{\bm{V}}^{\rm{T}},$
(23)

with
$\bm{R}=(\bm{r}_{1},...,\bm{r}_{p}),\;\bm{V}=(\bm{v}_{1},...,\bm{v}_{q}),\bm{\Sigma}=\left[{\begin{array}[]{*{20}c}{\bm{\Sigma}_{\kappa}}&0\\\
0&0\\\ \end{array}}\right]$, where $\bm{\Sigma}_{\kappa}={\text{
diag}}(\rho_{1},...,\rho_{\kappa})$, $\kappa$ denotes the number of nonzero
singular values and $\rho_{i}\;(i=1,...,\kappa)$ are canonical correlation
coefficients. $\bm{r}_{i}\;(i=1,...,p)$ and $\bm{v}_{i}\;(i=1,...,q)$ are
corresponding singular vectors.

The canonical correlation vectors can be derived as

${\bm{J}}=[{\bm{J}}_{1},...,{\bm{J}}_{p}]=\bm{\Sigma}_{\bm{U}}^{-1/2}{\bm{R}}\in\mathbb{R}^{p\times
d},$ (24)
${\bm{L}}=[{\bm{L}}_{1},...,{\bm{L}}_{d}]=\bm{\Sigma}_{\bm{Y}}^{-1/2}{\bm{V}}\in\mathbb{R}^{d\times
p}.$ (25)

Assume that PL $r$ and its both ends are corresponding to the physical path
${\text{PN}}_{m}\mathop{\to}\limits^{{\text{PL}}_{m,l}}{\text{PN}}_{l}$. To
detect the working state of PL $r$, i.e., ${\text{PL}}_{m,l}$ in a distributed
manner, we need to find each SFC, whose virtual path
${\text{VN}}_{m}\mathop{\to}\limits^{{\text{VL}}_{m,l}}{\text{VN}}_{l}$ is
embedded to this physical path. Then, the correlation of measurements between
${\text{VN}}_{m}$ and ${\text{VN}}_{l}$ needs to be analyzed in its Network
Slice Manager distributedly. Assume that the current measurements of
${\text{VN}}_{m}$ and ${\text{VN}}_{l}$ are represented by $\bm{u}$ and
$\bm{y}$, respectively, then the key step for analyzing the correlation of
measurements between neighbor VNs is to produce an anomaly detection residual,
which is shown as follows:

${\bm{r}}={\bm{J}}^{\rm{T}}\bm{u}-\bm{\Sigma}{\bm{L}}^{\rm{T}}\bm{y}.$ (26)

Then, the $T^{2}$ test for residual can be established as

$T_{r}^{2}={\bm{r}}^{\rm{T}}\bm{\Sigma}_{r}^{-1}{\bm{r}},\;\;{\text{where}}\;\bm{\Sigma}_{r}={\bm{I}}_{p}-\bm{\Sigma}\bm{\Sigma}^{\rm{T}}.$
(27)

Therefore, the correlation of measurements between neighbor VNs can be
determined based on the following decision logic as

$\left\\{{\begin{array}[]{*{20}c}{T_{r}^{2}\leqslant
T_{r,{\text{cl}}}^{2}\;\;\Rightarrow{\text{anomaly - free}}},\hfill\\\
{T_{r}^{2}>T_{r,{\text{cl}}}^{2}\;\Rightarrow{\text{anomaly}}},\hfill\\\
\end{array}}\right.$ (28)

where $T_{r,{\text{cl}}}^{2}$ is the control limit of $T_{r}^{2}$, which
represents the threshold that distinguishes between normal and anomalous
correlation patterns.

_Proposition 1_ : Because the residual vector
${\bm{r}}={\bm{J}}^{\rm{T}}{\bm{u}}-\bm{\Sigma}{\bm{L}}^{\rm{T}}{\bm{y}}$ has
the minimum covariance, it is an optimal residual for detecting the variation
of correlation between two sets of variables $\bm{u}$ and $\bm{y}$.

_Proof_ : See Appendix A.

Since new measurements will be generated in each period, to obtain the real-
time canonical correlation vectors ${\bm{J}}(t)$ and ${\bm{L}}(t)$ for the
online PL anomaly detection, we need to calculate covariances in formula (23)
with all the historical data at each time $t$. Its computational complexity
will become extremely high over time. To reduce the computational complexity
and the storage consumption, we derive a new method to update covariances in
formula (23), with which the Network Slice Manager only needs to store the
covariance matrices and mean vectors of current VN data instead of all the
historical data.

_Proposition 2_ : Assume that the covariance matrix of ${\bm{U}}(t)$ is
$\bm{\Sigma}_{\bm{U}(t)}=\frac{1}{{t-1}}\left[{\begin{array}[]{*{20}c}{a_{11}}&{a_{12}}&\ldots&{a_{1p}}\\\
{a_{21}}&{a_{22}}&\ldots&{a_{2p}}\\\ \vdots&\vdots&\ddots&\vdots\\\
{a_{p1}}&{a_{p2}}&\ldots&{a_{pp}}\\\ \end{array}}\right]$ and the mean vector
of ${\bm{U}}(t)$ is $[c_{1}(t),c_{2}(t),...,c_{p}(t)]$ at time $t$, if the
measurements of $\bm{u}$ is
${\bm{u}}_{t+1}=(u_{(t+1)1},u_{(t+1)2},...,u_{(t+1)p})$ at time $t+1$, the
covariance matrix of ${\bm{U}}(t+1)$ can be computed by

$\begin{gathered}\bm{\Sigma}_{\bm{U}(t+1)}=\hfill\\\
\frac{1}{t}\times\left[{\begin{array}[]{*{20}c}{a_{11}+m_{11}}&{a_{12}+m_{12}}&\ldots&{a_{1p}+m_{1p}}\\\
{a_{21}+m_{21}}&{a_{22}+m_{22}}&\ldots&{a_{2p}+m_{2p}}\\\
\vdots&\vdots&\ddots&\vdots\\\
{a_{p1}+m_{p1}}&{a_{p2}+m_{p2}}&\ldots&{a_{pp}+m_{pp}}\\\
\end{array}}\right],\\\ \end{gathered}$ (29)

where
$m_{i,j}=\frac{{t(c_{i}(t)-u_{(t+1)i})(c_{j}(t)-u_{(t+1)j})}}{{t+1}}(1\leq
i,j\leq p)$.

_Proof_ : See Appendix B.

Similarly, assume that the covariance matrix between $\bm{U}(t)$ and
$\bm{Y}(t)$ is
$\bm{\Sigma}_{\bm{U}(t)\bm{Y}(t)}=\frac{1}{{t-1}}\left[{\begin{array}[]{*{20}c}{e_{11}}&{e_{12}}&\ldots&{e_{1d}}\\\
{e_{21}}&{e_{22}}&\ldots&{e_{2d}}\\\ \vdots&\vdots&\ddots&\vdots\\\
{e_{p1}}&{e_{p2}}&\ldots&{e_{pd}}\\\ \end{array}}\right]$ and the mean vector
of $\bm{Y}(t)$ is $[f_{1}(t),f_{2}(t),...,f_{d}(t)]$ at time $t$, if the
measurements of ${\bm{y}}$ is
${\bm{y}}_{t+1}=(y_{(t+1)1},y_{(t+1)2},...,y_{(t+1)d})$ at time $t+1$, the
covariance matrix between $\bm{U}(t+1)$ and $\bm{Y}(t+1)$ can be computed by

$\begin{gathered}\bm{\Sigma}_{\bm{U}(t+1)\bm{Y}(t+1)}=\hfill\\\
\frac{1}{t}\times\left[{\begin{array}[]{*{20}c}{e_{11}+n_{11}}&{e_{12}+n_{12}}&\ldots&{e_{1d}+n_{1d}}\\\
{e_{21}+n_{21}}&{e_{22}+m_{22}}&\ldots&{e_{2d}+n_{2d}}\\\
\vdots&\vdots&\ddots&\vdots\\\
{e_{p1}+n_{p1}}&{e_{p2}+n_{p2}}&\ldots&{e_{pd}+n_{pd}}\\\
\end{array}}\right],\\\ \end{gathered}$ (30)

where
$n_{ij}=\frac{{t(c_{i}(t)-u_{(t+1)i})(d_{j}(t)-y_{(t+1)j})}}{{t+1}}(1\leq
i,(j)\leq p,(d))$.

Therefore, the covariance matrices at time $t+1$ can be calculated by only
keeping the covariance matrices and mean vectors at time $t$. Then, the
canonical correlation vectors ${\bm{J}}(t+1)$ and ${\bm{L}}(t+1)$ can be
derived to produce the optimal residual without storing all the historical
data, so large storage and computing resources can be saved.

The detailed steps of the CCA-based distributed online PL anomaly detection
algorithm are shown in Algorithm 2.

Algorithm 2 CCA-based distributed online PL anomaly detection algorithm

0: Initial number of labeled samples $t$, measurements $\bm{U}(t)$ and
$\bm{Y}(t)$ of ${\text{VN}}_{m}$ and ${\text{VN}}_{l}$ for each
${\text{VN}}_{m}\mathop{\to}\limits^{{\text{VL}}_{m,l}}{\text{VN}}_{l}$ mapped
to ${\text{PN}}_{m}\mathop{\to}\limits^{{\text{PL}}_{m,l}}{\text{PN}}_{m,l}$,
control limit $T_{r,{\text{cl}}}^{2}$ and number of iterations $T$

1: Compute the covariance matrices and mean vectors of $\bm{U}(t)$ and
$\bm{Y}(t)$: $\bm{\Sigma}_{\bm{U}(t)}$, $\bm{\Sigma}_{\bm{Y}(t)}$,
$\bm{\Sigma}_{\bm{U}(t)\bm{Y}(t)}$, $[c_{1}(t),c_{2}(t),...,c_{p}(t)]$ and
$[f_{1}(t),f_{2}(t),...,f_{d}(t)]$

2: for $t=t+1:T$ do

3: for each
${\text{VN}}_{l}\mathop{\to}\limits^{{\text{VL}}_{l,l+1}}{\text{VN}}_{l+1}$ do

4: Compute $\bm{\Sigma}_{\bm{U}(t)}$, $\bm{\Sigma}_{\bm{Y}(t)}$ and
$\bm{\Sigma}_{\bm{U}(t)\bm{Y}(t)}$ according to formulas (29) and (30)

5: Implement singular value decomposition on the matrix $\bm{K}(t)$ according
to formula (23)

6: Compute canonical correlation vectors $\bm{J}(t)$ and $\bm{L}(t)$ according
to formulas (24) and (25)

7: Produce the optimal anomaly detection residual ${\bm{r}}(t)$ according to
formula (26)

8: Establish the $T^{2}$ test:
$T_{r(t)}^{2}={\bm{r}}(t)^{\rm{T}}\bm{\Sigma}_{r(t)}^{-1}{\bm{r}}(t)$

9: end for

10: if $T_{r(t)}^{2}\leqslant T_{r,{\text{cl}}}^{2}$ for each
${\text{VN}}_{m}\mathop{\to}\limits^{{\text{VL}}_{m,l}}{\text{VN}}_{l}$ then

11: Determine that ${\text{PL}}_{m,l}$ (or PL $r$) is normal, and update
covariance matrices and mean vectors

12: else

13: Determine that ${\text{PL}}_{m,l}$ (or PL $r$) is anomalous, and reserve
covariance matrices and mean vectors at time $t-1$ and discard current ones

14: end if

15: end for

We usually deal with the anomaly detection as a special classification problem
by assuming that the entire dataset contains only normal samples. When the
number of anomalous data contained in training samples increases, the
performance of classical CCA method will decline greatly. By comparison, the
proposed algorithm is more robust for it can eliminate the impact of anomalous
data on the detection accuracy in each iteration (step 13).

The computation of CCA algorithm mainly focuses on the covariance calculation.
When the number of samples reaches $T$ and the number of VLs mapped to PL $r$
is $|J_{r}|$, the time and storage complexity of the classical CCA algorithm
are $O(T|J_{r}|(p^{2}+pd+d^{2}))$ and $O(T|J_{r}|(p+d+2))$, respectively. For
the distributed online mode, the time and storage complexity of the proposed
algorithm are $O(p^{2}+pd+d^{2})$ and $O(p+d)$, respectively.

##  VI Simulation Results

###  VI-A Evaluation Metrics

In this section, numerical simulations on both synthetic and real-world
network datasets are executed to evaluate the proposed distributed online PN
and PL anomaly detection algorithms. The performance analysis of the proposed
algorithms is evaluated using various metrics such as precision, recall and
f1-socre, which can be formulated by confusion matrix as presented in Table I.

TABLE I: Confusion matrix

| Predicted result
---|---
normal | anomalous
Actual result | normal | True Negative (TN) | False Positive (FP)
anomalous | False Negative (FN) | True Positive (TP)

TN, FN, FP and TP are defined as follows: TN represents the number of normal
working states which are identified correctly. FN indicates the number of
anomalies which are not identified. FP summarizes the normal working states
that have been judged as anomalies. TP represents the number of anomalies
which are identified correctly.

  * •

_Accuracy_ is the proportion of correctly predicted PN or PL working states to
the total ones.

$Accuracy=\frac{{TP+TN}}{{TP+TN+FP+FN}}$ (31)

  * •

_Precision_ represents the proportion of correctly predicted anomalous PN or
PL working states to the total predicted anomalous ones.

$Precision=\frac{{TP}}{{TP+FP}}$ (32)

  * •

_Recall_ represents the proportion of correctly predicted anomalous PN or PL
working states to the total actual anomalous ones.

$Recall=\frac{{TP}}{{TP+FN}}$ (33)

  * •

_F1-score_ is the weighted average of the precision and recall metrics.

$F1-score=2\times\frac{{precison\times recall}}{{precison+recall}}$ (34)

Accuracy is an effective evaluation metric when the datasets are balanced, but
for imbalanced ones it may give biased evaluation results. Compared to
accuracy, precision and recall are less biased metrics for the evaluation of
imbalanced datasets. Therefore, precision, recall and f1-score metrics are
mainly used to analyze the performance of the proposed anomaly detection
algorithms.

###  VI-B Synthetic Dataset

To generate synthetic data to validate the effectiveness of the proposed PN
and PL anomaly detection algorithms, a simple virtualized network slicing
scenario with 10 PNs and 6 SFCs has been established in MATLAB platform. Each
SFC is assumed to contain 4 $\sim$ 6 VNs [39]. Three different service
requests are contained to simulate the diverse service types in network
slicing. VNs and VLs in SFCs are randomly mapped to the substrate network to
provide end-to-end services for different user types. When PNs and PLs are in
normal states, their processing capacity is in normal range. To simulate the
anomalous cases, the loss rates of processing capacity for PNs and PLs, which
follows the Gaussian distribution $N(\mu,\sigma^{2})$ with the mean $\mu=0.5$
and the variance $\sigma^{2}=0.01$, are randomly injected into the networks.
Then, the simulated data can be acquired in each period, including processing
rate, data flow, queuing delay, processing delay, etc. Table II presents the
main parameters used for simulation.

TABLE II: Simulation parameters

Parameter | Value
---|---
PN number | 10
SFC number | 6
VN number in each SFC | 4$\sim$6
Service 1 (arrive rate, packet size) | (10packets/s, 200kbit/packets)
Service 2 (arrive rate, packet size) | (100packets/s, 10kbit/packets)
Service 3 (arrive rate, packet size) | (500packets/s, 1kbit/packets)
Random feature space | 100-dimensional

(a) Precision

(b) Recall

(c) F1-score

Figure 3: Performance comparison among different PN anomaly detection
algorithms on synthetic dataset.

For comparison, we have simulated the OCSVM based [12] and online OCSVM based
[29] anomaly detection methods. Fig. 3 shows the comparison results of
precision, recall and f1-score among the OCSVM, online OCSVM and the proposed
distributed online OCSVM (DO-OCSVM) algorithms. The presented results are
averaged total 500 Monte Carlo runs. To verify the impact of anomalous data
contained in training samples on the classical OCSVM algorithm, we have
introduced the concept of anomaly ratio in training data (ARTD). From Fig. 3,
compared with the classical OCSVM and the online OCSVM algorithms, the
proposed DO-OCSVM has a higher recall and f1-score, but its precision is a
little lower. This is because in the proposed DO-OCSVM algorithm, as long as
measurements of one VN mapped to a PN are detected as anomaly, this PN will be
determined to be anomalous. Therefore, the proposed algorithm can detect the
anomalous PNs more accurately, but the possibility of wrongly determining
normal PNs as anomalous ones also increases. Since we pay more attention to
the recall in anomaly detection, a slight drop in precision is acceptable.
Note that when the training samples contain anomalous data (ARTD=0.1), the
performance of the classical OCSVM in three metrics will be greatly degraded.
However, it is costly to obtain large amounts of accurate labeled data in
actual networks. Fortunately, the proposed DO-OCSVM algorithm can eliminate
the impact of anomalous data on the detection accuracy during the iteration
process without any labeled data, which makes it more economic and robust
anomaly detection method.

Fig. 4 shows the convergence process of $\bm{w}$ and $\rho$ in the proposed
DO-OCSVM algorithm. Since $\bm{w}$ is a 100-D vector, it is difficult to
depict all the components. Therefore, we randomly select 5 components for
illustration. From Fig. 4, we can see that $\bm{w}$ and $\rho$ in DO-OCSVM
algorithm nearly remain stable after 1000 iterations. Note that the
convergence process of parameters presented in Fig. 4 is basically consistent
with that of precision, recall and f1-score shown in Fig. 3.

Figure 4: Convergence process of $\bm{w}$ and $\rho$ in DO-OCSVM algorithm.

The performance comparison results between the CCA algorithm in [15], online
CCA algorithm and the proposed distributed online CCA (DO-CCA) algorithm for
PL anomaly detection are shown in Fig. 5. We set the initial number of labeled
samples $t=10$ and control limit $T_{r,{\rm{cl}}}^{2}=1$. The presented
results are averaged total 500 Monte Carlo runs. To verify the impact of
anomalous data in training samples on the classical CCA algorithm, we simulate
the condition when ARTD=0.1. Observing Fig. 5, the proposed DO-CCA has a
higher recall compared to the CCA and online CCA algorithms, which suggests
that proposed DO-CCA algorithm can predict the anomalous PLs more accurately.
Meanwhile, the precision of the proposed DO-CCA algorithm is a little lower,
implying that the possibility of wrongly determined normal PLs as anomalous
ones also has a little increase. Since we concentrate more on the performance
of recall in anomaly detection, a slightly small precision is acceptable. Note
that when the training samples contain anomalous data, the performance of the
classical CCA algorithm will decline greatly. Since the proposed DO-CCA
algorithm can eliminate the impact of anomalous data during the iteration
process, it is more adaptive.

(a) Precision

(b) Recall

(c) F1-score

Figure 5: Performance comparison among different PL anomaly detection
algorithms on synthetic dataset.

###  VI-C Real-World Network Dataset

Except validating the proposed distributed online anomaly detection algorithms
on the synthetic dataset from the simulated virtualized network slicing
scenario, we further evaluate their performance on the real-world network
dataset. The data rate performance measurements of NFV cloud native scaling
for a media application [40], which is selected from IEEEDataPort, include
CPU, memory, disk and network metrics. The system monitoring period, i.e., the
data sampling period was set to 15s. To evaluate the performance of different
PN anomaly detection algorithms, we introduce four cases to simulate the
anomalies in PNs, including endless loop in CPU, memory leak, Disk I/O fault
and network congestion.

For the purpose of further comparison, we also evaluate the performance of two
state-of-the-art classical anomaly detection methods on the four anomaly
cases, which are _k_ -nearest neighbor (_k_ -NN) method [41] and local outlier
factor (LOF) [42] algorithm. The comparison results of different anomaly
detection algorithms on four anomaly cases are summarized in Table III. The
results show that the online OCSVM and DO-OCSVM algorithms have comparable or
better detection performance than _k_ -NN and LOF. For example, under the case
of Disk I/O fault, the detection performance of LOF is significantly lower
than that of the other three algorithms. The reason is that anomalous disk
read or write rate is not much different from the normal one, so the local
density of monitored measurements has only a little fluctuation, which is
difficult to capture for the LOF algorithm. Besides, LOF has a better
performance on detecting endless loop in CPU, which indicates that LOF is more
sensitive to the fluctuant metrics. By contrast, our proposed DO-OCSVM based
PN anomaly detetion algorithm has the best detection performance on the four
anomaly cases.

TABLE III: The comparison results of different anomaly detection algorithms on
four anomaly cases

Anomalies |  _k_ -NN | LOF | Online OCSVM | Proposed DO-OCSVM
---|---|---|---|---
Recall | Precision | F1-score | Recall | Precision | F1-score | Recall | Precision | F1-score | Recall | Precision | F1-score
Endless loop in CPU | 0.951 | 0.975 | 0.963 | 0.905 | 0.960 | 0.932 | 0.984 | 0.961 | 0.972 | 1.00 | 1.00 | 1.00
Memory Leak | 0.903 | 0.939 | 0.923 | 0.894 | 0.849 | 0.871 | 0.976 | 0.934 | 0.955 | 0.988 | 0.969 | 0.978
Disk I/O fault | 0.979 | 0.977 | 0.976 | 0.588 | 0.600 | 0.594 | 0.976 | 0.931 | 0.953 | 0.979 | 0.960 | 0.969
Network congestion | 0.858 | 0.965 | 0.908 | 0.808 | 0.727 | 0.767 | 1.00 | 0.970 | 0.986 | 1.00 | 0.977 | 0.988

Figure 6: F1-score comparison among different PN anomaly detection algorithms
versus different sizes of real-world network dataset. Figure 7: CPU time
comparison among different PN anomaly detection algorithms versus different
sizes of real-world network dataset. Figure 8: Memory consumption comparison
among different PN anomaly detection algorithms versus different sizes of
real-world network dataset. Figure 9: F1-score comparison among different PL
anomaly detection algorithms versus different sizes of real-world network
dataset.

The f1-score, which is the weighted average of the precision and recall
metrics, is used to evaluate the overall performance of different anomaly
detection algorithms with different datatset sizes. The corresponding results
are shown in Fig. 6. Furthermore, we also record the CPU time and memory
consumption of each detection method with different dataset sizes, and the
results are shown in Figs. 7 and 8, respectively. It should be noted that for
the proposed distributed detection scheme, we only take the average CPU time
and memory consumption of every distributed manager into account. Observing
the simulation results depicted in Figs. 6-8, it is obvious that the online
OCSVM and DO-OCSVM algorithms take much less CPU time and memory consumption
than the classical _k_ -NN and LOF methods, and obtain better detection
performance. The f1-score of _k_ -NN and LOF methods increases as the dataset
size grows. By contrast, the online OCSVM and DO-OCSVM algorithms are not
sensitive to the dataset sizes. Besides, the required CPU time and memory
consumption of _k_ -NN and LOF methods increase exponentially as the dataset
size grows. For the online OCSVM and DO-OCSVM algorithms, the required memory
consumption remains nearly unchanged, and the required CPU time grows linearly
due to the online detection mode.

For validating the performance of the proposed DO-CCA based PL anomaly
detection algorithm on the real-world network dataset, in addition to the CCA
and online CCA algorithms, we also evaluate the performance of the principal
component analysis method (PCA)[43], which is the state-of-the-art classical
multivariate analysis method. The f1-score, the required CPU time and memory
consumption are used to evaluate the overall performance of different PL
anomaly detection algorithms. The corresponding results are shown in Figs.
9-11, respectively. Observing the simulation results depicted in Figs. 9-11,
we can see that the online CCA and DO-CCA based PL anomaly detection
algorithms take much less CPU time and memory consumption than the classical
PCA and CCA methods, and obtain better detection performance. The f1-score of
PCA and CCA methods increases, and the required CPU time and memory
consumption increase linearly as the dataset size grows. By contrast, for the
online CCA and DO-CCA algorithms, the f1-score, required CPU time and memory
consumption all remain nearly unchanged as the dataset size grows. Besides,
our proposed DO-CCA based PL anomaly detetion algorithm has the best detection
performance with the lowest CPU time and memory consumption for the
distributed online mode.

Figure 10: CPU time comparison among different PL anomaly detection algorithms
versus different sizes of real-world network dataset. Figure 11: Memory
consumption comparison among different PL anomaly detection algorithms versus
different sizes of real-world network dataset.

##  VII Conclusions

The accurate and rapid anomaly detection for PNs and PLs in substrate networks
is the prerequisite for ensuring the performance of virtualized network
slices. To realize the real-time anomaly detection for the substrate networks
with low communication and storage cost, a distributed online PN anomaly
detection algorithm was first proposed based on a decentralized OCSVM. This
detection algorithm could identify the working state of a PN through analyzing
the real-time measurements of each VN mapped to it in a distributed manner.
Then, we proposed a CCA-based distributed online algorithm to realize the PL
anomaly detection. This algorithm could infer the working state of a PL based
on the correlation of measurements between neighbor VNs, which were mapped to
both ends of the PL. Finally, the effectiveness and robustness of the proposed
distributed online anomaly detection algorithms were verified on both the
synthetic and real-world network datasets.

##  Appendix A Proof of Proposition 1

As ${\bm{J}}=\bm{\Sigma}_{\bm{U}}^{-1/2}{\bm{R}}$ and
$\bm{\Sigma}_{\bm{U}}^{-1/2}\bm{\Sigma}_{\bm{UY}}\bm{\Sigma}_{\bm{Y}}^{-1/2}={\bm{R}}\bm{\Sigma}{\bm{V}}^{\rm{T}}$,
then
${\bm{R}}^{\rm{T}}\bm{\Sigma}_{\bm{U}}^{-1/2}=\bm{\Sigma}{\bm{V}}^{\rm{T}}\bm{\Sigma}_{\bm{Y}}^{1/2}\bm{\Sigma}_{\bm{UY}}^{-1}$
and
${\bm{J}}^{\rm{T}}{\bm{u}}=\bm{\Sigma}{\bm{V}}^{\rm{T}}\bm{\Sigma}_{\bm{Y}}^{1/2}\bm{\Sigma}_{\bm{UY}}^{-1}{\bm{u}}$.
We can get that

$\begin{gathered}{\bm{r}}={\bm{J}}^{\rm{T}}{\bm{u}}-\bm{\Sigma}{\bm{L}}^{\rm{T}}{\bm{y}}\hfill\\\
\;\;\;=\bm{\Sigma}{\bm{V}}^{\rm{T}}\bm{\Sigma}_{\bm{Y}}^{1/2}\bm{\Sigma}_{\bm{UY}}^{-1}{\bm{u}}-\bm{\Sigma}{\bm{V}}^{\rm{T}}\bm{\Sigma}_{Y}^{-1/2}{\bm{y}}\hfill\\\
\;\;\;=\bm{\Sigma}{\bm{V}}^{\rm{T}}\bm{\Sigma}_{\bm{Y}}^{-1/2}(\bm{\Sigma}_{\bm{Y}}\bm{\Sigma}_{\bm{UY}}^{-1}{\bm{u}}-{\bm{y}}),\hfill\\\
\end{gathered}$ (35)

where ${\bm{\hat{y}}}=\bm{\Sigma}_{\bm{Y}}\bm{\Sigma}_{\bm{UY}}^{-1}{\bm{u}}$
is the least square estimation of $\bm{y}$ by $\bm{u}$. In this case, the
residual has the minimal covariance value, and the $T_{r}^{2}$ statistic,
which includes the inverse of the covariance matrix, has the best performance
for anomaly detection [38].

##  Appendix B Proof of Proposition 2

The covariance calculation processes of $\bm{U}(t)$ are as the following
steps.

Firstly, subtract average values from each column as

$\begin{gathered}\bar{\bm{U}}(t)=\left[{\begin{array}[]{*{20}c}{u_{11}-c_{1}(t)}&{u_{12}-c_{2}(t)}&\ldots&{u_{1p}-c_{p}(t)}\\\
{u_{21}-c_{1}(t)}&{u_{22}-c_{2}(t)}&\ldots&{u_{2p}-c_{p}(t)}\\\
\vdots&\vdots&\ddots&\vdots\\\
{u_{t1}-c_{1}(t)}&{u_{t2}-c_{2}(t)}&\ldots&{u_{tp}-c_{p}(t)}\\\
\end{array}}\right]\hfill\\\
\;\;\;\;\;\;\;\;=\left[{\begin{array}[]{*{20}c}{\bar{u}_{11}}&{\bar{u}_{12}}&\ldots&{\bar{u}_{1p}}\\\
{\bar{u}_{21}}&{\bar{u}_{22}}&\ldots&{\bar{u}_{2p}}\\\
\vdots&\vdots&\ddots&\vdots\\\
{\bar{u}_{t1}}&{\bar{u}_{t2}}&\ldots&{\bar{u}_{tp}}\end{array}}\right].\hfill\\\
\end{gathered}$ (36)

Then, the covariance of $\bm{U}(t)$ can be computed by

$\begin{gathered}\bm{\Sigma}_{\bm{U}(t)}=\frac{1}{{t-1}}(\bar{\bm{U}}(t))^{\rm{T}}\bar{\bm{U}}(t)\hfill\\\
\;\;\;\;\;\;\;\;\;\;=\frac{1}{{t-1}}\left[{\begin{array}[]{*{20}c}{a_{11}}&{a_{12}}&{\rm{\ldots}}&{a_{1p}}\\\
{a_{21}}&{a_{22}}&{\rm{\ldots}}&{a_{2p}}\\\ \vdots&\vdots&\ddots&\vdots\\\
{a_{p1}}&{a_{p2}}&{\rm{\ldots}}&{a_{pp}}\\\ \end{array}}\right].\\\
\end{gathered}$ (37)

The mean vector of $\bm{U}(t+1)$ is
$\left({c_{1}(t+1),...,c_{p}(t+1)}\right)=\left({\frac{{tc_{1}(t)+u_{(t+1)1}}}{{t+1}},...,\frac{{tc_{3}(t)+u_{(t+1)p}}}{{t+1}}}\right)$,
so

$\bar{\bm{U}}(t+1)=\left[{\begin{array}[]{*{20}c}{\bar{u}_{11}{\rm{+}}\frac{{c_{1}(t)-u_{(t+1)1}}}{{t+1}}}&\ldots&{\bar{u}_{1p}+\frac{{c_{p}(t)-u_{(t+1)p}}}{{t+1}}}\\\
{\;\;\;\;\;\;\;\;\;\;\;\vdots}&\ddots&{\;\;\;\;\;\;\;\;\;\;\;\vdots}\\\
{\frac{{t(u_{(t+1)1}-c_{1}(t))}}{{t+1}}}&\ldots&{\frac{{t(u_{(t+1)p}-c_{p}(t))}}{{t+1}}}\\\
\end{array}}\right].$ (38)

Assume that the covariance matrix of $\bm{U}(t+1)$ is
$\bm{\Sigma}_{\bm{U}(t+1)}=\frac{1}{t}\left[{\begin{array}[]{*{20}c}{b_{11}}&{b_{12}}&{\rm{\ldots}}&{b_{1p}}\\\
{b_{21}}&{b_{22}}&{\rm{\ldots}}&{b_{2p}}\\\ \vdots&\vdots&\ddots&\vdots\\\
{b_{p1}}&{b_{p2}}&{\rm{\ldots}}&{b_{pp}}\\\ \end{array}}\right]$, where

$\begin{gathered}b_{ij}=\sum\limits_{k=1}^{t}{\left({\bar{u}_{ki}{\rm{+}}\frac{{c_{i}(t)-u_{(t+1)i}}}{{t+1}}}\right)\left({\bar{u}_{kj}{\rm{+}}\frac{{c_{j}(t)-u_{(t+1)j}}}{{t+1}}}\right)}\\\
\;\;\;\;\;\;\;\;+\left({\frac{{t(u_{(t+1)i}-c_{i}(t))}}{{t+1}}}\right)\left({\frac{{t(u_{(t+1)j}-c_{j}(t))}}{{t+1}}}\right)\\\
\;\;\;=a_{ij}+\sum\limits_{k=1}^{t}{\bar{u}_{ki}}\frac{{c_{j}(t)-u_{(t+1)j}}}{{t+1}}+\sum\limits_{k=1}^{t}{\bar{u}_{kj}}\frac{{c_{i}(t)-u_{(t+1)i}}}{{t+1}}\\\
\;\;\;\;\;\;+t\left({\frac{{c_{i}(t)-u_{(t+1)i}}}{{t+1}}}\right)\left({\frac{{c_{j}(t)-u_{(t+1)j}}}{{t+1}}}\right)\hfill\\\
\;\;\;\;\;\;+\;\left({\frac{{t(u_{(t+1)i}-c_{i}(t))}}{{t+1}}}\right)\left({\frac{{t(u_{(t+1)j}-c_{j}(t))}}{{t+1}}}\right)\;\hfill\\\
\;\;\;=a_{ij}+\;\frac{{t(c_{i}(t)-u_{(t+1)i})(c_{j}(t)-u_{(t+1)j})}}{{t+1}}\hfill\\\
\;\;\;\;(\sum\limits_{k=1}^{t}{\bar{u}_{ki}}=0,\sum\limits_{k=1}^{t}{\bar{u}_{kj}}=0)(1\leq
i,j\leq p).\hfill\\\ \end{gathered}$ (39)

Therefore, the _Proposition_ 2 is proved.

## References

  * [1] J. Ordonez-Lucena _et al._ , “Network slicing for 5G with SDN/NFV: Concepts, architectures, and challenges,” _IEEE Commun. Mag._ , vol. 55, no. 5, pp. 80–87, May 2017.
  * [2] S. E. Elayoubi, S. B. Jemaa, Z. Altman, and A. Galindo-Serrano, “5G RAN slicing for verticals: Enablers and challenges,” _IEEE Commun. Mag._ , vol. 57, no. 1, pp. 28–34, Jan. 2019.
  * [3] A. Oi, D. Endou, T. Moriya, and H. Ohnishi, “Method for estimating locations of service problem causes in service function chaining,” in _Proc. IEEE Global Commun. Conf. (GLOBECOM)_ , Feb. 2016, pp. 1–6.
  * [4] F. Z. Yousaf, M. Bredel, S. Schaller, and F. Schneider, “NFV and SDN-Key technology enablers for 5G networks,” _IEEE J. Sel. Areas Commun._ , vol. 35, no. 11, pp. 2468–2478, Nov. 2017.
  * [5] Y. Yang _et al._ , “The stochastic-learning-based deployment scheme for service function chain in access network,” _IEEE Access_ , vol. 6, pp. 52 406–52 420, Sep. 2018.
  * [6] M. Miyazawa, M. Hayashi, and R. Stadler, “vNMF: Distributed fault detection using clustering approach for network function virtualization,” in _Proc. IFIP/IEEE Int. Symp. Integr. Netw. Manage. (IM)_ , Jul. 2015, pp. 640–645.
  * [7] V. Chandola, A. Banerjee, and V. Kumar, “Anomaly detection: A survey,” _ACM Comput. Surveys_ , vol. 41, no. 3, pp. 1–58, Jul. 2009.
  * [8] C. Callegari _et al._ , “A novel PCA-based network anomaly detection,” in _Proc. IEEE Int. Conf. Commun. (ICC)_ , Jul. 2011, pp. 1–5.
  * [9] K. Xie _et al._ , “Fast tensor factorization for accurate internet anomaly detection,” _IEEE/ACM Trans. Netw._ , vol. 25, no. 6, pp. 3794–3807, Dec. 2017.
  * [10] S. Liu and L. Xu, “An integrated model based on o-gan and density estimation for anomaly detection,” _IEEE Access_ , vol. 8, pp. 204 471–204 482, Nov. 2020.
  * [11] L. Li, J. Yan, H. Wang, and Y. Jin, “Anomaly detection of time series with smoothness-inducing sequential variational auto-encoder,” _IEEE Trans. Neural Netw. Learn. Syst._ , vol. 32, no. 3, pp. 1177–1191, Apr. 2021.
  * [12] J. Jiang and L. Yasakethu, “Anomaly detection via one class SVM for protection of SCADA systems,” in _Proc. Int. Conf Cyber-Enabled Distrib. Comput. Knowl. Discovery (CyberC)_ , Oct. 2013, pp. 82–88.
  * [13] F. Schuster, A. Paul, R. Rietz, and H. König, “Potentials of using one-class SVM for detecting protocol-specific anomalies in industrial networks,” in _Proc IEEE Symp. Ser. Comput. Intell._ , Dec. 2015, pp. 83–90.
  * [14] J. Cao _et al._ , “Maximum correntropy criterion-based hierarchical one-class classification,” _IEEE Trans. Neural Netw. Learn. Syst._ , vol. 32, no. 8, pp. 3748–3754, Aug. 2021.
  * [15] Z. Chen _et al._ , “Fault detection for non-gaussian processes using generalized canonical correlation analysis and randomized algorithms,” _IEEE Trans. Ind. Electron._ , vol. 65, no. 2, pp. 1559–1567, Feb. 2018.
  * [16] D. Cotroneo, R. Natella, and S. Rosiello, “A fault correlation approach to detect performance anomalies in virtual network function chains,” in _Proc. IEEE Int. Symp. Softw. Rel. Eng. (ISSRE)_ , Nov. 2017, pp. 90–100.
  * [17] C. O’Reilly, A. Gluhak, and M. A. Imran, “Distributed anomaly detection using minimum volume elliptical principal component analysis,” _IEEE Trans. Knowl. Data Eng._ , vol. 28, no. 9, pp. 2320–2333, Sep. 2016.
  * [18] Z. Chen _et al._ , “A distributed canonical correlation analysis-based fault detection method for plant-wide process monitoring,” _IEEE Trans. Ind. Inform._ , vol. 15, no. 5, pp. 2710–2720, May 2019.
  * [19] X. Li _et al._ , “Online internet anomaly detection with high accuracy: A fast tensor factorization solution,” in _Proc. IEEE Conf. Comput. Commun. (INFOCOM)_ , Jun. 2019, pp. 1900–1908.
  * [20] K. Xie _et al._ , “On-line anomaly detection with high accuracy,” _IEEE/ACM Trans. Netw._ , vol. 26, no. 3, pp. 1222–1235, Jun. 2018.
  * [21] J. Zhang _et al._ , “Boosting positive and unlabeled learning for anomaly detection with multi-features,” _IEEE Trans. Multimedia_ , vol. 21, no. 5, pp. 1332–1344, May 2019.
  * [22] Y. Liu _et al._ , “Generative adversarial active learning for unsupervised outlier detection,” _IEEE Trans. Knowl. Data Eng._ , vol. 32, no. 8, pp. 1517–1528, Aug. 2020.
  * [23] B. Schölkopf _et al._ , “Estimating the support of a high-dimensional distribution,” _Neural Comput._ , vol. 13, no. 7, pp. 1443–1471, Jul. 2001\.
  * [24] S. Garg, K. Kaur, N. Kumar, and J. J. Rodrigues, “Hybrid deep-learning-based anomaly detection scheme for suspicious flow detection in SDN: A social multimedia perspective,” _IEEE Trans. Multimedia_ , vol. 21, no. 3, pp. 566–578, Mar. 2019.
  * [25] Q. Jiang and X. Yan, “Multimode process monitoring using variational Bayesian inference and canonical correlation analysis,” _IEEE Trans. Autom. Sci. Eng._ , vol. 16, no. 4, pp. 1814–1824, Oct. 2019.
  * [26] L. Ma, J. Dong, K. Peng, and C. Zhang, “Hierarchical monitoring and root-cause diagnosis framework for key performance indicator-related multiple faults in process industries,” _IEEE Trans. Ind. Informatics_ , vol. 15, no. 4, pp. 2091–2100, Apr. 2019.
  * [27] K. Zhang _et al._ , “A correlation-based distributed fault detection method and its application to a hot tandem rolling mill process,” _IEEE Trans. Ind. Electron._ , vol. 67, no. 3, pp. 2380–2390, Mar. 2020.
  * [28] Z. Chen _et al._ , “A just-in-time-learning-aided canonical correlation analysis method for multimode process monitoring and fault detection,” _IEEE Trans. Ind. Electron._ , vol. 68, no. 6, pp. 5259–5270, Apr. 2021.
  * [29] V. Gómez-Verdejo, J. Arenas-García, M. Lázaro-Gredilla, and Á. Navia-Vázquez, “Adaptive one-class support vector machine,” _IEEE Trans. Signal Process._ , vol. 59, no. 6, pp. 2975–2981, Jun. 2011\.
  * [30] Y. He, Y. Peng, S. Wang, and D. Liu, “Admost: UAV flight data anomaly detection and mitigation via online subspace tracking,” _IEEE Trans. Instrum. Meas._ , vol. 68, no. 4, pp. 1035–1044, Apr. 2019.
  * [31] L. Zhang, J. Zhao, and W. Li, “Online and unsupervised anomaly detection for streaming data using an array of sliding windows and PDDs,” _IEEE Trans. Cybern._ , vol. 51, no. 4, pp. 2284–2289, Apr. 2021.
  * [32] C. Wang _et al._ , “A distributed anomaly detection system for in-vehicle network using HTM,” _IEEE Access_ , vol. 6, pp. 9091–9098, Mar. 2018.
  * [33] J. Shu _et al._ , “Collaborative intrusion detection for VANETs: A deep learning-based distributed SDN approach,” _IEEE Trans. Intell. Transp. Syst._ , vol. 22, no. 7, pp. 4519–4530, Jul. 2021.
  * [34] X. Miao, Y. Liu, H. Zhao, and C. Li, “Distributed online one-class support vector machine for anomaly detection over networks,” _IEEE Trans. Cybern._ , vol. 49, no. 4, pp. 1475–1488, Apr. 2019.
  * [35] A. De la Oliva _et al._ , “5G-TRANSFORMER: Slicing and orchestrating transport networks for industry verticals,” _IEEE Commun. Mag._ , vol. 56, no. 8, pp. 78–84, Aug. 2018.
  * [36] P. A. Forero, A. Cano, and G. B. Giannakis, “Consensus-based distributed support vector machines,” _J. Mach. Learn. Res._ , vol. 11, pp. 1663–1707, Mar. 2010.
  * [37] A. Rahimi and B. Recht, “Random features for large-scale kernel machines,” in _Proc. Adv. Neural Inf. Process. Syst._ , Dec. 2007, pp. 1177–1184.
  * [38] Q. Jiang, S. X. Ding, Y. Wang, and X. Yan, “Data-driven distributed local fault detection for large-scale processes based on the GA-regularized canonical correlation analysis,” _IEEE Trans. Ind. Electron._ , vol. 64, no. 10, pp. 8148–8157, Oct. 2017.
  * [39] X. Fu _et al._ , “Service function chain embedding for NFV-enabled IoT based on deep reinforcement learning,” _IEEE Commun. Mag._ , vol. 57, no. 11, pp. 102–108, Nov. 2019.
  * [40] P. Melas and S. Taylor, “Data rate performance measurements of NFV cloud native scaling for a media application,” IEEE Dataport, May 2020. [Online]. Available: https://dx.doi.org/10.21227/2e1c-rd87
  * [41] S. Ramaswamy, R. Rastogi, and K. Shim, “Efficient algorithms for mining outliers from large data sets,” in _Proc. ACM Int. Conf. Manag. Data (SIGMOD)_ , vol. 29, Jun. 2000, pp. 427–438.
  * [42] M. Breunig, H.-P. Kriegel, R. Ng, and J. Sander, “LOF: Identifying density-based local outliers,” in _Proc. ACM Int. Conf. Manag. Data (SIGMOD)_ , vol. 29, Jun. 2000, pp. 93–104.
  * [43] S. Ding _et al._ , “On the application of PCA technique to fault diagnosis,” _Tsinghua Sci. Technol._ , vol. 15, no. 2, pp. 138–144, Apr. 2010.