# Event2Graph: Event-driven Bipartite Graph for
Multivariate Time-series Anomaly Detection

Yuhang Wu, Mengting Gu, Lan Wang, Yusan Lin, Fei Wang, Hao Yang  Visa
Research, VisaPalo Alto, CAUnited States

###### Abstract.

Modeling inter-dependencies between time-series is the key to achieve high
performance in anomaly detection for multivariate time-series data. The de-
facto solution to model the dependencies is to feed the data into a recurrent
neural network (RNN). However, the fully connected network structure
underneath the RNN (either GRU or LSTM) assumes a static and complete
dependency graph between time-series, which may not hold in many real-world
applications. To alleviate this assumption, we propose a dynamic bipartite
graph structure to encode the inter-dependencies between time-series. More
concretely, we model time series as one type of nodes, and the time series
segments (regarded as event) as another type of nodes, where the edge between
two types of nodes describe a temporal pattern occurred on a specific time
series at a certain time. Based on this design, relations between time series
can be explicitly modelled via dynamic connections to event nodes, and the
multivariate time-series anomaly detection problem can be formulated as a
self-supervised, edge stream prediction problem in dynamic graphs. We
conducted extensive experiments to demonstrate the effectiveness of the
design.

Anomaly detection, multivariate time-series, graph neural networks

##  1\. Introduction

Detecting anomalies in time-series data has been an important problem in the
research community of data mining as well as the finance industry. In many
circumstances, anomaly patterns in multiple time-series need to be taken into
account together to disclose the full picture of the system. For example,
before a financial crisis taking place, multiple macro and micro economic
indicators can get aberrant in a sequential manner. It is important to analyze
their transition and anomaly patterns all together rather than treat each
individual signal separately. In another smaller scope example, to model a
merchant default risk, acquirer banks needs to collect multiple key business
indicators such as cash flow, asset, liabilities, etc. To comprehensively
monitor the financial health of the merchant, anomalies in multiple key
indicators need to form up a story-line to describe the default probability of
the merchant.

Previous works in multivariate anomaly detection mainly rely on recurrent
neural networks (RNNs). Malhotra et al.(DBLP:journals/corr/MalhotraRAVAS16, )
and Hundman et al.(DBLP:conf/kdd/HundmanCLCS18, ) employed LSTM models
(DBLP:journals/neco/HochreiterS97, ) to capture the temporal dependencies of
multivariate signals and adopted prediction and reconstruction errors,
respectively, to identify the anomaly patterns. Su et
al.(DBLP:conf/kdd/SuZNLSP19, ) improved the classical RNN model by modeling
the probability distribution of time series via incorporating stochastic
variables. To further model the correlations of time-series explicitly, Zhao
et al.(zhao2020multivariate, ) proposed a graph attention network to propagate
information from different time series and aggregate the information together
before feeding into a GRU (DBLP:journals/corr/ChungGCB14, ). However, because
of the underlying RNN structure, previous works assume a static, complete
dependency graph among time series. These approaches may not perform well
under regime change of time series where the underlying inter-dependencies are
different. Our work is built upon the recent wisdom of dynamic graph neural
network. We allow the connectivity of the graph changes dynamically in
different time stamps based on the patterns on the time series. We believe
this adaptive design is more realistic and flexible in the real-world.
Experiments also show that this design leads to improved performance in
multivariate anomaly detection.

Figure 1. Overview of the event-driven bipartite graph structure. This toy
example has six time-series nodes and three event nodes. Time-step
$\textbf{t}_{90}$ correspond to an anomaly pattern while other time stamps are
normal. The relation between time sries nodes and event nodes are highlighted,
green dashed line means the predicted relation while the red line means the
actual relation. Figure 2. The inter-dependency map of multi-variate time
series assumed by different time series forecasting models. (a) No inter-
dependency assumed by uni-variate regression (e.g., ARIMA
(DBLP:conf/imc/ZhangGGR05, )), (b) Static and complete inter-dependency
assumed by RNN (DBLP:journals/corr/MalhotraRAVAS16, ;
DBLP:conf/kdd/HundmanCLCS18, ; DBLP:journals/neco/HochreiterS97, ;
DBLP:conf/kdd/SuZNLSP19, ; zhao2020multivariate, ), (c) static and partial
inter-dependency assumed by graph structure learning (deng2021graph, ), (d)
dynamic inter-dependency assumed by the proposed approach.

In order to construct a dynamic graph on-the-fly, one important question is
how to determine the connectivity of the graph in each timestamp. Previous
works such as graph structural learning provide a feasible solution to
constructing a static latent graph from time-series, but to the best of our
knowledge, no previous work covers how to learn a dynamic and explainable
graph for multivariate anomaly detection. Our work is inspired by the recent
progress of evolutionary event graph (DBLP:conf/aaai/Cheng00HZS20, ;
DBLP:conf/wsdm/Hu0CYR21, ) where the nodes in the graph represent the time-
sequences segments (events) and directed links represent the transition of the
segments (events). Compare to previous works, this line of research naturally
models the time-varying relations among time-series states via dynamic
connections, and each state carries a physical meaning that is understandable
by human. However, one major limitation of (DBLP:conf/wsdm/Hu0CYR21, ) is that
the event nodes employed in this work capture the information across all the
time series. Assume there are $K$ segment patterns in each time series, and
the number of time series is $D$. The model would need $K^{D}$ number of event
nodes to represent the multivariate signals. This exponential number of event
nodes strongly limits the information processing capability of the
evolutionary state graph and therefore allows only a very small number of time
series ($D$) to be analyzed in practice (four in the dataset used in the paper
(DBLP:conf/wsdm/Hu0CYR21, )).

To address the problem of exponential number of combinations, we disentangle
the time-series nodes and the event nodes in our design, and model them as two
types of nodes in a dynamic bipartite graph (as depicted in Fig. 1). Each
event node only represents a time segment on one individual time series,
instead of integrating patterns across all time series. The undirectional
connection between two types of nodes indicates event $e$ happens on the
$d^{th}$ time series at time $t$. So the maximum number of edges in the graph
is $O(KD)$, which is much smaller than $O(K^{2D})$. To further improve the
efficiency and generalizability of the algorithm, we built upon the framework
based on the recent advances in edge streams (DBLP:conf/kdd/KumarZL19, ;
xu2020inductive, ; DBLP:journals/corr/abs-2006-10637, ), where connections
between nodes are modeled as incoming attributed edges instead of constructing
adjacency matrices. The complete system, with the name Event2Graph (Event-
driven bipartite Graph), outperforms previous state-of-the-art on three public
data-sets, and we summarize our main contributions as follows:

(i) We propose a bipartite event-graph-based system to analyze the
multivariate time series and model the interactions between time series and
event segments via edge streams. This design significantly reduces the
complexity of event encoding in multivariate time-series, and establishes
state-of-the-art results in multivariate time-series anomaly detection.
(ii) Events employed in the system is highly interpretable and easy for human
to interact with. Through backtracing the historical event series or
predicting future event series, the system can provide explainable root cause
analysis and scenario-based forecasting.

##  2\. Related work

Most of the existing works in anomaly detection focus on univariate time
series. In multivariate anomaly detection, based on how the methods model the
dependencies among time series, they can be classified into two major
categories, namely: “static inter-dependencies” approaches, and “dynamic
inter-dependencies” approaches. In this section, we introduce these two
categories of methods in detail. We also summarize the recent progress in
dynamic graph analysis, from where we derive the insight of the proposed
solution. A schematic view is presented in Fig. 2 to introduce these methods.

###  2.1. Static inter-dependency relation

We argue that most of the existing multivariate anomaly detection models
belong to this category of approach, where relations (e.g., correlation)
between multivariate time series are fixed once learned, and the model assumes
all the time series can influence each other (complete inter-dependency).
LSTM-based framework has been widely employed in this category of work. For
instance, Malhotra et al.(DBLP:journals/corr/MalhotraRAVAS16, ) proposed a
LSTM-based auto-encoder network to detect anomalies from multiple sensors, Su
et al.(DBLP:conf/kdd/SuZNLSP19, ) proposed a stochastic LSTM framework to
model the data distribution of time series in a global perspective. As RNN was
originally proposed to model the temporal dependencies between different
timestamps, and all the dimensions of the input of RNN are used to describe a
single concept (e.g.word) all together, it was not tailored to model the
inter-dependencies among variables by design. Recently, Zhao et
al.(zhao2020multivariate, ) proposed a graph attention network-based (GAT)
approach (DBLP:conf/iclr/VelickovicCCRLB18, ) where each time series is
regarded as an individual node, and information are aggregated based on the
underlying similarity of the signals. While this solution partially mitigates
the problem by taking into account of dynamic pairwise similarities between
time series, the design assumes a complete and static inter-dependency graph.
Also, the processed information after GAT is simply aggregated all together
and fed into a single GRU. Hence the module still suffers from similar
problems as previous designs.

###  2.2. Dynamic inter-dependency relation

A few recent works started to explore the partial inter-dependency relations
among time series. Deng et al.(deng2021graph, ) constructed a neighborhood
graph on-the-fly based on sensor embeddings to describe the dependencies
between different sensors. The sparsity level of this incomplete graph can be
customized by users, but the connectivity of the graph is fixed once
constructed. Hu et al.(DBLP:conf/wsdm/Hu0CYR21, ) proposed the first dynamic
dependency design in multivariate anomaly detection. In this design, human
interpretable time-series segments are used to create the nodes, and the
transition of segments are explicitly modeled by dynamic graph neural network.
This design allows the system to represent the time-varying relations among
time series. However, one major limitation of the method is it uses a single
event node to represent segments across all time series. This design brings
the combination explosion problem in pattern representation (mentioned in
Section I) and makes it hard to tackle a moderate scale problem.

###  2.3. Dynamic Graph Neural Networks

Our model is built upon the advances in dynamic graph neural networks. Most of
the literatures in dynamic graph networks assume discrete-time graph dynamics
where the graphs are represented as a sequence of snapshots
(DBLP:conf/ijcai/YuCACW17, ; DBLP:conf/kdd/YuCAZCW18, ;
DBLP:conf/wsdm/SankarWGZY20, ). Recent works start to explore more flexible
design and assumes edges can appear at any time (xu2020inductive, ;
DBLP:conf/bigdataconf/NguyenLRAKK18, ; DBLP:conf/kdd/KumarZL19, ;
DBLP:conf/iclr/TrivediFBZ19, ; DBLP:journals/corr/abs-2006-10637, ). Our
solution is derived from the state-of-the-art solution named Temporal Graph
Networks (TGN) proposed by Rossi et al.(DBLP:journals/corr/abs-2006-10637, ).
The model is built upon the temporal graph attention network (TGAT)
(xu2020inductive, ) but extended with an unique node-wise memory which attempt
to model the temporal dependency on the node level instead of the system level
as (zhao2020multivariate, ). This allow us to model the dynamic inter-
dependency in our bi-partite event graph accurately, and provide us the
flexibility to incorporate new nodes that have not seen in training.

Table 1. Notation. X | Input time series
---|---
$D$ | Total number of time series sequences
$d$ | The $d^{th}$ time series sequence
$K$ | Total number of events
$\beta$ | Stride size of sliding window
$\tau$ | Length of sliding window
$\dot{t}$ | A single time step
$t$ | A time window
$\textbf{v}_{e}$ | Event node
$\textbf{v}_{m}$ | Time series sequence node
$\mathbb{G}^{t}_{B}$ | Bipartite event graph at time step $t$
$A$ | Attributed edge between $\textbf{v}_{e}$ and $\textbf{v}_{m}$
|

##  3\. Problem definition

Problem definition:  Following (DBLP:conf/wsdm/Hu0CYR21, ), we assume the
multivariate time-series composes multiple univariate time-series from the
same entity, and the target problem can be defined as follows: A multivariate
time-series is defined as $\textbf{X}\in\mathbb{R}^{T\times d}$ where $T$ is
the maximum number of time stamps in the time-series and $d$ is the dimension.
The anomaly detection problem is to detect specific time stamps
$\dot{t}^{*}\in\mathbb{O}$ in a multivariate time series where the time series
behavior deviating from the normal patterns of the time series. Set
$\mathbb{O}$ contains the timestamps that marked as anomaly by domain expert.
In case $T$ is large, a common way to handle long time-series is to use a
sliding window approach. Let $\tau$ denotes the length of the sliding window.
We hereby formulate the problem as a binary classification problem with the
objective to identify time windows with
$\textbf{X}^{[\dot{t}-\tau:\dot{t}]\times d}$ that contain anomaly time-
stamps.

Background: One solution of solving the multivariate time-series anomaly
detection problem is to convert time series into a homogeneous graph
(DBLP:conf/aaai/Cheng00HZS20, ; DBLP:conf/wsdm/Hu0CYR21, ), where
$\mathbb{G}_{e}$ is defined on the representative patterns of a multivariate
time series. The node of the graph $\textbf{v}_{e}$ corresponds to a human
interpretable time series patterns $\textbf{p}\in\mathbb{R}^{\tau\times d}$.
The pattern should be representative enough so that the whole time series can
be approximated by transitions of symbolic events (states) in the graph. The
transitional relation between two sequential nodes (end with stamp $t$ and
$t+1$) is defined as an edge in the graph. Many existing works (e.g., Bag of
Patterns, Shaples, sequence clustering) can be used to distill representative
patterns (event) from time series, but most of them only effective on
univariate time series where $d=1$.

##  4\. METHODOLOGY

In this section, we present the proposed framework in detail. This section is
organized as follows: In subsection 4.1, we present the structure of our
Event2Graph system and introduce how to use it for anomaly detection. Then, in
subsection 4.2, we introduce how events nodes are created and added in the
graph given a multivariate time series. Subsection 4.3 describes the network
design to digest the information flow generated by Event2Graph. Finally, in
Section 4.4, the optimization and inference methods are detailed.

###  4.1. Bipartite Event Stream

The core intuition behind the dynamic bipartite design is to decouple three
key concepts in the traditional event graph, namely: where (which time
series), when (at what time), and which (the event category). This practice
helps us avoid the problem of exponential pattern combinations.

Specifically, the proposed algorithm is built upon a dynamic bipartite graph.
We define it as a sequence of undirected bipartite event graph
$\mathbb{G}^{t}_{B}=\\{(\textbf{v}_{m}^{t},\textbf{v}_{e}^{t},A(\textbf{v}_{m}^{t},\textbf{v}_{e}^{t}))\\}$.
Here $t$ represents when the time window that the event graph is formulated.
If the window size is set to be 1, then the time window $t$ is essentially the
same as the time step $\dot{t}$. A time sequence index node
$\textbf{v}_{m}^{t}$ indicates “where” an event $\textbf{v}_{e}^{t}$ is
happening. An attributed edge $A(\textbf{v}_{m}^{t},\textbf{v}_{e}^{t})$ that
connect event node $\textbf{v}_{e}^{t}$ and sequence index node
$\textbf{v}_{m}^{t}$ indicates an event $\textbf{v}_{e}^{t}$ happened on time
series $m$ at time window $t$. For simplicity, we also denote the edge as
$A_{m,e}$ in the following sections. We employ an edge stream representation
so an edge $A(\textbf{v}_{m}^{t},\textbf{v}_{e}^{t})$ is only constructed to
represent the relation that actually existed. A major benefit of using edge
stream representation over adjacency matrix is that it allows the graph
structure to be more scalable and flexible, which provides us the generality
of incorporating new events that have not appeared in training.

Based on this bipartite graph structure, we convert the multivariate anomaly
detection problem into a self-supervised edge prediction problem in dynamic
graph. Given the historical sequence of $\mathbb{G}_{B}^{(1:t)}$ and event
nodes $\textbf{v}_{e}^{t+1}$ in time window $t+1$, we are predicting edges
$\hat{A}(\textbf{v}_{m}^{t+1},\textbf{v}_{e}^{t+1}))$ in
$\mathbb{G}_{B}^{t+1}$. The anomaly is derived as a mismatching score from the
predicted edge set $\hat{A}(\textbf{v}_{m}^{t+1},\textbf{v}_{e}^{t+1}))$ and
the observed edge set $A(\textbf{v}_{m}^{t+1},\textbf{v}_{e}^{t+1}))$ with a
read-out function $r(\cdot)$.

So the procedure of solving the multivariate anomaly detection problem
becomes:

  1. (1)

Given a normal sequence $\textbf{X}_{tr}$, identify representative patterns on
each of the time sequence respectively in multivariate time series.

  2. (2)

With identified representative patterns events, merge similar events across
time-series to indicate affinity relation between time series.

  3. (3)

Build a sequence of bipartite event graph for multivariate time series (with a
stride $\beta$) based on event matching.

  4. (4)

Analyze the sequence of bi-partite event graphs and derive a model $\Psi$ to
describe the intra- and inter-dependency relations.

  5. (5)

Given a testing sequence $\textbf{X}_{te}$, repeat (3), apply the model
learned in (4) to predict $\hat{A}$ at time window $t^{*}$.

  6. (6)

Derive anomaly scores based on the predicted $\hat{A}$ and original time
series with a proposed readout function $r(\cdot)$.

###  4.2. Event Node Detection And Matching

We defined the event nodes on each time series separately and employed an
unsupervised algorithm to identify the events. As we have no prior knowledge
about what an anomaly pattern would look like in most real world applications,
the system is designed to learn the representative patterns from observed
normal data.

Representative segment detection: We adopted Matrix Profile
(DBLP:conf/icdm/YehZUBDDSMK16, ), a state-of-the-art unsupervised algorithm to
identify representative patterns from time sequences. Matrix profile is able
to identify the top-$K$ repeated patterns on time series with high accuracy
and low computation time. We employed a recent implementation of matrix
profile called SCRIMP++ (DBLP:conf/icdm/0014YZKK18, ) on each of the time
series in $\textbf{X}_{tr}$, and the algorithm yields a list of single
dimensional representative patterns $\textbf{p}_{m,k}$ ($m\in\\{1:D\\}$,
$k\in\\{1:K\\}$), each with size $\tau$.

Event generation: many time series share similar representative patterns with
each other. For example, many products sold in a supermarket may share similar
seasonality behaviors with other products during a year, even though the
actual amount may vary. Based on this observation, we argue that it is
reasonable to merge similar segment patterns together across different time
series to create the event nodes. Another advantage of merging detected
segments is it helps to build the affinity relation cross time series. If two
time series always share similar patterns, a spike happens in one time series
that break the pair-wise correlation may be highly possible due to a true
anomaly. To measure the similarity between time-series, we employed dynamic
time wrapping (DBLP:conf/kdd/BerndtC94, ), which provides us with a robust
distance matrix that is insensitive to small misalignments between two time
series segments. After that, we employed the H-DBSCAN
(DBLP:journals/jossw/McInnesHA17, ) to cluster the patterns $\textbf{p}_{m,k}$
into clusters. We select the centroid of each cluster to be an event, which is
a representative segment across all time-series within the cluster. Finally,
we obtain an event set $\mathbb{E}_{\tau}$ which contains all events each with
length $\tau$.

Figure 3. Three group of event bars. Events are represented by colorful grids
in the event bar. In each event bar, X axis corresponding to the time, each
row corresponds a time-series. One event bar summarizes all the event happened
on an individual time-series. For each group of event bar, two type event bars
corresponding to the events generated by pattern matching and residual
computing. The predicted event bars are generated by the proposed system via
forecasting the next connected event node for each time-series. The anomaly
score is generated based on the mismatch between the predicted events and the
actual events with a carefully designed score function. Figure 4. Overview of
the proposed approach. The time-series nodes are represented by blue circles,
two types of event nodes (results from time series matching and residual
computation) are represented by yellow and orange circle, respectively.

Event matching: after extracting $K$ events from time series segments, we
match each event with the original time series to identify where and when each
event is taking place. We employed the highly efficient C language-based
Dynamic Time Wrapping algorithm (FastestSimilaritySearch, ) to conduct the
event matching, which provides a $(T-\tau)\backslash\beta\times D\times K$
similarity tensor $\mathcal{L}_{\tau}$ to indicate the similarity of each
event to all the sliding window with size $\tau$ and stride $\beta$ in the
input time sequence.

Event nodes created by pattern matching: for a single time stamp $t$, given
$\mathcal{L}_{\tau}$ the event graph is created by selecting the best matched
events for each time series. The event node is then connected with the
corresponding time-series node $\textbf{v}_{m}^{t}$ in the event graph
$\mathbb{G}^{t}_{B}$ with an attribute edge
$A(\textbf{v}_{e}^{t},\textbf{v}_{m}^{t})$. By finding and linking the best
matched event for each time series, we added $d$ number of event nodes into
the graph.

Event nodes created by residual error: we noticed that even for the best
matched time series, the pattern matching still left with small residual
errors. We hereby define two general residual nodes which indicate whether the
residual error in a time-series is larger than a threshold $\theta$. One
residual event node denoted as $\textbf{v}_{e^{+}}^{t}$ indicates the residual
error is larger than $\theta$, another residual event node denoted as
$\textbf{v}_{e^{-}}^{t}$ indicates the residual error is equal or smaller than
$\theta$. We learn $\theta$ in a data-driven manner with SPOT
(siffer2017anomaly, ) algorithm, where we employ the whole training data-set
for initialization and testing for on-going adaptation. All time-series shared
these two residual nodes as shown in Fig. 3.

After generating a sequence of bipartite event graph $\mathbb{G}^{(1:t)}_{B}$
for time series $\textbf{X}_{tr}$, a model $\Psi$ is trained on
$\mathbb{G}^{(1:t)}_{B}$ so that given a testing sequence $\textbf{X}_{te}$,
it is able to predict the connectivity of an event graph in any specific time
stamp $\dot{t}$ by observing historical sequence observed on
$\textbf{X}_{te}$. In Fig. 5, we visualize the events forecasted/detected on
ten time-series. The difference between the top and the middle group result in
the bottom group where we can easily identify the anomaly locations.

###  4.3. Network Design With Node-wise Memory

For each pair of $\textbf{v}_{m}^{t}$ and $\textbf{v}_{e}^{t}$ in graph
$\mathbb{G}^{(t)}_{B}$, the node features of them are defined as $\bm{v}_{m}$
and $\bm{v}_{e}$, respectively. We also define $\bm{\epsilon}_{m,e}$ as the
edge features between $\textbf{v}_{m}$ and $\textbf{v}_{e}$. In order to
support the dynamic inter-dependency graph, we avoid explicitly defining a
static adjacency matrix to describe the relation between nodes. Instead, we
adopt a state-message framework to model the node interactions. For each node
$\textbf{v}_{m}$ at time stamp $t$ (here we use time series node $m$ as an
example; the same rule applies to the event node $e$), we define a state
vector $\bm{s}_{m}(t)$ to represent its interaction history with other
$\textbf{v}_{e}$ nodes before $t$ in a compressed format. By initiating
$\bm{s}_{m}(0)$ as an all zero vector, the interaction at time $t$ is encoded
with a message vector $\bm{\varrho}_{m}(t)$:

(1) $\displaystyle\bm{\varrho}_{m}(t)=[\bm{\epsilon}_{m,e}(t)||\Delta
t||\bm{s}_{m}(t^{-})||\bm{s}_{e}(t^{-})]$

where $\Delta t$ is the time elapse between the previous time stamp $t^{-}$
and $t$, symbol —— means concatenating operation. After aggregating all the
messages from neighbors, the state vector of $\textbf{v}_{m}$ is updated as:

(2)
$\displaystyle\bm{s}_{m}(t)=mem(agg\\{\bm{\varrho}_{m}(t_{1}),...,\bm{\varrho}_{m}(t_{b})\\},\bm{s}_{m}(t^{-}))$

Here agg($\cdot$) is a general aggregation operation (support learnable
parameters). For the sake of simplicity, we only compute the mean of the most
recent messages in aggregation. We use mem($\cdot$) to represent a trainable
update function (e.g., GRU).

Build upon the updated state vector $\bm{s}_{m}(t)$ and $\bm{s}_{m}(e)$, a
time-aware node embedding can be generated at any time $t$ as following:

(3) $\displaystyle\bm{z}_{m}(t)=\sum_{j\in
n_{m}^{k}([0,t])}TGA(\bm{s}_{m}(t),\bm{s}_{e}(t),e_{m,e},\bm{v}_{m}(t),\bm{v}_{e}(t))$

TGA represents the temporal graph attention module
(DBLP:journals/corr/abs-2006-10637, ), where $L$ graph attention layers
compute $m$’s embedding by aggregating information from its L-hop temporal
neighbors.

Time-encoding: A finite dimensional mapping function is used to encode the
time elapse between $t$ and $t_{0}$ as the functional time encoding:
$\phi(t-t_{0})$. The time encoding function allows us to encode time elapse
with other graph features in an end-to-end manner. Specifically, the generic
time encoding function used in (rossi2020temporal, ) is employed, which is
invariant to time rescaling, and can capture the periodical patterns in the
data.

Embedding update: The TGA embedding module consists of a series of $L$ graph
attention layers, which aggregates information $\tilde{\bm{h}}_{i}^{(l)}(t)$
from each node’s $L$-hop temporal neighborhood:

(4) $\displaystyle\bm{h}_{m}^{(0)}(t)=$
$\displaystyle\bm{s}_{m}(t)+\bm{v}_{m}(t)$ (5) $\displaystyle\bm{z}_{m}(t)=$
$\displaystyle
MLP^{(l)}(\bm{h}_{m}^{(l-1)}(t)||\tilde{\bm{h}}_{m}^{(l)}(t))=\bm{h}_{m}^{(l)}(t)$

In each attention layer, a multi-head-attention is performed where a node
attends to its neighboring nodes, generating key, query and values based on
neighboring nodes’ representation and the encoded time elapses. After temporal
graph attention, an MLP (showing above) is used to integrate the reference
node representation with the aggregated information:

(6) $\displaystyle\tilde{\bm{h}}_{m}^{(l)}(t)=$ $\displaystyle
MultiHeadAttention^{(l)}(\bm{q}^{(l)}(t),\bm{K}^{(l)}(t),\bm{V}^{(l)}(t))$ (7)
$\displaystyle\bm{q}^{(l)}(t)=$ $\displaystyle\bm{h}_{i}^{(l-1)}(t)||\phi(0)$
(8) $\displaystyle\bm{K}^{(l)}(t)=$ $\displaystyle\bm{V}^{(l)}(t))$ (9)
$\displaystyle=$
$\displaystyle[\bm{h}_{e_{1}}^{(l-1)}(t)||\epsilon_{m,e_{1}}(t_{1})||\phi(t-t_{1}),...,$
(10)
$\displaystyle\bm{h}_{{e}_{N}}^{(l-1)}(t)f||\epsilon_{m,e_{N}}(t_{N})||\phi(t-t_{N})]$

###  4.4. Optimization and inference

The aforementioned model is trained in a self-supervised fashion. As our goal
is to predict the events that might happen in the time sequences in the next
time step, which corresponds to the edges linking the time sequence nodes and
event nodes. Therefore, we train our model with the edge prediction task. We
use the cross entropy to capture the prediction loss, which is the standard
loss function in edge prediction.

Event Forecasting Score: To convert the predicted event edges $\hat{A}$ into
an anomaly score, for each time series node $\textbf{v}_{m,\tau}^{t}$ with
window size $\tau$, the event $\textbf{v}_{e,\tau}^{t}$ that has the highest
probability connect to $\textbf{v}_{m,\tau}^{t}$ is retrieved and its pattern
in the original signal space is denoted as $\textbf{s}^{t}_{e,\tau}$. We
project the event $\textbf{s}^{t}_{e,\tau}$’s pattern back to its original
signal space, we then compute anomaly score based on the dynamic time wrapping
distance as follows:

(11)
$\displaystyle\omega_{1,m}^{t}=DTW(\textbf{X}^{t}_{m,\tau},\textbf{s}^{t}_{e,\tau})$

Residual Score: For a positive residual event $\textbf{v}_{e^{+}}^{t}$ at time
stamp $t$ where the forecasted results is not a positive residual
$\hat{\textbf{v}}_{e^{-}}^{t}$, we calculate a changing point score to
quantify the surprisal level as follows:

(12)
$\displaystyle\omega_{2,m}^{t}=\psi_{NLG}(||\textbf{X}^{t}_{m,\tau}-\textbf{X}^{t-\tau}_{m,\tau}||)$

where $\psi_{NLG}$ is a standard function that maps a scalar into the negative
log likelihood, which indicates the sparsity of this changing point in the
training data. A frequent changing signal may results in small
$\omega_{2,m}^{t}$ after the mapping. The function is learned in a data-driven
manner based on the training data of time series $m$.

The final anomaly score at time stamp $t$ is calculated as one of the two
following equations:

(13)
$\displaystyle\omega^{t}_{max}=\max_{m}(\omega_{1,m}^{t}\cdot\omega_{2,m}^{t})$
(14)
$\displaystyle\omega^{t}_{sum}=\sum_{m}(\omega_{1,m}^{t}\cdot\omega_{2,m}^{t})$

The selection of the above equations is driven by the property of anomaly in a
testing dataset. An empirical study shows that if very few (even a single)
time-series determine the anomaly labels of the multivariate time series $t$,
Eq. (14) is preferred. Otherwise, if the anomaly score of time $t$ is
determined by multiple time series, in majority of cases, Eq. (15) would be a
better choice. We evaluated both settings with comprehensive experiments.

##  5\. Experiment

We performed experiments on three datasets to demonstrate the effectiveness of
our model on multivariate anomaly detection. We adopted three public datasets:
SMAP (Soil Moisture Active Passive satellite) (DBLP:conf/kdd/HundmanCLCS18, ),
MSL (Mars Science Laboratory rover) (DBLP:conf/kdd/HundmanCLCS18, ), and SMD
(Server Machine Dataset) (DBLP:conf/kdd/SuZNLSP19, ). To compare with other
SOTA models, we use the same set of metrics to evaluate our model, including
precision, recall, and F1-score.

We follow the evaluation protocol of (su2019robust, ) which assumes anomalies
occur continuously to form contiguous anomaly segments and an alert is
classified as a correct catch if it is triggered within any subset of a ground
truth anomaly segment.

Table 2. Dataset statistics and hyper-parameters

| SMAP | MSL | SMD
---|---|---|---
Number of Time series nodes | 25 | 55 | 38
Number of event nodes | 130 | 70 | 64
Time window length ($\tau$) | 20 | 20 | 50
Time window stride ($\beta$) | 5 | 5 | 10
|  |  |

###  5.1. Settings

The proposed bipartite event graph contains $D+K$ number of nodes, where $D$
is the number of time series and $K$ is the number of events. The detailed
number of nodes is shown in Table 2. The length ($\tau$) and stride ($\beta$)
of the time window of each dataset are also displayed in Table 2. Our model
can be easily extended to a multi-scale version by incorporating more event
nodes with varieties of window sizes and strides. In this work, to keep things
simple, we demonstrated our method on a single time-scale. For each time
series, we set the maximum number of motifs detected to be three, and the
minimum cluster size of H-DBSCAN to be three. For the temporal attention
model, we set the number of multi-head to be 2. The GRU is employed to model
the time encoding. The dimension of node embedding and message vector are both
set to 64 respectively. Each model is trained after 10 epochs with the
learning rate 0.0001. The node features of both time-series and event nodes
are randomly initialized, and edge features are the one-hot embedding of the
numerical ID of the nodes on both sides of the edge.

Figure 5. The multivariate time-series anomaly on each of the three datasets:
SMAP, MSL, and SMD. The training partition is highlighted in blue.

###  5.2. Compare with state-of-the-art

We compared our solution with multiple baselines on SMAP, SMD, and MSL
datasets: DAGMM (DBLP:conf/iclr/ZongSMCLCC18, ) \- an autoencoder-based
anomaly detection model without taking into account of temporal information;
LSTM-VAE (8279425, ), LSTM-NDT (DBLP:conf/kdd/HundmanCLCS18, ), two state-of-
the-art LSTM-based anomaly detection solutions; and the most recent stochastic
VAE-based approaches (e.g., Omni-Anomaly(su2019robust, )) and graph attention-
based method such as MTAD-GAT(zhao2020multivariate, ). We selected these
baselines mainly because: (i) they are self-supervised algorithms that do not
need any training labels (different from (DBLP:conf/wsdm/Hu0CYR21, )), (ii)
they rely on a single scale of time-window (instead of multi-scale
(DBLP:conf/nips/ShenLK20, )) so that the performances are directly comparable.

The results are reported in Table 3. From the results we observed that the
proposed Event2Graph achieves new state-of-the-art performance on SMAP
dataset, ranked 2nd on SMD dataset, and perform very competitive on MSL
dataset. In SMAP dataset, we observe that some of the data suffers from
significant regime change during both training and testing, and our dynamic
graph-based solution helps the algorithm adapt to the regime change faster
than MTAD-GAT(zhao2020multivariate, ) and Omni-Anomaly(su2019robust, ). We
also observe that our algorithm significantly outperforms simple LSTM-based
solution (LSTM-VAE and LSTM-NDT), which assumes a complete inter-dependency
graph. They did not model the dynamic inter-dependency among time series,
while our node-level model explicitly encode the temporal information along
with the attention, which helps to reduce the false alarms in anomaly
detection. Furthermore, since DAGMM assumes a completely static relationship
between time series, the algorithm lacks of the capability to adapt to any
temporal evolving pattern.

Table 3. Compare with current state-of-the-art approaches. The top three
performance algorithms are highlighted.

Dataset | Method | Precision | Recall | F1
---|---|---|---|---
| DAGMM | 59.51 | 88.82 | 70.94
| LSTM-VAE | 79.22 | 70.75 | 78.42
| LSTM-NDT | 56.84 | 64.38 | 60.37
SMD | OmniAnomaly | 83.34 | 94.49 | 88.57
| MTAD-GAT | - | - | -
| Event2Graph (max) | 85.35 | 83.71 | 83.47
| Event2Graph (sum) | 88.61 | 83.38 | 84.93
| DAGMM | 54.12 | 99.34 | 70.07
| LSTM-VAE | 52.57 | 95.46 | 67.80
| LSTM-NDT | 59.34 | 53.74 | 56.40
MSL | OmniAnomaly | 88.67 | 91.17 | 89.89
| MTAD-GAT | 87.54 | 94.40 | 90.84
| Event2Graph (max) | 81.21 | 89.69 | 85.24
| Event2Graph (sum) | 88.12 | 83.11 | 85.55
| DAGMM | 58.45 | 90.58 | 71.05
| LSTM-VAE | 85.51 | 63.66 | 72.98
| LSTM-NDT | 89.65 | 88.46 | 89.05
SMAP | OmniAnomaly | 74.16 | 97.76 | 84.34
| MTAD-GAT | 89.06 | 91.23 | 90.13
| Event2Graph (max) | 90.48 | 93.21 | 91.82
| Event2Graph (sum) | 86.54 | 94.33 | 90.27
|  |  |  |

We also observe that the summation-based aggregation (Eq. 14) generally
outperforms max-based aggregation (Eq. 13). However, on SMAP dataset, the
summation-based performs quite well. One reason is that the anomalies on SMAP
dataset are mainly caused by a few time series, so the summation may bring
extra noise if the algorithm takes into account all the time series together.
However, the anomaly in SMD and MSL dataset may occur in a concurrent manner
where the summation operation boosted the true positives that concurrently
happened, while suppress false-alarms that happen individually on a specific
time series. From now on, the default setting of Event2Graph would be with
summation-base aggregation.

###  5.3. Quantitative Ablation study

The objective of this experiment is to provide detailed analysis over the
effectiveness of each proposed module. The experiments are conducted on SMD
dataset. By removing each of the critical module of our model, the
performances are reported in Table 4.

Table 4. Ablation study on SMD dataset.

Model | Precision | Recall | F1
---|---|---|---
Event2Graph | 88.61 | 83.38 | 84.93
w/o TGA | 84.81 | 82.94 | 82.51
w/o score $\omega_{1}$ | 82.83 | 80.29 | 79.45
w/o score $\omega_{2}$ | 81.51 | 71.76 | 73.81
w/o TGA and w/o $\omega_{1}$ | 80.56 | 81.09 | 79.47
w/o TGA and w/o $\omega_{2}$ | 75.84 | 64.33 | 64.51
$\psi_{NLG}$ only | 85.89 | 75.22 | 78.48
|  |  |

####  5.3.1. Effectiveness of temporal graph attention

Through replacing the temporal graph attention with a simple MLP module, the
model’s F1 score is reduced by 2.42%. It demonstrated that the temporal
attention plays an important role in aggregating information from neighborhood
nodes.

####  5.3.2. Effectiveness of event forecasting score $\omega_{1}$:

We remove the event matching score (Eq. 11) ) and just use a residual score
(Eq. 11) for anomaly detection. We observe that the model performance drop by
5.48%. The experiment indicates the event forecasting module in Event2Graph is
essential for accurate anomaly detection. It helps the system understand the
normality of time series, model the context of changing point, and work
jointly with the residual score to achieve high model performance. From the
experiment, we also observe that TGA plays a critical role in guarantee the
quality of the forecasting score $\omega_{1}$. By removing TGA, the
forecasting-based model performance (w/o score $\omega_{2}$) degrade from
73.81% F1 score to 64.51%.

####  5.3.3. Effectiveness of residual score $\omega_{2}$

We remove the residual score (Eq. 12) and only use forecasting score to detect
anomaly. The model suffers more than 10% performance degradation. This results
shows that modeling the residual score is essential to allow the system to
identify the anomaly regions that may not well characterized by the event
matching scores. Since this score mainly captures the abrupt change of values
in time-series, we can infer that the anomaly pattern in SMD dataset are
mainly due to the changing points. We also observe that the TGA may not
contribute much to residual event modeling, a $\omega_{2}$ only model perform
competitively well with or without TGA. This phenomenon tells us the residual
event is something unexpected and an accurate changing point detection is
better than learning simple, repeatable patterns.

####  5.3.4. Effectiveness of log likelihood function $\psi_{NLG}$

To demonstrate our assumption, we compare the results with a baseline that
using a naive changing point detection score $\psi_{NLG}$, without adopting
any forecasting or event modeling model. Surprisingly, the changing-point only
model is able to achieve 78.48% F1-score. Outperforms DAGMM, LSTM-VAE, as well
as LSTM-NDT. This results further confirm that it is hard to learn any static
pattern by using LSTM or autoencoder for anomaly detection on a challenge
dataset, a dynamic model that is able to adapt to the changing of time-series
patterns is preferred.

###  5.4. Qualitative Ablation study

Figure 6. Visualization of the event forecasting score $\omega_{1}$ and
residual scores $\omega_{2}$ on a randomly selected machine in the SMD
dataset. The top 10 time-series of this machine is visualized via (a), and the
ground-truth anomaly in all sub-figures are highlighted with red. Subfigure
(b,c) corresponding to the estimation of $\omega_{1}$ (green line) and
$\omega_{2}$ (blue line) under different settings of Table 4. In (d), the
change point score generated by $\psi_{NLG}$ is marked in blue, green line
corresponding to the minimal residual error by matching all possible event
nodes.

To better understand the behavior of our Event2Graph approach, we compared the
propose approach with other variations via visualization. We present the two
scores ($\omega_{1}$ and $\omega_{2}$) generated by the system on a randomly
selected machine from SMD dataset. The original time series of this machine is
shown in Fig. 6(a). The $\omega_{1}$ and $\omega_{2}$ are visualized in Fig.
6(b-d) under different experiment protocol. From the experiment, we observe
that the temporal graph attention plays a very important role in terms of
stabilizing the system. By removing this module, the whole system seems to
overreact to small noises in the time series. From Fig. 6(d) we can also have
a better understanding of the system’s behavior without modeling the event
transition patterns with Event2Graph. Both scores become very similar across
times. This phenomenon indicates the system have no understanding about the
current state of the time series and cannot make a reasonable guess on what
will happen next. In this case, the anomaly are purely reflected by the value
change, instead of jointly considering the context of the value change like
Event2Graph.

Finally, we wish to highlight that our Event2Graph system is able to capture
the anomalies that have not been well marked by the ground-truth labels. For
example in Fig. 7, we visualize the anomaly scores generated by Event2Graph on
a testing data in the SMD dataset. From the outputs of our system, it is
clearly observed that there is an potential anomaly region after the ground-
truth anomaly (red peak in Fig. 7(a)). Our pattern matching score $\omega_{1}$
is significantly increased near the end of testing period, and also the
residual error scores $\omega_{2}$ are frequently spiked. This phenomenon
aligned well with human eye’s inspection where there is an increasing, noisy
trend emerging on three out of ten time series close to the end of testing
period. We wish future works are able to re-label this dataset so that the
ground-truth is able to reflect all types of the anomalies that corresponding
to human’s intuition.

Figure 7. Case study of a potential false negative. Testing region is
highlighted in the left, and outputs of Event2Graph: $\omega_{1}$ (green line)
and $\omega_{2}$ (blue line) are visualized on the right.

##  6\. Conclusion

We proposed an event-driven bipartite graph solution for multivariate time
series anomaly detection. The solution does not assume any inter-dependency on
time series, and all the relations are learned in a dynamic, data-driven
manner. Our design is based on edge-stream so no adjacency matrix of the graph
is required as input. As the system’s memory is defined on the node level, our
design left plenty space for future extensions such as inductive learning and
parallel computation. Our solution achieved very competitive results on three
anomaly detection datasets, and we encourage future works to explore further
using bipartite event graph for multivariate anomaly detection.

## References

  * (1) Berndt, D. J., and Clifford, J. Using dynamic time warping to find patterns in time series.  In Knowledge Discovery in Databases: Papers from the 1994 AAAI Workshop, Seattle, Washington, USA, July 1994. Technical Report WS-94-03 (1994), U. M. Fayyad and R. Uthurusamy, Eds., AAAI Press, pp. 359–370.
  * (2) Cheng, Z., Yang, Y., Wang, W., Hu, W., Zhuang, Y., and Song, G. Time2graph: Revisiting time series modeling with dynamic shapelets.  In The Thirty-Fourth AAAI Conference on Artificial Intelligence, AAAI 2020, The Thirty-Second Innovative Applications of Artificial Intelligence Conference, IAAI 2020, The Tenth AAAI Symposium on Educational Advances in Artificial Intelligence, EAAI 2020, New York, NY, USA, February 7-12, 2020 (2020), AAAI Press, pp. 3617–3624.
  * (3) Chung, J., Gülçehre, Ç., Cho, K., and Bengio, Y. Empirical evaluation of gated recurrent neural networks on sequence modeling.  CoRR abs/1412.3555 (2014).
  * (4) Deng, A., and Hooi, B. Graph neural network-based anomaly detection in multivariate time series.  In Proceedings of the 35th AAAI Conference on Artificial Intelligence, Vancouver, BC, Canada (2021), pp. 2–9.
  * (5) Hochreiter, S., and Schmidhuber, J. Long short-term memory.  Neural Comput. 9, 8 (1997), 1735–1780.
  * (6) Hu, W., Yang, Y., Cheng, Z., Yang, C., and Ren, X. Time-series event prediction with evolutionary state graph.  In WSDM ’21, The Fourteenth ACM International Conference on Web Search and Data Mining, Virtual Event, Israel, March 8-12, 2021 (2021), L. Lewin-Eytan, D. Carmel, E. Yom-Tov, E. Agichtein, and E. Gabrilovich, Eds., ACM, pp. 580–588.
  * (7) Hundman, K., Constantinou, V., Laporte, C., Colwell, I., and Söderström, T. Detecting spacecraft anomalies using lstms and nonparametric dynamic thresholding.  In Proceedings of the 24th ACM SIGKDD International Conference on Knowledge Discovery & Data Mining, KDD 2018, London, UK, August 19-23, 2018 (2018), Y. Guo and F. Farooq, Eds., ACM, pp. 387–395.
  * (8) Kumar, S., Zhang, X., and Leskovec, J. Predicting dynamic embedding trajectory in temporal interaction networks.  In Proceedings of the 25th ACM SIGKDD International Conference on Knowledge Discovery & Data Mining, KDD 2019, Anchorage, AK, USA, August 4-8, 2019 (2019), A. Teredesai, V. Kumar, Y. Li, R. Rosales, E. Terzi, and G. Karypis, Eds., ACM, pp. 1269–1278.
  * (9) Malhotra, P., Ramakrishnan, A., Anand, G., Vig, L., Agarwal, P., and Shroff, G. Lstm-based encoder-decoder for multi-sensor anomaly detection.  CoRR abs/1607.00148 (2016).
  * (10) McInnes, L., Healy, J., and Astels, S. hdbscan: Hierarchical density based clustering.  J. Open Source Softw. 2, 11 (2017), 205.
  * (11) Mueen, A., Zhu, Y., Yeh, M., Kamgar, K., Viswanathan, K., Gupta, C., and Keogh, E. The fastest similarity search algorithm for time series subsequences under euclidean distance.  http://www.cs.unm.edu/~mueen/FastestSimilaritySearch.html.
  * (12) Nguyen, G. H., Lee, J. B., Rossi, R. A., Ahmed, N. K., Koh, E., and Kim, S. Dynamic network embeddings: From random walks to temporal random walks.  In IEEE International Conference on Big Data, Big Data 2018, Seattle, WA, USA, December 10-13, 2018 (2018), N. Abe, H. Liu, C. Pu, X. Hu, N. K. Ahmed, M. Qiao, Y. Song, D. Kossmann, B. Liu, K. Lee, J. Tang, J. He, and J. S. Saltz, Eds., IEEE, pp. 1085–1092.
  * (13) Park, D., Hoshi, Y., and Kemp, C. C. A multimodal anomaly detector for robot-assisted feeding using an lstm-based variational autoencoder.  IEEE Robotics and Automation Letters 3, 3 (2018), 1544–1551.
  * (14) Rossi, E., Chamberlain, B., Frasca, F., Eynard, D., Monti, F., and Bronstein, M. Temporal graph networks for deep learning on dynamic graphs.  arXiv preprint arXiv:2006.10637 (2020).
  * (15) Rossi, E., Chamberlain, B., Frasca, F., Eynard, D., Monti, F., and Bronstein, M. M. Temporal graph networks for deep learning on dynamic graphs.  CoRR abs/2006.10637 (2020).
  * (16) Sankar, A., Wu, Y., Gou, L., Zhang, W., and Yang, H. Dysat: Deep neural representation learning on dynamic graphs via self-attention networks.  In WSDM ’20: The Thirteenth ACM International Conference on Web Search and Data Mining, Houston, TX, USA, February 3-7, 2020 (2020), J. Caverlee, X. B. Hu, M. Lalmas, and W. Wang, Eds., ACM, pp. 519–527.
  * (17) Shen, L., Li, Z., and Kwok, J. T. Timeseries anomaly detection using temporal hierarchical one-class network.  In Advances in Neural Information Processing Systems 33: Annual Conference on Neural Information Processing Systems 2020, NeurIPS 2020, December 6-12, 2020, virtual (2020), H. Larochelle, M. Ranzato, R. Hadsell, M. Balcan, and H. Lin, Eds.
  * (18) Siffer, A., Fouque, P.-A., Termier, A., and Largouet, C. Anomaly detection in streams with extreme value theory.  In Proceedings of the 23rd ACM SIGKDD International Conference on Knowledge Discovery and Data Mining (2017), pp. 1067–1075.
  * (19) Su, Y., Zhao, Y., Niu, C., Liu, R., Sun, W., and Pei, D. Robust anomaly detection for multivariate time series through stochastic recurrent neural network.  In Proceedings of the 25th ACM SIGKDD International Conference on Knowledge Discovery & Data Mining, KDD 2019, Anchorage, AK, USA, August 4-8, 2019 (2019), A. Teredesai, V. Kumar, Y. Li, R. Rosales, E. Terzi, and G. Karypis, Eds., ACM, pp. 2828–2837.
  * (20) Su, Y., Zhao, Y., Niu, C., Liu, R., Sun, W., and Pei, D. Robust anomaly detection for multivariate time series through stochastic recurrent neural network.  In Proceedings of the 25th ACM SIGKDD International Conference on Knowledge Discovery & Data Mining (2019), pp. 2828–2837.
  * (21) Trivedi, R., Farajtabar, M., Biswal, P., and Zha, H. Dyrep: Learning representations over dynamic graphs.  In 7th International Conference on Learning Representations, ICLR 2019, New Orleans, LA, USA, May 6-9, 2019 (2019), OpenReview.net.
  * (22) Velickovic, P., Cucurull, G., Casanova, A., Romero, A., Liò, P., and Bengio, Y. Graph attention networks.  In 6th International Conference on Learning Representations, ICLR 2018, Vancouver, BC, Canada, April 30 - May 3, 2018, Conference Track Proceedings (2018), OpenReview.net.
  * (23) Xu, D., Ruan, C., Körpeoglu, E., Kumar, S., and Achan, K. Inductive representation learning on temporal graphs.  In 8th International Conference on Learning Representations, ICLR 2020, Addis Ababa, Ethiopia, April 26-30, 2020 (2020), OpenReview.net.
  * (24) Yeh, C. M., Zhu, Y., Ulanova, L., Begum, N., Ding, Y., Dau, H. A., Silva, D. F., Mueen, A., and Keogh, E. J. Matrix profile I: all pairs similarity joins for time series: A unifying view that includes motifs, discords and shapelets.  In IEEE 16th International Conference on Data Mining, ICDM 2016, December 12-15, 2016, Barcelona, Spain (2016), F. Bonchi, J. Domingo-Ferrer, R. Baeza-Yates, Z. Zhou, and X. Wu, Eds., IEEE Computer Society, pp. 1317–1322.
  * (25) Yu, W., Cheng, W., Aggarwal, C. C., Chen, H., and Wang, W. Link prediction with spatial and temporal consistency in dynamic networks.  In Proceedings of the Twenty-Sixth International Joint Conference on Artificial Intelligence, IJCAI 2017, Melbourne, Australia, August 19-25, 2017 (2017), C. Sierra, Ed., ijcai.org, pp. 3343–3349.
  * (26) Yu, W., Cheng, W., Aggarwal, C. C., Zhang, K., Chen, H., and Wang, W. Netwalk: A flexible deep embedding approach for anomaly detection in dynamic networks.  In Proceedings of the 24th ACM SIGKDD International Conference on Knowledge Discovery & Data Mining, KDD 2018, London, UK, August 19-23, 2018 (2018), Y. Guo and F. Farooq, Eds., ACM, pp. 2672–2681.
  * (27) Zhang, Y., Ge, Z., Greenberg, A. G., and Roughan, M. Network anomography.  In Proceedings of the 5th Internet Measurement Conference, IMC 2005, Berkeley, California, USA, October 19-21, 2005 (2005), USENIX Association, pp. 317–330.
  * (28) Zhao, H., Wang, Y., Duan, J., Huang, C., Cao, D., Tong, Y., Xu, B., Bai, J., Tong, J., and Zhang, Q. Multivariate time-series anomaly detection via graph attention network.  In 20th IEEE International Conference on Data Mining, ICDM 2020, Sorrento, Italy, November 17-20, 2020 (2020), C. Plant, H. Wang, A. Cuzzocrea, C. Zaniolo, and X. Wu, Eds., IEEE, pp. 841–850.
  * (29) Zhu, Y., Yeh, C. M., Zimmerman, Z., Kamgar, K., and Keogh, E. J. Matrix profile XI: SCRIMP++: time series motif discovery at interactive speeds.  In IEEE International Conference on Data Mining, ICDM 2018, Singapore, November 17-20, 2018 (2018), IEEE Computer Society, pp. 837–846.
  * (30) Zong, B., Song, Q., Min, M. R., Cheng, W., Lumezanu, C., Cho, D., and Chen, H. Deep autoencoding gaussian mixture model for unsupervised anomaly detection.  In 6th International Conference on Learning Representations, ICLR 2018, Vancouver, BC, Canada, April 30 - May 3, 2018, Conference Track Proceedings (2018), OpenReview.net.


{\put(0.0,0.0){}}\ignorespaces{\color[rgb]{0.75,0.75,0.75}\definecolor[named]{pgfstrokecolor}{rgb}{0.75,0.75,0.75}\pgfsys@color@gray@stroke{0.75}\pgfsys@color@gray@fill{0.75}\put(0.0,0.0){}
\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\put(0.0,0.0){}}\ignorespaces{\color[rgb]{0.75,0.75,0.75}\definecolor[named]{pgfstrokecolor}{rgb}{0.75,0.75,0.75}\pgfsys@color@gray@stroke{0.75}\pgfsys@color@gray@fill{0.75}\put(0.0,0.0){}
\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\put(0.0,0.0){}}\ignorespaces\end{picture}
$\displaystyle=2\times(g)^{{\mbox{\tiny\rm{BB}}}}$
$\displaystyle=-2g_{{\mbox{\tiny\rm{3d}}}}^{2}C_{\mbox{\tiny\rm{A}}}^{2}C_{\mbox{\tiny\rm{F}}}\int_{{\bf
pk}}e^{ip_{z}L}\frac{1}{p^{2}k^{2}(p-k)^{2}}$
$\displaystyle\hphantom{=-2g_{{\mbox{\tiny\rm{3d}}}}^{2}C_{\mbox{\tiny\rm{A}}}^{2}C_{\mbox{\tiny\rm{F}}}}\times\frac{1}{2}\Bigl{[}3p^{2}+k^{2}-(p-k)^{2}+2(d-2)p_{z}(p_{z}+k_{z})\Bigr{]}$
$\displaystyle=-2g_{{\mbox{\tiny\rm{3d}}}}^{2}C_{\mbox{\tiny\rm{A}}}^{2}C_{\mbox{\tiny\rm{F}}}\frac{3(d-3)(d-1)}{(d-4)}\int_{{\bf
pk}}e^{ip_{z}L}\frac{1}{k^{2}(p-k)^{2}}\;,$ (103)

where the last line holds in dimensional regularization and vanishes for $d=3$
since the integrand is finite. It was derived using integration-by-parts
identities at two-loop level similar to eq. (67). The relevant relation is

$\displaystyle\int_{{\bf p}_{\perp}{\bf
k}}\frac{p_{z}(p_{z}+k_{z})}{p^{2}k^{2}(p-k)^{2}}$
$\displaystyle=\frac{3}{2}\frac{(2d-5)}{(d-4)}\int_{{\bf p}_{\perp}{\bf
k}}\frac{1}{k^{2}(p-k)^{2}}\;.$ (104)

##  Appendix B Simulation details

In order to carry out non-perturbative computations, we discretize the
continuum action of EQCD (14) on a three-dimensional numerical grid with
lattice spacing $a$,

$\displaystyle S_{\mathrm{EQCD},L}=$
$\displaystyle\;\beta\sum_{x,i>j}\Bigl{(}1-\frac{1}{3}\Box_{x,ij}\Bigr{)}$
$\displaystyle+2\sum_{x,i}\mathrm{Tr}\,\Bigl{(}\Phi_{\mbox{\tiny\rm{L}}}^{2}(x)-\Phi_{\mbox{\tiny\rm{L}}}(x)U_{i}(x)\Phi_{\mbox{\tiny\rm{L}}}(x+a\hat{i})U^{\dagger}_{i}(x)\Bigr{)}$
$\displaystyle+\sum_{x}\left[Z_{4}(x+\delta
x)\left(\mathrm{Tr}\,\Phi_{\mbox{\tiny\rm{L}}}^{2}(x)\right)^{2}+Z_{2}(y+\delta
y)\mathrm{Tr}\,\Phi^{2}_{\mbox{\tiny\rm{L}}}(x)\right]$ (105)
$\displaystyle\Box_{x,ij}\equiv$
$\displaystyle\;U_{i}(x)U_{j}(x+a\hat{i})U^{\dagger}_{i}(x+a\hat{j})U^{\dagger}_{j}(x)\;,$
(106)

where $U_{i}(x)$ is the gauge link in spatial direction $i$, connecting
lattice sites $x$ and $x+a\hat{i}$. We rescaled the adjoint scalar field to
its lattice version $\Phi_{\mbox{\tiny\rm{L}}}$ such that its wave-function
normalization $Z_{\Phi}=1$ is always enforced. The subscript L will be dropped
in the context of lattice calculations. Analytical calculations in EQCD
lattice perturbation theory yield expressions for the counterterm of the
inverse gauge coupling $\beta$, the quartic coupling $\delta x$, and the
multiplicative mass ($Z_{2}$) and quartic ($Z_{4}$) renormalization that
compensate for discretization errors up to $\mathcal{O}(a)$ Moore:1997np . A
semi-analytical computation of $\delta y$ completes the
$\mathcal{O}(a)$-improvement at the Lagrangian level Moore:2019lua . At the
operator level, we use the lattice implementation of the modified Wilson line
(16) described in Panero:2013pla , with the $\mathcal{O}(a)$-improvement
delivered in DOnofrio:2014mld , to connect the respective pairs of operators.
The ‘color-electric’ field operator is discretized using gauge-invariant
central derivative (see Moore:2020wvy ), whereas the ‘color-magnetic’ field
operator is calculated using the clover discretization Weisz:1982zw ;
Weisz:1983bn . For the $\langle BB\rangle$-correlator, only the single-
plaquette correlations of the pair of clover operators are used, as elaborated
in sec. 5. Consequently, the only remaining sources of errors at
$\mathcal{O}(a)$ are of the form (47) and can be determined in an overall fit,
as outlined in sec. 5.

Our numerical implementation is based on openQCD-1.6 by Martin Lüscher openQCD
. Using a combined update of one heatbath sweep succeeded by four over-
relaxation sweeps through the volume, we update the lattice sites in a
checkerboard ordering. We use the multi-level algorithm proposed by Lüscher
and Weisz Luscher:2001up to reduce the noise in relation to the signal for
non-local operators in lattice gauge theories. In particular, we divide our
volume in four sub-volumes along the $z$-axis and freeze the surfaces between
them for 80 combined heatbath/over-relaxation sweeps, before we allow a single
combined sweep through the entire volume. For further details, see
Moore:2019lua ; Moore:2019lgw ; Moore:2020wvy . We give the raw data at finite
lattice spacings $a$ in an online repository rawData .

We provide the parameters of our simulations in tab. 4. As argued in
Hietanen:2008tv , EQCD possesses a mass gap and therefore finite volume
effects are exponentially suppressed. Thus, one can keep finite volume effects
under control by maintaining a sufficiently large volume along the rules of
thumb given in Hietanen:2008tv . We measured all three correlators on the same
lattice configurations introducing correlations among different
$g_{{\mbox{\tiny\rm{3d}}}}^{2}L$ and correlators that had to be accounted for
in the jackknife error analysis by keeping the binned jackknife data until the
numerical integration and calculating the error thereafter. Equal
computational resources were spent on all four scenarios
$(T,N_{\mbox{\scriptsize f}})$; the difference in statistics can be explained
by a slight increase in the acceptance rate for smaller $x$, since the quartic
self-coupling of the scalars is taken into account via a Metropolis step in
the scalar heatbath-update.

$g_{{\mbox{\tiny\rm{3d}}}}^{2}a$ | $x_{\mathrm{cont}}$ | $y_{\mathrm{cont}}$ | $N_{\mbox{\tiny\rm{x}}}N_{\mbox{\tiny\rm{y}}}N_{\mbox{\tiny\rm{z}}}$ | $L/a$ | statistics
---|---|---|---|---|---
$1/4$ | $0.08896$ | $0.452423$ | $24^{2}\times 48$ | $4,6,8,10,12$ | $14260$
$1/6$ | $0.08896$ | $0.452423$ | $36^{2}\times 72$ | $6,12,18$ | $8440$
$1/8$ | $0.08896$ | $0.452423$ | $48^{2}\times 96$ | $4,6,8,12,16,20,24$ | $8709$
$1/12$ | $0.08896$ | $0.452423$ | $72^{2}\times 144$ | $6,12,18,24,30,36$ | $5720$
$1/16$ | $0.08896$ | $0.452423$ | $96^{2}\times 192$ | $4,8,12,16,24,32,40$ | $6060$
$1/24$ | $0.08896$ | $0.452423$ | $144^{2}\times 288$ | $6,12,18,24$ | $780$
$1/32$ | $0.08896$ | $0.452423$ | $192^{2}\times 384$ | $4,8,16,24,32$ | $440$
$1/4$ | $0.0677528$ | $0.586204$ | $24^{2}\times 48$ | $4,6,8,10,12$ | $17940$
$1/6$ | $0.0677528$ | $0.586204$ | $36^{2}\times 72$ | $6,12,18$ | $8700$
$1/8$ | $0.0677528$ | $0.586204$ | $48^{2}\times 96$ | $4,6,8,12,16,20,24$ | $8620$
$1/12$ | $0.0677528$ | $0.586204$ | $72^{2}\times 144$ | $6,12,18,24,30,36$ | $6820$
$1/16$ | $0.0677528$ | $0.586204$ | $96^{2}\times 192$ | $4,8,12,16,24,32,40$ | $6060$
$1/24$ | $0.0677528$ | $0.586204$ | $144^{2}\times 288$ | $6,12,18,24$ | $1400$
$1/32$ | $0.0677528$ | $0.586204$ | $192^{2}\times 384$ | $4,8,16,24,32$ | $820$
$1/4$ | $0.0463597$ | $0.823449$ | $24^{2}\times 48$ | $4,6,8,10,12$ | $18000$
$1/6$ | $0.0463597$ | $0.823449$ | $36^{2}\times 72$ | $6,12,18$ | $10680$
$1/8$ | $0.0463597$ | $0.823449$ | $48^{2}\times 96$ | $4,6,8,12,16,20,24$ | $8620$
$1/12$ | $0.0463597$ | $0.823449$ | $72^{2}\times 144$ | $6,12,18,24,30,36$ | $6120$
$1/16$ | $0.0463597$ | $0.823449$ | $96^{2}\times 192$ | $4,8,12,16,24,32,40$ | $5160$
$1/24$ | $0.0463597$ | $0.823449$ | $144^{2}\times 288$ | $6,12,18,24$ | $1420$
$1/32$ | $0.0463597$ | $0.823449$ | $192^{2}\times 384$ | $4,8,16,24,32$ | $640$
$1/4$ | $0.0178626$ | $1.64668$ | $24^{2}\times 48$ | $4,6,8,10,12$ | $19380$
$1/6$ | $0.0178626$ | $1.64668$ | $36^{2}\times 72$ | $6,12,18$ | $10680$
$1/8$ | $0.0178626$ | $1.64668$ | $48^{2}\times 96$ | $4,6,8,12,16,20,24$ | $8600$
$1/12$ | $0.0178626$ | $1.64668$ | $72^{2}\times 144$ | $6,12,18,24,30,36$ | $6500$
$1/16$ | $0.0178626$ | $1.64668$ | $96^{2}\times 192$ | $4,8,12,16,24,32,40$ | $3400$
$1/24$ | $0.0178626$ | $1.64668$ | $144^{2}\times 288$ | $6,12,18,24$ | $1300$
$1/32$ | $0.0178626$ | $1.64668$ | $192^{2}\times 384$ | $4,8,16,24,32$ | $680$
Table 4: Parameters for all EQCD multi-level simulations. We give the
continuum expressions for the rescaled quartic scalar self-coupling
$x_{\mbox{\tiny\rm{cont}}}$ and the rescaled screening mass
$y_{\mbox{\tiny\rm{cont}}}$. A conversion into lattice units is done using the
respective lattice-continuum relations Moore:1997np ; Moore:2019lua .

## References

  * (1) M. Connors, C. Nattrass, R. Reed and S. Salur, _Jet measurements in heavy ion physics_ , _Rev. Mod. Phys._ 90 (2018) 025005 [1705.01974].
  * (2) L. Cunqueiro and A.M. Sickles, _Studying the QGP with Jets at the LHC and RHIC_ , 2110.14490.
  * (3) V. Klimov, _Spectrum of Elementary Fermi Excitations in Quark Gluon Plasma. (In Russian)_ , _Sov. J. Nucl. Phys._ 33 (1981) 934\.
  * (4) V. Klimov, _Collective Excitations in a Hot Quark Gluon Plasma_ , _Sov. Phys. JETP_ 55 (1982) 199.
  * (5) H. Weldon, _Covariant Calculations at Finite Temperature: The Relativistic Plasma_ , _Phys. Rev. D_ 26 (1982) 1394.
  * (6) H. Weldon, _Effective Fermion Masses of Order gT in High Temperature Gauge Theories with Exact Chiral Invariance_ , _Phys. Rev. D_ 26 (1982) 2789.
  * (7) E. Braaten and R.D. Pisarski, _Simple effective Lagrangian for hard thermal loops_ , _Phys. Rev. D_ 45 (1992) 1827.
  * (8) S. Caron-Huot, _On supersymmetry at finite temperature_ , _Phys. Rev. D_ 79 (2009) 125002 [0808.0155].
  * (9) A.D. Linde, _Infrared Problem in Thermodynamics of the Yang-Mills Gas_ , _Phys. Lett._ 96B (1980) 289.
  * (10) S. Caron-Huot, _O(g) plasma effects in jet quenching_ , _Phys. Rev._ D79 (2009) 065039 [0811.1603].
  * (11) E. Braaten, _Solution to the perturbative infrared catastrophe of hot gauge theories_ , _Phys.Rev.Lett._ 74 (1995) 2164 [hep-ph/9409434].
  * (12) E. Braaten and A. Nieto, _Effective field theory approach to high temperature thermodynamics_ , _Phys. Rev._ D51 (1995) 6990 [hep-ph/9501375].
  * (13) E. Braaten and A. Nieto, _Free energy of QCD at high temperature_ , _Phys.Rev._ D53 (1996) 3421 [hep-ph/9510408].
  * (14) K. Kajantie, M. Laine, K. Rummukainen and M.E. Shaposhnikov, _Generic rules for high temperature dimensional reduction and their application to the standard model_ , _Nucl.Phys._ B458 (1996) 90 [hep-ph/9508379].
  * (15) K. Kajantie, M. Laine, K. Rummukainen and M.E. Shaposhnikov, _3-D SU(N) + adjoint Higgs theory and finite temperature QCD_ , _Nucl. Phys._ B503 (1997) 357 [hep-ph/9704416].
  * (16) G.D. Moore and N. Schlusser, _The nonperturbative contribution to asymptotic masses_ , _Phys. Rev. D_ 102 (2020) 094512 [2009.06614].
  * (17) T. Appelquist and R.D. Pisarski, _High-Temperature Yang-Mills Theories and Three-Dimensional Quantum Chromodynamics_ , _Phys. Rev._ D23 (1981) 2305.
  * (18) S. Nadkarni, _Dimensional Reduction in Hot QCD_ , _Phys. Rev._ D27 (1983) 917.
  * (19) M. Laine and Y. Schröder, _Two-loop QCD gauge coupling at high temperatures_ , _JHEP_ 03 (2005) 067 [hep-ph/0503061].
  * (20) I. Ghişoiu, J. Möller and Y. Schröder, _Debye screening mass of hot yang-mills theory to three-loop order_ , _JHEP_ 2015 (2015) 121 [1509.08727].
  * (21) M. Laine, P. Schicho and Y. Schröder, _Soft thermal contributions to 3-loop gauge coupling_ , _JHEP_ 05 (2018) 037 [1803.08689].
  * (22) M. Laine, P. Schicho and Y. Schröder, _A QCD Debye mass in a broad temperature range_ , _Phys. Rev. D_ 101 (2020) 023532 [1911.09123].
  * (23) ALPHA collaboration, _QCD Coupling from a Nonperturbative Determination of the Three-Flavor $\Lambda$ Parameter_, _Phys. Rev. Lett._ 119 (2017) 102001 [1706.03821].
  * (24) A. Vuorinen and L.G. Yaffe, _Z(3)-symmetric effective theory for SU(3) Yang-Mills theory at high temperature_ , _Phys. Rev. D_ 74 (2006) 025011 [hep-ph/0604100].
  * (25) J. Casalderrey-Solana and D. Teaney, _Transverse Momentum Broadening of a Fast Quark in a N=4 Yang Mills Plasma_ , _JHEP_ 04 (2007) 039 [hep-th/0701123].
  * (26) F. D’Eramo, H. Liu and K. Rajagopal, _Transverse Momentum Broadening and the Jet Quenching Parameter, Redux_ , _Phys. Rev. D_ 84 (2011) 065015 [1006.1367].
  * (27) M. Benzke, N. Brambilla, M.A. Escobedo and A. Vairo, _Gauge invariant definition of the jet quenching parameter_ , _JHEP_ 02 (2013) 129 [1208.4253].
  * (28) P. Aurenche, F. Gelis, G.D. Moore and H. Zaraket, _Landau-Pomeranchuk-Migdal resummation for dilepton production_ , _JHEP_ 12 (2002) 006 [hep-ph/0211036].
  * (29) P. Giovannangeli, _Two loop renormalization of the magnetic coupling and non-perturbative sector in hot QCD_ , _Nucl. Phys. B_ 738 (2006) 23 [hep-ph/0506318].
  * (30) M. Laine and O. Philipsen, _Gauge invariant scalar and field strength correlators in 3-D_ , _Nucl. Phys. B_ 523 (1998) 267 [hep-lat/9711022].
  * (31) J. Ghiglieri, P. Schicho, N. Schlusser and E. Weitz _, in preparation_ (2023) .
  * (32) S. Nadkarni, _Nonabelian Debye Screening. 1. The Color Averaged Potential_ , _Phys. Rev. D_ 33 (1986) 3738.
  * (33) Y. Burnier, M. Laine and M. Vepsäläinen, _Dimensionally regularized Polyakov loop correlators in hot QCD_ , _JHEP_ 01 (2010) 054 [0911.3480].
  * (34) G.D. Moore, _O(a) errors in 3-D SU(N) Higgs theories_ , _Nucl. Phys._ B523 (1998) 569 [hep-lat/9709053].
  * (35) G.D. Moore and N. Schlusser, _Full O(a) improvement in electrostatic QCD_ , _Phys. Rev._ D100 (2019) 034510 [1905.09708].
  * (36) M. Panero, K. Rummukainen and A. Schäfer, _Lattice Study of the Jet Quenching Parameter_ , _Phys. Rev. Lett._ 112 (2014) 162001 [1307.5850].
  * (37) M. D’Onofrio, A. Kurkela and G.D. Moore, _Renormalization of Null Wilson Lines in EQCD_ , _JHEP_ 03 (2014) 125 [1401.7951].
  * (38) P. Weisz, _Continuum Limit Improved Lattice Action for Pure Yang-Mills Theory. 1._ , _Nucl. Phys. B_ 212 (1983) 1.
  * (39) P. Weisz and R. Wohlert, _Continuum Limit Improved Lattice Action for Pure Yang-Mills Theory. 2._ , _Nucl. Phys. B_ 236 (1984) 397.
  * (40) http://luscher.web.cern.ch/luscher/openQCD/index.html.
  * (41) M. Lüscher and P. Weisz, _Locality and exponential error reduction in numerical lattice gauge theory_ , _JHEP_ 09 (2001) 010 [hep-lat/0108014].
  * (42) G.D. Moore and N. Schlusser, _Transverse momentum broadening from the lattice_ , _Phys. Rev. D_ 101 (2020) 014505 [1911.13127].
  * (43) J. Ghiglieri, G.D. Moore, P. Schicho and N. Schlusser, _ff-eqcd: v1.0_ , Jan., 2022.  10.5281/zenodo.5884631.
  * (44) A. Hietanen, K. Kajantie, M. Laine, K. Rummukainen and Y. Schröder, _Three-dimensional physics and the pressure of hot QCD_ , _Phys. Rev._ D79 (2009) 045018 [0811.4664].


# Isometry theorem of gradient Shrinking Ricci solitons

Absos Ali Shaikh1 and Chandan Kumar Mondal2
12Department of Mathematics,
University of Burdwan, Golapbag,
Burdwan-713104,
West Bengal, India<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract.

In this paper, we have proved that if a complete conformally flat gradient
shrinking Ricci soliton has linear volume growth or the scalar curvature is
finitely integrable and also the reciprocal of the potential function is
subharmonic, then the manifold is isometric to the Euclidean sphere. As a
consequence, we have showed that a four dimensional gradient shrinking Ricci
soliton satisfying some conditions is isometric to $\mathbb{S}^{4}$ or
$\mathbb{RP}^{4}$ or $\mathbb{CP}^{2}$. We have also deduced a condition for
the shrinking Ricci soliton to be compact with quadratic volume growth.

††footnotetext: $\mathbf{2020}$ Mathematics Subject Classification: 53C20;
53C21
Key words and phrases: Shrinking Ricci soliton, fundamental group, simply
connected, volume growth, Riemannian manifold.

##  1\. Introduction and results

A complete Riemannian manifold $(M,g)$, of dimension $n\geq 2$, is called a
Ricci soliton if there exists a vector field $X$ such that the following
condition holds:

(1) $Ric+\frac{1}{2}\pounds_{X}g=\lambda g,$

where $\lambda$ is a constant and $\pounds$ denotes the Lie derivative. The
vector field $X$ is called potential vector field. The Ricci solitons are
self-similar solutions to the Ricci flow, which is developed by Hamilton [13,
15]. If $X=\nabla f$, for some function $f\in C^{\infty}(M)$, then (1) takes
the form

(2) $\nabla f^{2}+Ric=\lambda g.$

It is called gradient Ricci soliton and the function $f$ is called potential
function. If $f$ is constant then the gradient Ricci soliton is simply the
Einstein manifold. An Einstein manifold with constant $f$ is called trivial
Ricci soliton. The Ricci soliton is called shrinking, steady or expanding if
$\lambda>0$, $\lambda=0$ or $\lambda<0$, respectively. If the gradient Ricci
soliton is shrinking, then after rescaling the metric $g$ we assume that
$\lambda=\frac{1}{2}$. Then (2) takes the form

(3) $\nabla f^{2}+Ric=\frac{1}{2}g.$

Shrinking Ricci soliton has great impact on the topology of the manifold. One
of the first result in this connection is the Myers’s classical theorem which
says that a compact Riemannian manifold with positive Ricci curvature has
finite fundamental group. Fernández-López and García-Río [12] proved, using
universal cover, that the fundamental group is finite for the compact
shrinking Ricci soliton. Zhang [30] also gave an alternative proof. The same
result has been proved earlier by Lott [19] for gradient Ricci soliton. Naber
[21] generalized this result and proved for bounded Ricci curvature. Wylie
[27] proved the strongest result that any complete shrinking Ricci soliton has
finite fundamental group. Locally conformally flat Ricci soliton have been
investigated immensely in the last few years. Cao et al. [5] Petersan and
Wylie [20] and Zhang [28] classified locally conformally flat gradient Ricci
soliton and they showed that it is isometric to
$\mathbb{S}^{n},\mathbb{R}^{n},\mathbb{R}\times\mathbb{S}^{n-1}$ or one of
their quotients. For more results on Ricci soliton see [22, 28].
In this paper, we have showed that a complete conformally flat gradient
shrinking Ricci soliton with reciprocal of the potential functional being
subharmonic and the manifold satisfying linear volume growth is isometric to
the Euclidean sphere. As a result of this, we have derived that the
conformally flat restriction condition can be neglected in $3$-dimension
gradient shrinking Ricci soliton and in case of $4$-dimension, half-
conformally flat condition is needed. We have also deduced a compactness
criteria for gradient shrinking Ricci soliton with quadratic volume growth.
The set of all finitely integrable function on the manifold $M$ is denoted by
$L^{1}(M)$, i.e.,

$L^{1}(M)=\Big{\\{}f|M\rightarrow\mathbb{R}:\int_{M}f<+\infty\Big{\\}}.$

######  Theorem 1.1.

Let $(M,g,f)$ be an $n$-dimensional complete non-flat gradient shrinking Ricci
soliton with scalar curvature $R$ satisfying $R\in L^{1}(M)$. If $\frac{1}{f}$
is subharmonic, then $M$ is compact Einstein manifold with positive Ricci
curvature and scalar curvature $R=\frac{n}{2}$. Furthermore, if $M$ is
conformally flat, then $M$ is a space form. Furthermore, the manifold is
isometric to the $n$-dimensional Euclidean sphere $\mathbb{S}^{n}$.

Since, in $3$-dimension, the Weyl conformal curvature vanishes, the following
result can be stated:

######  Corollary 1.1.1.

A $3$-dimensional gradient shrinking Ricci soliton $(M,g,f)$ with $R\in
L^{1}(M)$ and $\Delta\frac{1}{f}\geq 0$ is simply connected.

A Riemannian manifold $M$ is said to have linear volume growth, if

$\limsup_{r\rightarrow\infty}\frac{Vol(B(p,r))}{r}<+\infty,$

for every $p\in M$. Various properties of Riemannian manifold with linear
volume growth have been studied by many authors, for example see [23, 24]. If
we use the linear volume growth condition, then the curvature restricted
condition of Theorem 1.1 can be omitted.

######  Theorem 1.2.

Suppose $(M,g,f)$ is an $n$-dimensional complete conformally flat gradient
shrinking Ricci soliton with linear volume growth. If $\frac{1}{f}$ is
subharmonic, then $M$ is isometric to $\mathbb{S}^{n}$.

The Weyl curvature tensor $W$ in $4$-dimensional oriented Riemannian manifold
can be decomposed into two irreducible components, $W^{+}$ and $W^{-}$ under
the action of special orthogonal group [3]. If $W^{+}$ or $W^{-}$ vanishes,
then the manifold is called half-conformally flat. Hitchin [16] classified
$4$-dimensional compact half-conformally flat Einstein manifolds.

######  Theorem 1.3.

[3, Theorem 13.30] Let $M$ be a compact half-conformally flat four dimensional
Einstein manifold with positive scalar curvature. Then $M$ is isometric to
$\mathbb{S}^{4}$ or $\mathbb{C}P^{2}$.

Now by using the Theorem 1.1, Theorem 1.2 and the Theorem 1.3, we conclude the
following:

######  Corollary 1.3.1.

Suppose $(M,g,f)$ is a $4$-dimensional complete gradient shrinking Ricci
soliton with $\frac{1}{f}$ is subharmonic. If $M$ is half-conformally flat and
any one of the following conditions holds
(i) $R\in L^{1}(M)$
(ii) $M$ has linear volume growth,
then $M$ is isometric to $\mathbb{S}^{4}$ or $\mathbb{C}P^{2}$. Therefore, $M$
has trivial fundamental group.

Recently, Zhang [31] proved the following result:

######  Theorem 1.4.

[31] Suppose $(M,g)$ is a compact four-dimensional Einstein manifold with
$Ric=g$ and the sectional curvature $K\leq\frac{1}{12}+\frac{\sqrt{30}}{8}$.
Then $(M,g)$ is isometric to either $\mathbb{S}^{4}$ or $\mathbb{RP}^{4}$ or
$\mathbb{CP}^{2}$.

By using the Theorem 1.4 and proof of Theorem 1.2, the following can be
stated:

######  Corollary 1.4.1.

Let $(M,g,f)$ be a $4$-dimensional complete gradient shrinking Ricci soliton
with linear volume growth and sectional curvature
$K\leq\frac{1}{12}+\frac{\sqrt{30}}{8}$. If $\frac{1}{f}$ is subharmonic, then
$(M,g)$ is isometric to either $\mathbb{S}^{4}$ or $\mathbb{RP}^{4}$ or
$\mathbb{CP}^{2}$.

Next using the proof of the Theorem 1.1, we will prove the following
compactness criteria:

######  Theorem 1.5.

If an $n$-dimensional gradient shrinking Ricci soliton $(M,g,f)$ realizes the
following conditions
(i) $Vol(B_{r}(p))\leq Cr^{2}$ for all $p\in M$ and $r>0$, where $C>0$ is a
constant,
(ii) $f$ satisfies the inequality $\Delta f^{2}\leq|\nabla f|^{2}$,
then $M$ is compact.

##  2\. Proof of the results

To prove the Theorem 1.1, we need the following result:

######  Theorem 2.1.

[32] A conformally flat compact manifold with metric of positive Ricci
curvature and constant scalar curvature is a space form.

###### Proof of Theorem 1.1.

Taking trace of (3), we get

(4) $R+\Delta f=\frac{n}{2}.$

Again Hamilton [14] proved that the potential function $f$ in a shrinking
gradient Ricci soliton, after adding some constant if necessary, satisfies

(5) $R+|\nabla f|^{2}=f.$

Since $M$ is shrinking Ricci soliton, the scalar curvature is non-negative
[28] and also $f$ satisfies the asymptotic behavior [4], i.e., there exist
positive constants $C$ and $C^{\prime}$ such that $f$ obeys the estimation

(6) $\frac{1}{4}\Big{[}(d(x,p)-C)_{+}\Big{]}^{2}\leq
f(x)\leq\frac{1}{4}\Big{(}d(x,p)+C^{\prime}\Big{)}^{2},$

where $C_{+}=max(C,0)$. Again since $M$ is non-flat [8], $R>0$. Thus (5),
yields $f>0$. Therefore, in normal coordinate we calculate

$\Big{(}\frac{1}{f}\Big{)}_{i}=-\frac{1}{f^{2}}f_{i},\
\nabla\Big{(}\frac{1}{f}\Big{)}=-\frac{1}{f^{2}}\nabla f,$

and

$\Big{(}\frac{1}{f}\Big{)}_{ii}=\Big{(}-\frac{1}{f^{2}}f_{i}\Big{)}_{i}=\frac{2}{f^{3}}f_{i}^{2}-\frac{1}{f^{2}}f_{ii}.$

Therefore

(7) $\Delta\Big{(}\frac{1}{f}\Big{)}=\frac{2}{f^{3}}|\nabla
f|^{2}-\frac{1}{f^{2}}\Delta f.$

Similarly, we calculate

$\Big{(}\frac{1}{f^{2}}\Big{)}_{i}=-\frac{2}{f^{3}}f_{i},\
\nabla\Big{(}\frac{1}{f^{2}}\Big{)}=-\frac{2}{f^{3}}\nabla f,$

and

$\Big{(}\frac{1}{f^{2}}\Big{)}_{ii}=\Big{(}-\frac{2}{f^{3}}f_{i}\Big{)}_{i}=\frac{6}{f^{4}}f_{i}^{2}-\frac{2}{f^{3}}f_{ii}.$

Hence

(8) $\Delta\Big{(}\frac{1}{f^{2}}\Big{)}=\frac{6}{f^{4}}|\nabla
f|^{2}-\frac{2}{f^{3}}\Delta f.$

For a fixed constant $c>0$, as introduced in [9], we define the function $u$
in $M$ by

$u(x)=R(x)-\frac{c}{f(x)}-\frac{nc}{f^{2}(x)}\ \text{ for }x\in M.$

Hence, using weighted Laplacian, (4) and (5), we obtain

(9) $\displaystyle\Delta_{f}(f^{-1})$ $\displaystyle=$ $\displaystyle\Delta
f^{-1}-g(\nabla f,\nabla f^{-1})$ $\displaystyle=$
$\displaystyle\frac{2}{f^{3}}|\nabla f|^{2}-\frac{1}{f^{2}}\Delta f-g(\nabla
f,-\frac{1}{f^{2}}\nabla f)$ $\displaystyle=$
$\displaystyle\frac{2}{f^{3}}|\nabla f|^{2}-\frac{1}{f^{2}}\Delta
f+\frac{1}{f^{2}}|\nabla f|^{2}$ $\displaystyle=$
$\displaystyle\frac{2}{f^{3}}|\nabla
f|^{2}-\frac{1}{f^{2}}(\frac{n}{2}-R)+\frac{1}{f^{2}}(f-R)$ $\displaystyle=$
$\displaystyle\frac{2}{f^{3}}|\nabla f|^{2}-\frac{n}{2f^{2}}+\frac{1}{f}$
$\displaystyle=$
$\displaystyle\frac{1}{f}-\frac{1}{f^{2}}\Big{(}\frac{n}{2}-\frac{2|\nabla
f|^{2}}{f}\Big{)}$

and

(10) $\displaystyle\Delta_{f}(f^{-2})$ $\displaystyle=$ $\displaystyle\Delta
f^{-2}-g(\nabla f,\nabla f^{-2})$ $\displaystyle=$
$\displaystyle\frac{6}{f^{4}}|\nabla f|^{2}-\frac{2}{f^{3}}\Delta
f+\frac{2}{f^{3}}|\nabla f|^{2}$ $\displaystyle=$
$\displaystyle\frac{6}{f^{4}}|\nabla
f|^{2}-\frac{2}{f^{3}}(\frac{n}{2}-R)+\frac{2}{f^{3}}(f-R)$ $\displaystyle=$
$\displaystyle\frac{6}{f^{4}}|\nabla f|^{2}-\frac{n}{f^{3}}+\frac{2}{f^{2}}$
$\displaystyle=$
$\displaystyle\frac{2}{f^{2}}-\frac{1}{f^{3}}\Big{(}n-\frac{6|\nabla
f|^{2}}{f}\Big{)}.$

From [11], we get

$\Delta R=g(\nabla f,\nabla R)+R-2|Ric|^{2}.$

Then, we conclude

(11) $\Delta_{f}R=R-2|Ric|^{2}.$

Now, (9) and (11) together imply that for any constant $c>0$

(12) $\displaystyle\Delta_{f}\Big{(}R-\frac{c}{f}\Big{)}$ $\displaystyle=$
$\displaystyle
R-2|Ric|^{2}-\frac{c}{f}+\frac{c}{f^{2}}\Big{(}\frac{n}{2}-\frac{2|\nabla
f|^{2}}{f}\Big{)}$ $\displaystyle\leq$ $\displaystyle
R-\frac{c}{f}+\frac{c}{f^{2}}\Big{(}\frac{n}{2}-\frac{2|\nabla
f|^{2}}{f}\Big{)}.$

Then, from above inequality and (10), we obtain

(13) $\displaystyle\Delta_{f}u$ $\displaystyle=$
$\displaystyle\Delta_{f}R-\Delta_{f}\Big{(}\frac{c}{f}\Big{)}-\Delta_{f}\Big{(}\frac{nc}{f^{2}}\Big{)}$
$\displaystyle\leq$ $\displaystyle
R-\frac{c}{f}+\frac{c}{f^{2}}\Big{(}\frac{n}{2}-\frac{2|\nabla
f|^{2}}{f}\Big{)}-\Big{(}\frac{2}{f^{2}}-\frac{1}{f^{3}}\Big{(}n-\frac{6|\nabla
f|^{2}}{f}\Big{)}\Big{)}$ $\displaystyle=$ $\displaystyle
u-\frac{cn}{f^{3}}\Big{(}\frac{f}{2}-n\Big{)}-\frac{c}{f^{4}}\Big{(}2f+6n\Big{)}|\nabla
f|^{2}.$

Now for sufficiently small $c>0$ we can consider $u$ as positive inside
$B(p,C+3n)$. Let us suppose that there is a point $x_{0}\in M-B(p,C+3n)$ where
$u$ attains its negative minimum. Then evaluation (13) at the point $x_{0}$
and using maximum principle, we get $\frac{f(x_{0})}{2}\leq n$. Again (6)
shows that $f(x_{0})\geq\frac{9n^{2}}{4}$, which leads to a contradiction.
Hence, we conclude that $u\geq 0$ in $M$. Therefore for sufficiently small
$c>0$, we obtain

(14) $R\geq\frac{c}{f}+\frac{cn}{f^{2}}\geq\frac{c}{f}.$

Now construct a cut-off function, introduced in [7], $\varphi_{r}\in
C^{\infty}_{0}(B(p,2r))$ in $M$ with the property

$\begin{cases}0\leq\varphi_{r}\leq 1&\text{ in }B(p,2r)\\\
\varphi_{r}=1&\text{ in }B(p,r)\\\ |\nabla\varphi_{r}|\leq\frac{1}{2r}&\text{
in }B(p,2r)\\\ \varphi_{r}=0&\text{ in }\partial B(p,2r).\end{cases}$

Since $\frac{1}{f}$ is subharmonic, it follows that

$\displaystyle 0$ $\displaystyle\leq$
$\displaystyle\int_{B_{2r}(p)}\varphi_{r}^{2}\frac{1}{f}\Delta\Big{(}\frac{1}{f}\Big{)}$
$\displaystyle=$
$\displaystyle-\int_{B_{2r}(p)}\varphi_{r}^{2}|\nabla\frac{1}{f}|^{2}-2\int_{B_{2r}(p)}\varphi_{r}\frac{1}{f}g\Big{(}\nabla\varphi_{r},\nabla\frac{1}{f}\Big{)}.$

Then we have

$\displaystyle\int_{B_{2r}(p)}\varphi_{r}^{2}|\nabla\frac{1}{f}|^{2}\leq-2\int_{B_{2r}(p)}\varphi_{r}\frac{1}{f}g\Big{(}\nabla\varphi_{r},\nabla\frac{1}{f}\Big{)}$
$\displaystyle\leq
2\Big{(}\int_{B_{2r}(p)}\varphi_{r}^{2}|\nabla\frac{1}{f}|^{2}\Big{)}^{1/2}\Big{(}\int_{B_{2r}(p)}\big{(}\frac{1}{f}\big{)}^{2}|\nabla\varphi_{r}|^{2}\Big{)}^{1/2}.$

Therefore, we obtain

(15) $\displaystyle\int_{B_{r}(p)}|\nabla\frac{1}{f}|^{2}$ $\displaystyle\leq$
$\displaystyle\int_{B_{2r}(p)}\varphi_{r}^{2}|\nabla\frac{1}{f}|^{2}\leq
4\int_{B_{2r}(p)}\big{(}\frac{1}{f}\big{)}^{2}|\nabla\varphi_{r}|^{2}$
$\displaystyle\leq$
$\displaystyle\frac{1}{r^{2}}\int_{B_{2r}(p)}\big{(}\frac{1}{f}\big{)}^{2}.$

Now from (14), we get

$R\geq\frac{nc}{f^{2}}\text{ in }M.$

Then the inequality (15) reduces to

(16)
$\int_{B_{r}(p)}|\nabla\frac{1}{f}|^{2}\leq\frac{1}{ncr^{2}}\int_{B_{2r}(p)}R.$

Since $R\in L^{1}(M)$, taking limit $r\rightarrow\infty$ in both sides we get
$|\nabla\frac{1}{f}|^{2}=0$, which implies that $f$ is constant. Therefore, we
conclude, using Myers’s theorem, from (3) that $M$ is compact with constant
scalar curvature $R=\frac{n}{2}$. Also from (3) we state that $M$ is Einstein.
Now for the second part, if $M$ is conformally flat, then we conclude from
Theorem 2.1, that $M$ is the space form. Again, the Ricci curvature of $M$ is
positive, therefore, $M$ has constant positive sectional curvature. Thus, $M$
is isometric to $\mathbb{S}^{n}$. ∎

###### Proof of Theorem 1.2.

The scalar curvature $R$ of gradient shrinking Ricci soliton satisfies the
following inequality [4]:

$\int_{B_{r}(p)}R\leq\frac{n}{2}Vol(B(p,r)),$

for any $p\in M$. Therefore, (16) implies that

(17)
$\int_{B_{r}(p)}|\nabla\frac{1}{f}|^{2}\leq\frac{1}{2cr^{2}}Vol(B(p,2r)).$

On the other hand, the manifold $M$ satisfies the linear volume growth
condition. Hence, taking the limit $r\rightarrow\infty$ of the above
inequality, we get $|\nabla\frac{1}{f}|^{2}=0$. Therefore, using the same
argument of the proof of Theorem 1.1, we conclude that the result. ∎

For the proof of Theorem 1.5, the following results are needed:

######  Lemma 2.2.

If $f>0$ and $\Delta f^{2}\leq|\nabla f|^{2}$, then $(\frac{1}{f})$ is
subharmonic.

###### Proof.

Given that $\Delta f^{2}\leq|\nabla f|^{2}.$ Then

$\displaystyle 0$ $\displaystyle\geq$ $\displaystyle\Delta f^{2}-6|\nabla
f|^{2}$ $\displaystyle=$ $\displaystyle(\Delta f^{2}-2|\nabla f|^{2})-4|\nabla
f|^{2}$ $\displaystyle=$ $\displaystyle f\Delta f-2|\nabla f|^{2}.$

Therefore

$2|\nabla f|^{2}-f\Delta f\geq 0.$

Since $f>0$, the above inequality implies that

$\frac{2}{f^{3}}|\nabla f|^{2}-\frac{1}{f^{2}}\Delta f\geq 0.$

The right side of the above inequality is $\Delta(\frac{1}{f})$, which is non-
negative. Therefore, $\frac{1}{f}$ is subharmonic. ∎

######  Theorem 2.3.

[17] Let $f\in C^{\infty}(M)$ be a nonnegative subharmonic function. Then
either $f$ is constant or

$\liminf_{r\rightarrow\infty}\frac{1}{r^{2}}\int_{B_{r}(p)}f^{p}dV=+\infty,$

for every $p\in M$ and $p>1.$

###### Proof of Theorem 1.5.

Since $M$ is non-flat shrinking Ricci soliton, the scalar curvature is
positive [8, 28]. Thus (5) implies that $f>0$. Then from (14), we have

(18) $\frac{nc}{f^{2}}\leq R.$

Furthermore, Lemma 2.2 implies $\frac{1}{f}$ is subharmonic. Again

$\Delta\Big{(}\frac{1}{f^{2}}\Big{)}=2\frac{1}{f}\Delta\frac{1}{f}+2|\nabla\frac{1}{f}|^{2}.$

Hence, $\frac{1}{f^{2}}$ is also subharmonic. Now, the scalar curvature
estimation in [4] implies that

(19)
$\int_{B_{r}(p)}\frac{nc}{f^{2}}dV\leq\int_{B_{r}(p)}RdV\leq\frac{n}{2}Vol(B(p,r)).$

Since $M$ has quadratic volume growth, i.e., $Vol(B(p,r))\leq C_{1}r^{2}$,
where $C_{1}$ is a positive constant, the above inequality implies that

$c\int_{B_{r}(p)}\frac{1}{f^{2}}dV\leq\frac{C_{1}}{2}r^{2},$

i.e.,

(20) $\frac{1}{r^{2}}\int_{B_{r}(p)}\frac{1}{f^{2}}dV\leq C^{\prime},$

where $C^{\prime}=C_{1}/2c$ is a positive constant. This implies that

$\liminf_{r\rightarrow\infty}\frac{1}{r^{2}}\int_{B_{r}(p)}\frac{1}{f^{2}}dV<+\infty,$

and it contradicts Theorem 2.3. Thus we conclude that $f$ must be constant.
Hence (3) implies $Ric=\frac{1}{2}g$. Therefore Myers’s Theorem [1] implies
that the manifold $M$ is compact. ∎

## References

  * [1] Ambrose, W., A Theorem of Myers, Duke Math. J., 24 (1957), 345–348.
  * [2] Aubin, T., Some nonlinear problems in Riemannian geometry, Springer, 2013.
  * [3] Besse, A. L., Einstein manifolds, Springer, 1987.
  * [4] Cao, H. D. and Zhou, D., On complete gradient shrinking Ricci soliton, J. Diff. Geom., 85 (2010), 175–185.
  * [5] Cao, X., Wang, D. and Zhang, Z., On locally conformally flat gradient shrinking Ricci solitons, Commun. Contemp. Math., 13 (2011), 269–282.
  * [6] Carrillo, J. A and Ni, L., Sharp logarithmic Sobolev inequalities on gradient solitons and applications, Comm. Anal. Geom., 17(4) (2010), 1–33.
  * [7] Cheeger, J. and Colding, T. H., Lower bounds on Ricci curvature and the almost rigidity of warped products, Ann. Math., 144(1) (1996), 189–237.
  * [8] Chen, B. L., Strong uniqueness of the Ricci flow, J. Diff. Geom., 82(2) (2009), 362–382.
  * [9] Chow, B., Lu, P. and Yang, B., Lower bounds for the scalar curvature of noncompact gradient Ricci soliton, C. R. Math. Acad. Sci. Paris, 349(23-24) (2011), 1265–1267.
  * [10] Chow, B. and Lu, P., A bound for the order of the fundamental group of a complete noncompact Ricci shrinker, Proc. Amer. Math. Soc., 144 (2016), 2623–2625.
  * [11] Eminenti M, La Nave G and Mantegazza C., Ricci solitons: the equation point of view, Manu. Math., 127(1) (2008), 345–367.
  * [12] Fernández-López, M. and García-Río, E., A remark on compact Ricci solitons, Math. Ann., 340 (2008), 893–896.
  * [13] Hamilton, R. S., The Ricci flow on surfaces, Mathematics and general relativity (Santa Cruz, CA, 1986), Contemp. Math., 71 (1988), 237-262.
  * [14] Hamilton, R. S., The formation of singularities in the Ricci flow, Surveys in differential geom., 2 (1995), 7–136.
  * [15] Hamilton, R. S., Three-manifolds with positive Ricci curvature, J. Differ. Geom., 17 (1982), 255–306.
  * [16] Hitchin, N., Compact four-dimensional Einstein manifolds, J. Differ. Geom., 9(3) (1974), 435–441.
  * [17] Karp, L., Subharmonic functions on real and complex manifolds, Math. Z., 179(4) (1982), 535–554.
  * [18] Lohkamp, J., The space of negative scalar curvature metrics, Invent. Math., 110(1) (1992), 57–95.
  * [19] Lott, J., Some geometric properties of the Bakry-Émery-Ricci tensor, Comment. Math. Helv., 78(4) (2003), 865–883.
  * [20] Petersen, P. and Wylie, W., On the classification of gradient Ricci solitons, Geom. Topol., 14 (2010), 2277–2300.
  * [21] Naber, A., Some Geometry and Analysis on Ricci Solitons, arXiv:math/0612532.
  * [22] Shaikh, A. A. and Mondal, C. K.,  On Ricci solitons whose potential is convex, Accepted in Proc. Indian Acad. Sci. (Math. Sci.).
  * [23] Sormani, C.,  Busemann functions on manifolds with lower Ricci curvature bounds and minimal volume growth, J. Diff. Geom., 48 (1998), 557–585.
  * [24] Sormani, C., Harmonic functions on manifolds with nonnegative Ricci curvature and Linear volume Growth, Pacific J. Math., 192 (2000), 183–189.
  * [25] Wan, J., On the fundamental group of complete manifolds with almost Euclidean volume growth, Proc. Amer. Math. Soc., 147(10) (2019), 4493–4498.
  * [26] Wu, P., On the potential function of gradient steady Ricci solitons, J. Geom. Anal., 23 (2013), 221–228.
  * [27] Wylie, W., Complete shrinking Ricci solitons have finite fundamental group, Proc. Amer. Math. Soc., 136(5) (2007), 1803–1806.
  * [28] Zhang, Z. H., On the completeness of gradient Ricci solitons, Proc. Amer. Math. Soc., 137(8) (2009), 2755–2759.
  * [29] Zhang, Z. H., Gradient shrinking solitons with vanishing Weyl tensor, Pacific J. Math., 242(1) (2009), 189–200.
  * [30] Zhang, Z., On the finiteness of the fundamental group of a compact shrinking Ricci soliton, Collq. Math., 107(2) (2007), 297–299.
  * [31] Zhang, Z., Einstein four-manifolds with sectional curvature bounded from above, J. Geom. Anal, 30(1) (2020), 182–196.
  * [32] Goldberg, S. I., on conformally flat spaces with definite Ricci curvature, Kodai. Math. J., 21 (1969), 226–232.