# Asymptotic theory for quadratic variation of
harmonizable fractional stable processes

Andreas Basse-O’Connor Department of Mathematics, University of Aarhus,
E-mail<EMAIL_ADDRESS>Mark Podolskij Department of Mathematics, University
of Luxembourg, E-mail<EMAIL_ADDRESS>Podolskij gratefully
acknowledges financial support of ERC Consolidator Grant 815703 “STAMFORD:
Statistical Methods for High Dimensional Diffusions”.

###### Abstract

In this paper we study the asymptotic theory for quadratic variation of a
harmonizable fractional $\alpha$-stable process. We show a law of large
numbers with a non-ergodic limit and obtain weak convergence towards a Lévy-
driven Rosenblatt random variable when the Hurst parameter satisfies
$H\in(1/2,1)$ and $\alpha(1-H)<1/2$. This result complements the asymptotic
theory for fractional stable processes investigated in e.g. [4, 5, 6, 7, 17,
21].

Key words: fractional processes, harmonizable processes, limit theorems,
quadratic variation, stable Lévy motion.

AMS 2010 subject classifications: 60F05, 60F15, 60G22, 60G48, 60H05.

##  1 Introduction

Since the seminal work of Mandelbrot and Van Ness [20] fractional stochastic
processes and fields have gained a lot of attention in the probabilistic and
statistical literature. One of the key results of [20] was the introduction of
the fractional Brownian motion, the most famous fractional stochastic process.
The authors show that, up to scaling, fractional Brownian motion with Hurst
parameter $H\in(0,1)$ is the unique self-similar centered Gaussian process
with stationary increments. When the Gaussianity assumption is replaced by
$\alpha$-stable distributions the aforementioned class becomes much richer:
Rosinski [24] shows that it consists of (mixed) moving average $\alpha$-stable
processes, $\alpha$-stable harmonizable processes and processes of a third
kind (which lack an explicit description). The most well-known representative
of the first class is the linear fractional stable motion (lfsm), whose
probabilistic and statistical properties have been studied in numerous
articles. We refer to [4, 5, 6, 7, 17, 22, 23] for a detailed exposition of
limit theorems for functionals of lfsm, and to [1, 10, 12, 18, 19, 21] for
statistical estimation of lfsm and related processes.

In the paper [9] the authors introduced the $\alpha$-stable harmonizable
fractional motion, defined as

$X_{t}=\int_{\mathbb{R}}\frac{\exp(its)-1}{is}|s|^{1-H-1/\alpha}\,dL_{s},$
(1.1)

where $i:=\sqrt{-1}$, $H\in(0,1)$ is the Hurst parameter and
$L_{t}:=L_{t}^{1}+iL_{t}^{2}$, $t\in\mathbb{R}$, denotes an isotropic complex-
valued $\alpha$-stable Lévy motion with $\alpha\in(0,2)$. The processes
$(X_{t})_{t\in\mathbb{R}}$ is self-similar with index $H$, i.e.
$(X_{t})_{t\in\mathbb{R}}\stackrel{{\scriptstyle
d}}{{=}}(a^{-H}X_{at})_{t\in\mathbb{R}}$, and possesses stationary increments.
Probabilistic properties of model (1.1) and their extensions have been studied
in numerous articles, see e.g. [2, 8, 9, 11, 15] among others. However, to the
best of our knowledge, asymptotic theory for statistics of harmonizable
fractional $\alpha$-stable processes has not yet being investigated in the
literature. In contrast to lfsm, harmonizable fractional $\alpha$-stable
processes are not ergodic, which indicates potential difficulties when
studying limit theorems.

This paper aims to derive asymptotic results for the quadratic statistic

$\sum_{j=0}^{n-1}\|X_{j}-X_{j-1}\|^{2},$ (1.2)

and our main statement is the following theorem:

######  Theorem 1.1.

Let $(X_{t})_{t\geq 0}$ be the $\alpha$-stable harmonizable fractional motion
introduced in (1.1) with parameters $\alpha\in(0,2)$ and $H\in(0,1)$.

  1. (i)

As $n\to\infty$,

$\frac{1}{n}\sum_{j=0}^{n-1}\|X_{j}-X_{j-1}\|^{2}\stackrel{{\scriptstyle\mathbb{P}}}{{\longrightarrow}}U:=2\int_{\mathbb{R}}|s|^{-2H-2/\alpha}\left(1-\cos(s)\right)d([L^{1}]_{s}+[L^{2}]_{s}),$
(1.3)

where $[L^{k}]$ denotes the quadratic variation of $L^{k}$ for $k=1,2$.

  2. (ii)

For $H>1/2$ and $\alpha(1-H)<1/2$ we obtain the weak convergence

$\displaystyle
n^{2-2H}\left(\frac{1}{n}\sum_{j=0}^{n-1}\|X_{j}-X_{j-1}\|^{2}-U\right)$ (1.4)
$\displaystyle\qquad\quad\stackrel{{\scriptstyle\raisebox{-1.0pt}{{\tiny
d}}}}{{\longrightarrow}}2{\textrm{Re}}\Big{(}\int_{\mathbb{R}^{2}}\frac{1-\exp(i(s-u))}{i(s-u)}|su|^{1-H-1/\alpha}1_{\\{u<s\\}}\,d\overline{L}_{u}\,dL_{s}\Big{)}.$
(1.5)

######  Remark 1.2.

(A) The result of Theorem 1.1(i) is in sharp contrast to the setting of short
memory processes. To make this comparison concrete, let us consider a sequence
$(Y_{i})_{i\in\mathbb{N}}$ of independent $\alpha$-stable isotropic random
variables with $\alpha\in(0,2)$. According to the Stable Central Limit
Theorem, cf. [13, Theorem 5.25],

$\frac{1}{n^{2/\alpha}}\sum_{j=0}^{n-1}\|Y_{j}\|^{2}\stackrel{{\scriptstyle\raisebox{-1.0pt}{{\tiny
d}}}}{{\longrightarrow}}Z,$ (1.6)

where $Z$ is an $(\alpha/2)$-stable random variable. That is, for independent
$\alpha$-stable random variables, the quadratic variation
$\sum_{j=0}^{n}\|Y_{j}\|^{2}$ is of the order $n^{\alpha/2}$, whereas
$\sum_{j=0}^{n-1}\|X_{j}-X_{j-1}\|^{2}$ is of the order $n$, cf. Theorem
1.1(i). This difference may be viewed as an indication of the strong
dependence between the random variables $(X_{j}-X_{j-1})_{j\in\mathbb{N}}$.

(B) We believe that the randomness of the limit in (1.3) is connected to the
non-ergodicity of the process $X$.

(C) The weak convergence result of (1.4) is similar in spirit to the Gaussian
case. Indeed, the quadratic variation of the fractional Brownian motion with
Hurst parameter $H$ is asymptotically distributed according to the Rosenblatt
random variable when $H\in(3/4,1)$ with convergence rate $n^{2H-2}$; the
classical Rosenblatt random variable possesses the same representation as the
right hand side of (1.4) where $L$ is replaced by a complex-valued Gaussian
measure. We also recall that there is no convergence regime of the form (1.4)
in the setting of lfsm (cf. [4, 5]).

(D) The proof of Theorem 1.1 heavily relies on the quadratic form of the
statistic. Extensions to more general functional forms and statistical
estimation of the process $X$ are therefore far from obvious. $\diamond$

The rest of the paper is structured as follows. In Chapter 2, we discuss the
integration theory and a dominated convergence result we will use in our
proof. Chapter 3 is devoted to the proof of Theorem 1.1.

### Notation

Throughout the paper all random variables are defined on a probability space
$(\Omega,\mathcal{F},\mathbb{P})$. By
$\stackrel{{\scriptstyle\mathbb{P}}}{{\longrightarrow}}$ we denote convergence
in probability, and $\stackrel{{\scriptstyle\raisebox{-1.0pt}{{\tiny
d}}}}{{\longrightarrow}}$ denotes convergence in distribution. For a complex
number $x\in\mathbb{C}$, $\overline{x}$ denotes its complex conjugate and
$\|x\|$ its norm. Furthermore, we denote by $\text{Re}(x)$ (resp.
$\text{Im}(x)$) the real (resp. imaginary) part of $x$.

##  2 Some background

In this section we will collect some important results on stochastic
integration theory. A complex-valued random variable $Z$ is called isotropic
if it holds that $\exp(i\theta)Z\stackrel{{\scriptstyle d}}{{=}}Z$ for any
$\theta\in(0,2\pi)$. Similarly, a stochastic process
$(Z_{t})_{t\in\mathbb{R}}$ is isotropic when any linear combination
$\sum_{j=1}^{d}a_{j}Z_{t_{j}}$ is isotropic. It is easy to see that a complex-
valued Lévy process $L$ with characteristic triplet $(0,0,\nu)$, where $\nu$
denotes the Lévy measure, is isotropic if and only if

$\nu\circ(x\mapsto x\exp(i\theta))^{-1}=\nu\qquad\forall\theta\in(0,2\pi).$
(2.1)

In particular, $L^{1}$ and $L^{2}$ both have the same symmetric
$\alpha$-stable distribution. For any measurable function
$f:\mathbb{R}\to\mathbb{C}$ we define

$\displaystyle\int_{\mathbb{R}}f(s)\,dL_{s}$ (2.2)
$\displaystyle\quad:=\Big{(}\int_{\mathbb{R}}\text{Re}(f(s))\,dL^{1}_{s}-\int_{\mathbb{R}}\text{Im}(f(s))\,dL^{2}_{s}\Big{)}+i\Big{(}\int_{\mathbb{R}}\text{Re}(f(s))\,dL^{2}_{s}+\int_{\mathbb{R}}\text{Im}(f(s))\,dL^{1}_{s}\Big{)},$
(2.3)

whenever the four real-valued integrals on the right-hand side exist. We note
that the integral $\int_{\mathbb{R}}f(s)\,dL_{s}$ is well defined if and only
if $\int_{\mathbb{R}}\|f(s)\|\,dL^{1}_{s}$ is well defined, which holds if
$\int_{\mathbb{R}}\|f(s)\|^{\alpha}\,ds<\infty$. The next result, which
follows by [25, Theorem 2], gives a criterion for the existence of the double
integral:

######  Proposition 2.1.

Let $\psi:\mathbb{R}\to(0,\infty)$ be a measure function satisfying
$\int_{\mathbb{R}}\psi(s)^{\alpha}\,ds=1$, and suppose that
$(Z_{t})_{t\in\mathbb{R}}$ is a two-dimensional symmetric $\alpha$-stable Lévy
process with $Z_{t}=(Z^{1}_{t},Z^{2}_{t})$. The double integral
$\int_{\mathbb{R}^{2}}f(u,v)\,dZ^{1}_{u}\,dZ^{2}_{v}$ of a measurable function
$f:\mathbb{R}^{2}\to\mathbb{R}$ vanishing outside
$\\{(u,v)\in\mathbb{R}:u<v\\}$ exists when the following condition holds:

$\int_{\mathbb{R}}\int_{\mathbb{R}}|f(s,t)|^{\alpha}\Big{[}1+\log_{+}\Big{(}\frac{|f(s,t)|}{\psi(s)\psi(t)}\Big{)}\Big{]}\,ds\,dt<\infty,$
(2.4)

where $\log_{+}(x):=\log(x)$ for $x\geq 1$ and $\log_{+}(x):=0$ for $x<1$.

We remark that [25, Theorem 2] is only stated in the case where $Z^{1}=Z^{2}$;
however, the extension to $Z^{1}\neq Z^{2}$ follows directly using a
polarization techniques similar to the polarization equality for real numbers
$xy=((x+y)^{2}-(x-y)^{2})/4$. Proposition 2.4 extends directly to complex-
valued functions and processes as well.

######  Remark 2.2.

Let $g:\mathbb{R}\to\mathbb{C}$ is measurable function with
$\|g\|_{\alpha}^{\alpha}:=\int_{\mathbb{R}}\|g(x)\|^{\alpha}dx<\infty$ and set
$f(u,v)=\overline{g}(u)g(v)1_{\\{u<v\\}}$ for all $u,v\in\mathbb{R}$. Then we
have the following statements.

(A) The function $f$ satisfies condition (2.4), and hence the double integral
$\int_{\mathbb{R}^{2}}f(u,v)\,d\overline{L}_{u}\,dL_{v}$ is well defined, cf.
Proposition 2.1. Indeed, (2.4) follows by choosing
$\psi(s)=\|g(s)\|/\|g\|_{\alpha}$.

(B) The following identity holds

$\Big{\|}\int_{\mathbb{R}}g(s)\,dL_{s}\Big{\|}^{2}=2\text{Re}\Big{(}\int_{\mathbb{R}^{2}}\overline{g}(u)g(s)1_{\\{u<s\\}}\,d\overline{L}_{u}\,dL_{s}\Big{)}+\int_{\mathbb{R}}\|g(s)\|^{2}\,d([L^{1}]_{s}+[L^{2}]_{s}).$
(2.5)

To show (2.5) we note that $Y_{t}:=\int_{-\infty}^{t}g(u)\,dL_{u}$ exists for
all $t\in\mathbb{R}\cup\\{\infty\\}$. For all $w<t$ in $\mathbb{R}$ we have by
integration by parts that

$\displaystyle Y_{t}\overline{Y}_{t}-Y_{w}\overline{Y}_{w}={}$
$\displaystyle\int_{w}^{t}Y_{s-}\,d\overline{Y}_{s}+\int_{w}^{t}\overline{Y}_{s-}\,dY_{s}+[Y,\overline{Y}]_{t}-[Y,\overline{Y}]_{w}$
(2.6) $\displaystyle={}$ $\displaystyle
2\text{Re}\Big{(}\int_{u}^{t}\overline{Y}_{s-}\,dY_{s}\Big{)}+[Y,\overline{Y}]_{t}-[Y,\overline{Y}]_{w},$
(2.7)

where $[Z,\tilde{Z}]$ denotes the quadratic variation of two complex-valued
semimartingales $Z$ and $\tilde{Z}$. If we let $t\to\infty$ and $w\to-\infty$
in (2.7) we obtain

$Y_{\infty}\overline{Y}_{\infty}=2\text{Re}\Big{(}\int_{\mathbb{R}}\overline{Y}_{s-}\,dY_{s}\Big{)}+[Y,\overline{Y}]_{\infty}.$
(2.8)

By reconizing the integral on the right-hand side of (2.8) as a double
integral and calculating the quadratic variation $[Y,\overline{Y}]$ we derive
the identity

$\|Y_{\infty}\|^{2}=Y_{\infty}\overline{Y}_{\infty}=2\text{Re}\Big{(}\int_{\mathbb{R}^{2}}\overline{g}(u)g(v)1_{\\{u<v\\}}\,d\overline{L}_{u}\,dL_{v}\Big{)}+\int_{\mathbb{R}}\|g(s)\|^{2}\,d([L^{1}]_{s}+[L^{2}]_{s}),$
(2.9)

where we have used that $[L,\overline{L}]_{t}=[L^{1}]_{t}+[L^{2}]_{t}$.
Equation (2.9) implies (2.5). $\diamond$

We will also use the dominated convergence result [14, Theorem 6.2(c)], which
reads as follows.

######  Proposition 2.3.

Assume that $(h_{n})_{n\geq 1}$, $h$ and $g$ are deterministic measurable
functions $\mathbb{R}^{2}\to\mathbb{R}$ supported on
$\\{(u,v)\in\mathbb{R}^{2}\\!:u<v\\}$ all satisfying condition (2.4). Suppose
moreover that $(Z_{t})_{t\in\mathbb{R}}$ is a two-dimensional symmetric
$\alpha$-stable Lévy process with $Z_{t}=(Z^{1}_{t},Z^{2}_{t})$. If
$h_{n}(u,v)\to h(u,v)$ for Lebesgue almost all points $(u,v)\in\mathbb{R}^{2}$
and $\sup_{n}\|h_{n}\|\leq g$, then

$\int_{\mathbb{R}^{2}}h_{n}(u,v)\,dZ^{1}_{u}\,dZ^{2}_{v}\stackrel{{\scriptstyle\mathbb{P}}}{{\longrightarrow}}\int_{\mathbb{R}^{2}}h(u,v)\,dZ_{u}^{1}\,dZ^{2}_{v}\qquad\text{as
}n\to\infty.$ (2.10)

##  3 Proof of Theorem 1.1

Throughout this section all positive constants are denoted by $C$, although
they may change from line to line. We start with the identity

$X_{j}-X_{j-1}=\int_{\mathbb{R}}\exp(ijs)r(s)\,dL_{s}\qquad\text{with}\qquad
r(s):=\frac{1-\exp(-is)}{is}|s|^{1-H-1/\alpha}.$ (3.1)

In particular, we have that

$\|r(s)\|\leq C\left(|s|^{1-H-1/\alpha}1_{\\{|s|\leq
1\\}}+|s|^{-H-1/\alpha}1_{\\{|s|>1\\}}\right),$ (3.2)

which implies the statement

$\|r\|_{\alpha}<\infty.$ (3.3)

By using the representations (3.1) and (2.5) on the function
$g:\mathbb{R}\to\mathbb{C}$ given by $g(s)=\exp(ijs)r(s)$ we deduce the
identity

$\|X_{j}-X_{j-1}\|^{2}=V_{j}+U,$ (3.4)

where

$\displaystyle U$
$\displaystyle:=2\int_{\mathbb{R}}|s|^{-2H-2/\alpha}\left(1-\cos(s)\right)d([L^{1}]_{s}+[L^{2}]_{s})\qquad\text{and}$
$\displaystyle V_{j}$
$\displaystyle:=2\text{Re}\Big{(}\int_{\mathbb{R}^{2}}\exp(ij(s-u))r(s)\overline{r}(u)1_{\\{u<s\\}}\,d\overline{L}_{u}\,dL_{s}\Big{)}.$
(3.5)

###  3.1 Proof of convergence in (1.3)

By the decomposition (3.4) we have

$\frac{1}{n}\sum_{j=0}^{n-1}\|X_{j}-X_{j-1}\|^{2}=U+\frac{1}{n}\sum_{j=0}^{n-1}V_{j},$
(3.6)

and thus we only need to prove that
$n^{-1}\sum_{j=0}^{n-1}V_{j}\stackrel{{\scriptstyle\mathbb{P}}}{{\longrightarrow}}0$.
We first note that

$\frac{1}{n}\sum_{j=0}^{n-1}\exp(ij(s-u))=\frac{1-\exp(in(s-u))}{n(1-\exp(i(s-u)))}\to
0\qquad\text{as }n\to\infty,$ (3.7)

for any $s-u\not\in\\{2\pi m:~{}m\in\mathbb{Z}\\}$, using the geometric
series. On the other hand, we have

$\|\exp(ij(s-u))r(s)r(u)\|\leq\|r(s)\|\|r(u)\|.$ (3.8)

In view of (3.3) and Remark 2.2(A), we can apply Proposition 2.3 to conclude
that

$n^{-1}\sum_{j=0}^{n-1}V_{j}\stackrel{{\scriptstyle\mathbb{P}}}{{\longrightarrow}}0,$
(3.9)

which implies (1.3).

###  3.2 Proof of weak convergence

We will divide the proof into several steps. We write
$g_{n}(x):=(1-\exp(inx))/(1-\exp(ix))$ and define the following complex-valued
functions

$\displaystyle h_{n}(s,u)$
$\displaystyle:=n^{1-2H}g_{n}(s-u)r(s)\overline{r}(u)1_{\\{u<s\\}},$ (3.10)
$\displaystyle h(s,u)$
$\displaystyle:=\frac{1-\exp(i(s-u))}{i(s-u)}|su|^{1-H-1/\alpha}1_{\\{u<s\\}}.$
(3.11)

By the representation (3.4) we have

$n^{2-2H}\left(\frac{1}{n}\sum_{j=0}^{n-1}\|X_{j}-X_{j-1}\|^{2}-U\right)=n^{1-2H}\sum_{j=0}^{n-1}V_{j}=2\textrm{Re}\Big{(}\int_{\mathbb{R}^{2}}h_{n}(s,u)\,d\overline{L}_{u}\,dL_{s}\Big{)},$
(3.12)

where the last equality follows from the identity in (3.7) and the definition
of $h_{n}$. Thus, to show the limit theorem (1.4) it is enough to show the
convergence

$\int_{\mathbb{R}^{2}}h_{n}(s,u)\,d\overline{L}_{u}\,dL_{s}\stackrel{{\scriptstyle\raisebox{-1.0pt}{{\tiny
d}}}}{{\longrightarrow}}\int_{\mathbb{R}^{2}}h(s,u)\,d\overline{L}_{u}\,dL_{s}.$
(3.13)

The following subsections will be denoted to the proof of (3.13), and to this
aim we define

$\gamma:=1-H-1/\alpha.$

####  3.2.1 Existence of the limit

In this subsection we will show that the limiting integral
$\int_{\mathbb{R}^{2}}h(s,u)\,d\overline{L}_{u}\,dL_{s}$ from (3.13) is well-
defined by proving that $h$ satisfies condition (2.4). Set

$\psi(s):=c\left(|s|^{-r}\mathbbm{1}_{\\{|s|>1\\}}+\mathbbm{1}_{\\{|s|\leq
1\\}}\right)$ (3.14)

for all $s\in\mathbb{R}$, where $r$ is any real number satisfying
$r>\alpha^{-1}$ and $c=c_{r}>0$ is chosen such that
$\int_{\mathbb{R}}\psi(s)^{\alpha}\,ds=1$. To show condition (2.4) with
respect to the function $\psi$ we will use the following lemma.

######  Lemma 3.1.

Let $f:\mathbb{R}^{2}\to\mathbb{C}$ be a complex-valued function satisfying
the inequality

$\|f(s,u)\|\leq
C|su|^{-r_{1}}\left(1_{[s-1,s)}(u)+|s-u|^{-r_{2}}1_{(-\infty,s-1)}(u)\right),$
(3.15)

where the numbers $r_{1},r_{2}$ satisfy the conditions $r_{1}\in(1/2,1)$ and
$2r_{1}+r_{2}>2$. Then it holds that

$\int_{\mathbb{R}^{2}}\|f(s,u)\|\,du\,ds<\infty.$ (3.16)

###### Proof.

According to condition (3.15) we have that $\|f(s,u)\|\leq
f_{1}(s,u)+f_{2}(s,u)$ and hence we need to show the integrability of the
functions $f_{1}$ and $f_{2}$. We have that

$\int_{s-1}^{s}|u|^{-r_{1}}du\leq
C\left(1_{[-1,1]}(s)+s^{-r_{1}}1_{[-1,1]^{c}}(s)\right)$

since $r_{1}\in(1/2,1)$. Hence, we conclude

$\int_{\mathbb{R}^{2}}f_{1}(s,u)\,du\,ds\leq
C\int_{\mathbb{R}}|s|^{-r_{1}}\left(1_{[-1,1]}(s)+s^{-r_{1}}1_{[-1,1]^{c}}(s)\right)ds$

and the latter is finite as $r_{1}\in(1/2,1)$.

Now we turn our attention to the function $f_{2}$. We use the substitution
$s=xu$ to obtain the bound

$\displaystyle\int_{u+1}^{\infty}|s|^{-r_{1}}|s-u|^{-r_{2}}ds$
$\displaystyle=|u|^{1-r_{1}-r_{2}}\int_{\frac{u+1}{u}}^{\infty}|x|^{-r_{1}}|x-1|^{-r_{2}}dx$
(3.17) $\displaystyle\leq
C\left(1_{[-1,1]}(u)+\left(|u|^{-r_{1}}+|u|^{1-r_{1}-r_{2}}\right)1_{[-1,1]^{c}}(u)\right)$
(3.18)

Indeed, the latter integral is finite since $r_{1}+r_{2}>2r_{1}+r_{2}-1>1$.
Thus, we deduce

$\int_{\mathbb{R}^{2}}f_{2}(s,u)\,ds\,du\leq
C\int_{\mathbb{R}}|u|^{-r_{1}}\left(1_{[-1,1]}(u)+\left(|u|^{-r_{1}}+|u|^{1-r_{1}-r_{2}}\right)1_{[-1,1]^{c}}(u)\right)ds<\infty$

due to assumptions $r_{1}\in(1/2,1)$ and $2r_{1}+r_{2}>2$. This completes the
proof. ∎

Now, we proceed with the proof of existence. First of all, we observe that the
function $\|h(s,u)\|^{\alpha}$ satisfies the condition (3.15) with

$r_{1}=\alpha\gamma\qquad\text{and}\qquad r_{2}=\alpha.$

It is easy to see that assumptions $r_{1}\in(1/2,1)$ and $2r_{1}+r_{2}>2$ are
satisfied due to $H>1/2$ and $\alpha(1-H)<1/2$. On the other hand, for any
$\varepsilon>0$ we readily deduce that

$\log_{+}\left(\frac{\|h(s,u)\|}{\psi(s)\psi(u)}\right)\leq
C\left(|s|^{\varepsilon}+|s|^{-\varepsilon}\right)\left(|u|^{\varepsilon}+|u|^{-\varepsilon}\right).$
(3.19)

When we choose $\varepsilon>0$ small enough, we conclude that the function

$\|h(s,u)\|^{\alpha}\left(1+\log_{+}\left(\frac{\|h(s,u)\|}{\psi(s)\psi(u)}\right)\right)$

also satisfies the conditions of Lemma 3.1. This proves the existence of the
double integral $\int_{\mathbb{R}^{2}}h(s,u)\,d\overline{L}_{u}\,dL_{s}$, cf.
Proposition 2.1.

####  3.2.2 Main decomposition

Now, we show the weak limit theorem in (1.4). For this purpose, we introduce
the sets

$\mathcal{Z}:=\\{2\pi j:~{}j\geq
0\\}\qquad\text{and}\qquad\mathcal{Z}_{\delta}:=\\{x\in\mathbb{R}_{+}:~{}|x-\mathcal{Z}|<\delta\\}$

for some fixed $\delta\in(0,1)$. For $x\in\mathbb{R}_{+}$ we also use the
notation $[x]:=2\pi j$ where $2\pi j\in\mathcal{Z}$ is the closest number to
$x$. We obtain the following estimate for the function $g_{n}$:

$\displaystyle\|g_{n}(x)\|$ $\displaystyle\leq
C\left(\left\|\frac{1-\exp(in(x-[x]))}{i(x-[x])}\right\|1_{\mathcal{Z}_{\delta}}(x)+\delta^{-1}1_{\mathcal{Z}^{c}_{\delta}}(x)\right)$
$\displaystyle\leq
C\left(n1_{\mathcal{Z}_{\delta}\cap\mathcal{Z}_{1/n}}+\frac{1}{|x-[x]|}1_{\mathcal{Z}_{\delta}\cap\mathcal{Z}_{1/n}^{c}}+\delta^{-1}1_{\mathcal{Z}^{c}_{\delta}}(x)\right)$
$\displaystyle\leq
C\left(n^{1-a}|x-[x]|^{-a}1_{\mathcal{Z}_{\delta}}+\delta^{-1}1_{\mathcal{Z}^{c}_{\delta}}(x)\right)$
(3.20)

for any $a\in(0,1)$. We decomposition the double integral on the right-hand
side of (3.12) as follows

$\int_{\mathbb{R}^{2}}h_{n}(s,u)\,d\overline{L}_{u}\,dL_{s}=Z_{n}^{(1)}+Z_{n}^{(2)}+Z_{n}^{(3)},$
(3.21)

with

$\displaystyle Z_{n}^{(1)}$
$\displaystyle:=\int_{(-\delta,\delta)^{2}}h_{n}(s,u)\,d\overline{L}_{u}\,dL_{s}$
(3.22) $\displaystyle Z_{n}^{(2)}$
$\displaystyle:=\int_{(-\delta,\delta)^{c}\times(-\delta,\delta)^{c}}h_{n}(s,u)\,d\overline{L}_{u}\,dL_{s}$
(3.23) $\displaystyle Z_{n}^{(3)}$
$\displaystyle:=\int_{\big{(}(-\delta,\delta)^{c}\times(-\delta,\delta)\big{)}\cup\big{(}(-\delta,\delta)\times(-\delta,\delta)^{c}\big{)}}h_{n}(s,u)\,d\overline{L}_{u}\,dL_{s}.$
(3.24)

In the following, we will prove that $Z_{n}^{(1)}$ is the dominating term and
show its convergence, while we will prove that
$\|Z_{n}^{(2)}\|+\|Z_{n}^{(3)}\|\stackrel{{\scriptstyle\mathbb{P}}}{{\longrightarrow}}0$.
The following identity in distribution

$\int_{A}f(s,u)\,d\overline{L}_{u}\,dL_{s}\stackrel{{\scriptstyle
d}}{{=}}n^{-2/\alpha}\int_{nA}f(s/n,u/n)\,d\overline{L}_{u}\,dL_{s},$ (3.25)

where $nA:=\\{(ns,nu)\\!:~{}(s,u)\in A\\}$, will be used in our proof. The
identity (3.25) follows from the self-similarity of $L$; more precisely, from
the fact that $(L_{t/n})_{t\in\mathbb{R}}$ equals
$(n^{-1/\alpha}L_{t})_{t\in\mathbb{R}}$ in finite dimensional distributions.

####  3.2.3 Term $Z_{n}^{(1)}$

We apply the distributional identity (3.25) to deduce that

$Z_{n}^{(1)}\stackrel{{\scriptstyle d}}{{=}}n^{-2/\alpha}\int_{(-\delta
n,\delta n)^{2}}h_{n}(s/n,u/n)\,d\overline{L}_{u}\,dL_{s}.$ (3.26)

For any fixed $(s,u)$ with $s-u\not\in\mathcal{Z}$ we obtain the convergence

$\lim_{n\to\infty}n^{\gamma}r(s/n)=|s|^{\gamma}\qquad\text{and}\qquad\lim_{n\to\infty}n^{-1}g_{n}(s-u)=\frac{1-\exp(i(s-u))}{i(s-u)}.$
(3.27)

Hence, we conclude that

$\lim_{n\to\infty}n^{-2/\alpha}h_{n}(s/n,u/n)=h(s,u).$ (3.28)

On the other hand, we have that

$\|n^{-2/\alpha}h_{n}(s/n,u/n)\|^{\alpha}\leq
C|su|^{\gamma}\left(\mathbbm{1}_{\\{s-1\leq
u<s\\}}+\mathbbm{1}_{\\{u<s-1\\}}|s-u|^{-\alpha}\right)=:f(s,u).$ (3.29)

As in the proof of existence we deduce that the function $f$ satisfies the
condition (2.4). Thus,

$Z_{n}^{(1)}\stackrel{{\scriptstyle\raisebox{-1.0pt}{{\tiny
d}}}}{{\longrightarrow}}\int_{\mathbb{R}^{2}}\frac{1-\exp(i(s-u))}{i(s-u)}|su|^{1-H-1/\alpha}1_{\\{u<s\\}}\,d\overline{L}_{u}\,dL_{s}$
(3.30)

due to Proposition 2.3.

####  3.2.4 Term $Z_{n}^{(2)}$

We will use the approximation (3.20). First of all, we have that

$\|h_{n}(s,u)\|1_{\mathcal{Z}_{\delta}^{c}}(s-u)\leq
C\delta^{-1}n^{1-2H}|su|^{\gamma-1}1_{\\{u<s\\}}.$ (3.31)

Note that $n^{1-2H}\to 0$ since $H>1/2$. On the other hand, the function
$f(s,u):=|su|^{\gamma-1}1_{\\{u<s\\}}$ satisfies the condition

$\int_{(-\delta,\delta)^{c}\times(-\delta,\delta)^{c}}f(s,u)^{\alpha}\left(1+\log_{+}\left(\frac{f(s,u)}{\psi(s)\psi(u)}\right)\right)du\,ds<\infty$

where we use an estimate of the type (3.19) for the log term. For all
$a\in(0,1)$ we have

$\displaystyle\|h_{n}(s,u)\|1_{\mathcal{Z}_{\delta}}(s-u)$ $\displaystyle\leq
Cn^{2-2H-a}|su|^{\gamma-1}\left|(s-u)-[s-u]\right|^{-a}1_{\\{u<s\\}}.$ (3.32)

It is easy to show that the function
$(s,u)\mapsto|su|^{\gamma-1}|(s-u)-[s-u]|^{-a}1_{((-\delta,\delta)^{c})^{2}}(s,u)1_{\\{u<s\\}}$
satisfies the condition (2.4) when $a\alpha<1$. Choosing
$a=1/\alpha-\varepsilon$ we deduce that

$n^{2-2H-a}=n^{2-2H-1/\alpha+\varepsilon}\to 0\qquad\text{as }n\to\infty,$

if we choose $\varepsilon>0$ small enough since $\alpha(1-H)<1/2$. We,
therefore, conclude that there exists a measurable function $\phi$ satisfying
(2.4) such that $|h_{n}(u,v)|\leq\phi(u,v)$ for all
$(u,v)\in(-\delta,\delta)^{c}\times(-\delta,\delta)^{c}$. Since $h_{n}(u,v)\to
0$ Lebesgue almost surely, cf. (3.31)–(3.32), Proposition 2.3 yields
$Z_{n}^{(2)}\stackrel{{\scriptstyle\mathbb{P}}}{{\longrightarrow}}0$.

####  3.2.5 Term $Z_{n}^{(3)}$

The proof is similar to the one in Subsection 3.2.4. For simplicity we only
consider the integration region $(-\delta,\delta)\times(-\delta,\delta)^{c}$;
the $(-\delta,\delta)^{c}\times(-\delta,\delta)$ case is analogue. For all
$(u,s)\in(-\delta,\delta)\times(-\delta,\delta)^{c}$:

$\|h_{n}(s,u)\|1_{\mathcal{Z}_{\delta}^{c}}(s-u)\leq
C\delta^{-1}n^{1-2H}|u|^{\gamma}|s|^{\gamma-1}1_{\\{u<s\\}}.$ (3.33)

We have $n^{1-2H}\to 0$ and the function
$f(s,u):=|u|^{\gamma}|s|^{\gamma-1}1_{(-\delta,\delta)}(u)1_{(-\delta,\delta)^{c}}(s)1_{\\{u<s\\}}$
satisfies the condition (2.4). On the other hand, we have

$\|h_{n}(s,u)\|1_{\mathcal{Z}_{\delta}}(s-u)\leq
Cn^{2-2H-a}|u|^{\gamma}|s|^{\gamma-1}\left|(s-u)-[s-u]\right|^{-a}1_{\\{u<s\\}}.$
(3.34)

The bounds (3.33) and (3.34) imply the existence of a measurable function $f$
satisfying (2.4) such that $|h_{n}(u,v)|\leq f(u,v)$ for all
$(u,v)\in(-\delta,\delta)\times(-\delta,\delta)^{c}$. Since a similar bound
holds on $(-\delta,\delta)^{c}\times(-\delta,\delta)$ we conclude, cf.
Proposition 2.3, that
$Z_{n}^{(3)}\stackrel{{\scriptstyle\mathbb{P}}}{{\longrightarrow}}0$ since
$h_{n}\to 0$ (cf. (3.33)–(3.34)).

By the decomposition (3.21) and the convergence results of Subsections
3.2.3–3.2.5 we obtain the limit theorem (3.13). This completes the proof of
(1.4) and of Theorem 1.1. ∎

## References

  * [1] A. Ayachea and J. Hamoniera (2012): Linear fractional stable motion: A wavelet estimator of the $\alpha$ parameter. Statistics and Probability Letters 82, 1569–1575.
  * [2] A. Ayachea and Y. Xiao (2016): Harmonizable fractional stable fields: Local nondeterminism and joint continuity of the local times. Stochastic Processes and Their Applications 126(1), 171–185.
  * [3] A. Basse-O’Connor, T. Grønbæk and M. Podolskij (2021): Local asymptotic self-similarity for heavy-tailed harmonizable fractional Lévy motions. ESAIM: Probability and Statistics 25, 286–297.
  * [4] A. Basse-O’Connor, C. Heinrich and M. Podolskij (2019): On limit theory for functionals of stationary increments Lévy driven moving averages. Electronic Journal of Probability 24(79), 1–42.
  * [5] A. Basse-O’Connor, R. Lachièze-Rey and M. Podolskij (2017): Power variation for a class of stationary increments Lévy driven moving averages. Annals of Probability, 45(6B), 4477–4528.
  * [6] A. Basse-O’Connor and M. Podolskij (2017): On critical cases in limit theory for stationary increments Lévy driven moving averages. Festschrift for Bernt Øksendal, Stochastics 81(1), 360–383.
  * [7] A. Basse-O’Connor, M. Podolskij and C. Thäle (2020): A Berry-Esseén theorem for partial sums of functionals of heavy-tailed moving averages. Electronic Journal of Probability 25(31), 1–31.
  * [8] A. Benassi, S. Cohen and J. Istas (2002): Identification and properties of real harmonizable fractional Lévy motions. Bernoulli 8(1), 97–115.
  * [9] S. Cambanis and M. Maejima (1989): Two classes of self-similar stable processes with stationary increments. Stochastic Processes and Their Applications 32(2), 305–329.
  * [10] T.T.N. Dang and J. Istas (2017): Estimation of the Hurst and the stability indices of a $H$-self-similar stable process. Electronic Journal of Statistics 11, 4103–4150.
  * [11] M. Dozzi and G. Shevchenko (2011): Real harmonizable multifractional stable process and its local properties. Stochastic Processes and Their Applications 121(7), 1509–1523.
  * [12] D. Grahovac, N.N. Leonenko and M.S. Taqqu (2015): Scaling properties of the empirical structure function of linear fractional stable motion and estimation of its parameters. J. Stat. Phys. 158 105–119.
  * [13] J. Hoffmann-Jørgensen (1994): Probability with a View Toward Statistics Vol. I. Chapman & Hall, New York.
  * [14] O. Kallenberg and J. Szulga (1989): Multiple integration with respect to Poisson and Lévy processes. Probability Theory and Related Fields 83(1-2), 101–134.
  * [15] N. Kôno and M. Maejima (1991): Hölder continuity of sample paths of some self-similar stable processes. Tokyo J. Math. 14(1), 93–100.
  * [16] S. Kwapien and W.A. Woyczynski (1987): Double stochastic integrals, random quadratic forms and random series in Orlicz spaces. Annals of Probability 15(3), 1072–1096.
  * [17] M.M. Ljungdahl and M. Podolskij (2018): A limit theorem for a class of stationary increments Lévy moving average process with multiple singularities. Modern Stochastics: Theory and Applications 5(3), 297–316.
  * [18] M.M. Ljungdahl and M. Podolskij (2020): A minimal contrast estimator for the linear fractional stable motion. Statistical Inference for Stochastic Processes 23, 381–413.
  * [19] M.M. Ljungdahl and M. Podolskij (2022): Multidimensional parameter estimation of heavy-tailed moving averages. Scandinavian Journal of Statistics 49(2), 593–624.
  * [20] B.B. Mandelbrot and J. W. Van Ness (1968): Fractional Brownian motions, fractional noises and applications. SIAM Rev. 10, 422–437.
  * [21] S. Mazur, D. Otryakhin and M. Podolskij (2020): Estimation of the linear fractional stable motion. Bernoulli 26(1), 226–252.
  * [22] V. Pipiras and M.S. Taqqu (2003): Central limit theorems for partial sums of bounded functionals of infinite-variance moving averages. Bernoulli 9, 833–855.
  * [23] V. Pipiras, M.S. Taqqu and P. Abry (2007): Bounds for the covariance of functions of infinite variance stable random variables with applications to central limit theorems and wavelet-based estimation. Bernoulli 13(4), 1091–1123.
  * [24] J. Rosinski (1995): On the structure of stationary stable processes. Annals of Probability 23(3), 1163–1187.
  * [25] G. Samorodnitsky and M. Taqqu (1990): Multiple stable integrals of Banach-valued functions. Journal of Theoretical Probability 3, 267–287.


# DeepFIB: Self-Imputation for Time Series
Anomaly Detection

Minhao Liu, Zhijian Xu, Qiang Xu
{mhliu, zjxu<EMAIL_ADDRESS>

###### Abstract

Time series (TS) anomaly detection (AD) plays an essential role in various
applications, e.g., fraud detection in finance and healthcare monitoring. Due
to the inherently unpredictable and highly varied nature of anomalies and the
lack of anomaly labels in historical data, the AD problem is typically
formulated as an unsupervised learning problem. The performance of existing
solutions is often not satisfactory, especially in data-scarce scenarios. To
tackle this problem, we propose a novel self-supervised learning technique for
AD in time series, namely _DeepFIB_. We model the problem as a _Fill In the
Blank_ game by masking some elements in the TS and imputing them with the
rest. Considering the two common anomaly shapes (point- or sequence-outliers)
in TS data, we implement two masking strategies with many self-generated
training samples. The corresponding self-imputation networks can extract more
robust temporal relations than existing AD solutions and effectively
facilitate identifying the two types of anomalies. For continuous outliers, we
also propose an anomaly localization algorithm that dramatically reduces AD
errors. Experiments on various real-world TS datasets demonstrate that DeepFIB
outperforms state-of-the-art methods by a large margin, achieving up to
$65.2\%$ relative improvement in F1-score.

##  1 Introduction

Anomaly detection (AD) in time series (TS) data has numerous applications
across various domains. Examples include fault and damage detection in
industry (Hundman et al., 2018), intrusion detection in cybersecurity (Feng &
Tian, 2021), and fraud detection in finance (Zheng et al., 2018) or healthcare
(Zhou et al., 2019), to name a few.

Generally speaking, an anomaly/outlier is an observation that deviates
considerably from some concept of normality (Ruff et al., 2021). The somewhat
“vague” definition itself tells the challenges of the AD problem arising from
the rare and unpredictable nature of anomalies. With the lack of anomaly
labels in historical data, most AD approaches try to learn the expected values
of time-series data in an unsupervised manner (Bl’azquez-Garc’ia et al.,
2021). Various techniques use different means (e.g., distance-based methods
(Angiulli & Pizzuti, 2002), predictive methods (Holt, 2004; Yu et al., 2016;
Deng & Hooi, 2021) or reconstruction-based methods (Shyu et al., 2003;
Malhotra et al., 2016; Zhang et al., 2019; Shen et al., 2021)) to obtain this
expected value, and then compute how far it is from the actual observation to
decide whether or not it is an anomaly.

While existing solutions have shown superior performance on some time series
AD tasks, they are still far from satisfactory. For example, for the six ECG
datasets in (Keogh et al., 2005), the average F1-score of state-of-the-art
solutions (Kieu et al., 2019; Shen et al., 2021) with model ensembles are
barely over $40\%$. Other than the TS data’ complexity issues, one primary
reason is that the available data is often scarce while deep learning
algorithms are notoriously data-hungry.

Recently, self-supervised learning (SSL) that enlarges the training dataset
without manual labels has attracted lots of attention, and it has achieved
great success in representation learning in computer vision (Zhang et al.,
2016; Pathak et al., 2016; Chen et al., 2020), natural language processing
(Devlin et al., 2019), and graph learning (Hu et al., 2020) areas. There are
also a few SSL techniques for time series analysis proposed in the literature.
Most of them (Falck et al., 2020; Saeed et al., 2021; Fan et al., 2020) craft
contrastive TS examples for classification tasks. (Deldari et al., 2021) also
leverages contrastive learning for change point detection in time series.

(a) point-wise (b) sequence-wise

Figure 1: Anomalies in time series.

While interesting, the above SSL techniques do not apply to the AD task
because detecting anomalies in time series requires fine-grained models at the
element level. In this work, inspired by the context encoder for visual
feature learning (Pathak et al., 2016) and the BERT model for language
representation learning (Devlin et al., 2019), we propose a novel self-
supervised learning technique for time series anomaly detection, namely
_DeepFIB_. To be specific, we model the problem as a _Fill In the Blank_ game
by masking some elements in the TS and imputing them with other elements. This
is achieved by revising the TS forecasting model _SCINet_ (Liu et al., 2021)
for the TS imputation task, in which the masked elements are regarded as
missing values for imputation. Such self-imputation strategies facilitate
generating a large amount of training samples for temporal relation
extraction. As anomalies in time series manifest themselves as either discrete
points or subsequences (see Fig. 1), correspondingly, we propose two kinds of
masking strategies and use them to generate two pre-trained models. They are
biased towards recovering from _point-wise_ anomalies (_DeepFIB-p_ model for
_point outliers_) and _sequence-wise_ anomalies (_DeepFIB-s_ model for
_continuous outliers_), respectively. To the best of our knowledge, this is
the first SSL work for time series anomaly detection.

Generally speaking, AD solutions have difficulty detecting sequence-wise
anomalies because it is hard to tell the real outliers from their neighboring
normal elements due to their interplay. To tackle this problem, we propose a
novel anomaly localization algorithm to locate the precise start and end
positions of continuous outliers. As a post-processing step, we conduct a
local search after determining the existence of sequence-wise anomalies within
a timing window with our _DeepFIB-s_ model. By doing so, the detection
accuracy for continuous outliers is significantly improved.

We conduct experiments on several commonly-used time series benchmarks, and
results show that DeepFIB consistently outperforms state-of-the-art solutions.
In particular, the average F1-score of DeepFIB for the six ECG datasets is
more than $62\%$, achieving nearly $50\%$ relative improvement.

##  2 Related Work

In this section, we mainly discuss recent deep learning-based time series AD
approaches. A comprehensive survey on the traditional techniques can be found
in (Gupta et al., 2014).

Figure 2: Existing time series anomaly detection architectures.

Existing anomaly detection approaches can be broadly categorized into three
types (see Fig. 2): (i) _Density-based_ methods consider the normal instances
compact in the latent space and identify anomalies with one-class classifiers
or likelihood measurements (Su et al., 2019; Shen & Kwok, 2020; Feng & Tian,
2021). (ii) _Reconstruction-based_ methods use recurrent auto-encoders (RAE)
(Malhotra et al., 2016; Yoo et al., 2021; Kieu et al., 2019; Shen et al.,
2021; Zhang et al., 2019) or deep generative models such as recurrent VAEs
(Park et al., 2018) or GANs (Li et al., 2019; Zhou et al., 2019) for
reconstruction. The reconstruction errors are used as anomaly scores. (iii)
_Prediction-based_ methods rely on predictive models (Bontemps et al., 2016;
Deng & Hooi, 2021; Chen et al., 2021) and use the prediction errors as anomaly
scores.

While the above methods have been successfully used in many real-world
applications, practical AD tasks still have lots of room for improvement,
especially in data-scarce scenarios. Unlike existing AD approaches, the
proposed mask-and-impute method in _DeepFIB_ exploits the unique property of
TS data that missing values can be effectively imputed (Fang & Wang, 2020). By
constructing many training samples via self-imputation, _DeepFIB_ extracts
robust temporal relations of TS data and improves AD accuracy dramatically.
Moreover, for the more challenging sequence-wise anomalies, most prior work
assumes a user-defined fixed-length for anomaly subsequences (Cook et al.,
2020) or simplifies the problem by stating all the continuous outliers have
been correctly detected as long as one of the points is detected (Su et al.,
2019; Shen & Kwok, 2020). In _DeepFIB_ , we lift these assumptions and try to
locate the exact location of sequence-wise anomalies.

##  3 Method

In this section, we first introduce the overall self-imputation framework in
DeepFIB and then discuss the separate AD models for detecting point- and
sequence-wise anomalies with different mask-and-impute strategies, namely
_DeepFIB-p_ and _DeepFIB-s_ , respectively. Next, we describe the TS
imputation method used in _DeepFIB_ , based on an existing TS forecasting
approach _SCINet_ (Liu et al., 2021). Finally, we present our anomaly
localization algorithm for continuous outliers.

###  3.1 Self-Imputation for Anomaly Detection

Given a set of multivariate time series wherein
$X_{s}=\left\\{x_{1},x_{2},...,x_{T_{s}}\right\\}\epsilon\mathbb{R}^{d\times
T_{s}}$ ($T_{s}$ is the length of the $s_{th}$ time series $X_{s}$), the
objective of the AD task is to find all anomalous points
$x_{t}\in\mathbb{R}^{d}$ ($d$ is the number of variates) and anomalous
subsequences $X_{t,\tau}=\left\\{x_{t-\tau+1},...,x_{t}\right\\}$.

The critical issue to solve the above problem is obtaining an expected value
for each element in the TS, which requires a large amount of training data to
learn from, especially for deep learning-based solutions. However, time-series
data are often scarce, significantly restricting the effectiveness of
learning-based AD solutions.

DeepFIB is a simple yet effective SSL technique to tackle the above problem.
We model this problem as a _Fill In the Blank_ game by randomly masking some
elements in the TS and imputing them with the rest. Such self-imputation
strategies generate many training samples from every time series and hence
dramatically improve temporal learning capabilities.

Figure 3: Self-imputation strategies of _DeepFIB-p_ and _DeepFIB-s_.

In particular, we propose to train two self-imputation models (Fig. 3), biased
towards point- and sequence-wise anomalies in the TS data, respectively.

  * •

_DeepFIB-p model_ targets point outliers, as shown in Fig. 3(a), in which we
mask discrete elements and rely on the local temporal relations extracted from
neighboring elements for reconstruction. For each time series $X_{s}$, we
generate $M$ training samples by masking it $M$ times with randomly-selected
yet _non-overlapping_ $\frac{d\times T_{s}}{M}$ elements.

  * •

_DeepFIB-s model_ targets continuous outliers, as shown in Fig. 3(b), in which
we mask continuous elements and rely on predictive models for reconstruction.
For each time series $X_{s}$, we evenly divide it into $N$ _non-overlapping_
sub-sequences as
$\left\\{X^{d\times{\frac{T_{s}}{N}}}_{s,i},i\in\left[0,N-1\right]\right\\}$
and generate $N$ training samples by masking one of them each time.

During training, for each time series $X_{s}$, we obtain a set of non-
overlapped imputed data with the above model and integrate them together
results in a reconstructed time series $\widehat{X_{s}}$ (i.e.,
$\widehat{X_{s}}$_-p_ for _DeepFIB-p_ model and $\widehat{X_{s}}$_-s_ for
_DeepFIB-s_ model). The training loss for both models are defined as the
reconstruction errors between the input time series and the reconstructed one:

$\mathcal{L}=\frac{1}{T_{s}}\sum_{t=1}^{T_{s}}\left\|x_{t}-\widehat{x_{t}}\right\|$
$None$

where $x_{t}$ is the original input value at time step $t$ and the
$\widehat{x_{t}}$ denotes the reconstructed value from the corresponding
model, and $\left\|\cdot\right\|$ is the L1-norm of a vector.

During testing, to detect point outliers with the _DeepFIB-p_ model, we simply
use the residual error as the anomaly score, defined as
$e_{t}=\sum_{i=0}^{d}\left|\widehat{x_{t}}^{i}-x_{t}^{i}\right|$, and when
$e_{t}$ is larger than a threshold value $\lambda_{p}$, time step $t$ is
regarded as an outlier. In contrast, for continuous outliers, we use dynamic
time warping (DTW) (Sakoe & Chiba, 1978) distance metrics as our anomaly
scoring mechanism, which measures the similarity between the input time series
$X$ and reconstructed sequence $\widehat{X}$. If $DTW(X,\widehat{X})$ is above
a threshold value $\lambda_{s}$, a sequence-wise anomaly is detected.

###  3.2 Time Series Imputation in DeepFIB

While the time-series data imputation problem has been investigated for
decades (Fang & Wang, 2020), there are still lots of rooms for improvement and
various deep learning models are proposed recently (Cao et al., 2018; Liu et
al., 2019; Luo et al., 2019).

Figure 4: The structure of the SCINet.

SCINet (Liu et al., 2021)(Fig. 4) is an encoder-decoder architecture motivated
by the unique characteristics of time series data. It incorporates a series of
SCI-Blocks that conduct down-sampled convolutions and interactive learning to
capture temporal features at various resolutions and effectively blend them in
a hierarchical manner. Considering the highly-effective temporal relation
extraction capability of SCINet when compared to other sequence models, we
propose to revise it for the TS imputation task. More details about _SCINet_
can be found in (Liu et al., 2021).

To impute the missing elements from the two masking strategies with _DeepFIB-
p_ and _DeepFIB-s_ models, we simply change the supervisions for the decoder
part accordingly. For point imputation, we use the original input sequence as
the supervision of our _DeepFIB-p_ model, making it a reconstruction
structure. By doing so, the model concentrates more on the local temporal
relations inside the timing window for imputing discrete missing data, as
shown in Fig. 5(a). As for continuous imputation, we propose to change SCINet
as a bidirectional forecasting structure in our _DeepFIB-s_ model, with the
masked sub-sequence as supervision. As shown in Fig. 5(b), the two sub-models,
namely _F-SCINet_ and _B-SCINet_ , are used to conduct forecasting in the
forward and backward directions, respectively. By doing so, the model can
aggregate the temporal features from both directions and learn a robust long-
term temporal relations for imputing continuous missing data.

Figure 5: Time series imputation in DeepFIB.

###  3.3 Anomaly Localization Algorithm

During inference, we use a sliding window with stride $\mu$ to walk through
the time series and find anomalies in each window. For sequence-wise
anomalies, without knowing their positions a priori, we could mask some normal
elements in the window and use those unmasked outliers for prediction (see
Fig. 5(b)), thereby leading to mispredictions. To tackle this problem, we
propose to conduct a local search for the precise locations of the sequence-
wise anomalies.

Figure 6: Anomaly localization algorithm.

As shown in Fig. 6, the _Active window_ are the current input sequence to the
DeepFIB-s model with length $\omega$ ($\omega>\mu$), i.e.,
$X_{t}=\left\\{x_{t},x_{t+1},...,x_{t+\omega-1}\right\\}$ at time step $t$.
When the _DTW_ distance between the original time series in the _Active
window_ and the imputed sequence is above the threshold $\lambda_{s}$, a
sequence-wise anomaly is detected in the current window, and the localization
mechanism is triggered. As the sliding window is moving along the data stream
with stride $\mu$, if no outliers are detected in the previous window, the
start position of the sequence-wise anomaly can only exist at the end of
$X_{t}$ in the window
$\left\\{x_{t+\omega-\mu},...,x_{t+\omega-1},x_{t+\omega-1}\right\\}$ with
length $\mu$. Consequently, by gradually shifting the _Active window_ backward
to include one more element in the _Buffer_ window (see Fig. 6) at a time and
calculating the corresponding _DTW_ distances as
$\left\\{e_{1},...,e_{i},...,e_{\mu}\right\\}$, we can find the maximum $i$
with $e_{i}<\lambda_{s}$, indicating the following element after the _Active
window_ starting with $i$ is the start of the anomaly subsequence. The
_Anomaly flag_ is then activated from this position. Similarly, to determine
the ending position of the anomaly subsequence, we keep sliding the _Active
windows_ until we find a window with _DTW_ distance smaller than
$\lambda_{s}$, indicating that the ending position is within
$\left\\{x_{t-\mu},...,x_{t-2},x_{t-1}\right\\}$. Again, we shift the _Active
window_ backwardly one-by-one to include one element of the above window at a
time and calculate the corresponding _DTW_ distance, until we find the ending
position with its _DTW_ distance larger than $\lambda_{s}$.

##  4 Experiments

In this section, we conduct extensive experiments to answer the following two
questions: _Whether DeepFIB outperforms state-of-the-art AD methods_ (Q1)?
_How does each component of DeepFIB affect its performance_ (Q2)?

Table 1: Datasets used in experiments

Datasets |  $\\#$Dim |  $\\#$Train |  $\\#$Test | Anomaly
---|---|---|---|---
2d-gesture | 2 | 8590 | 2420 | 24.63$\%$
Power demand | 1 | 18145 | 14786 | 11.44$\%$
ECG | (A)chfdb_chf01_275 | 2 | 2888 | 1772 | 14.61 $\%$
(B)chfdb_chf13_45590 | 2 | 2439 | 1287 | 12.35 $\%$
(C)chfdbchf15 | 2 | 10863 | 3348 | 4.45 $\%$
(D)ltstdb_20221_43 | 2 | 2610 | 1121 | 11.51 $\%$
(E)ltstdb_20321_240 | 2 | 2011 | 1447 | 9.61 $\%$
(F)mitdb_100_180 | 2 | 2943 | 2255 | 8.38 $\%$
Credit Card | 3 | 142403 | 142404 | 0.173 $\%$

Experiments are conducted on a number of commonly-used benchmark TS datasets,
namely _2d-gesture_ , _Power demand_ , _ECG_ and _Credit Card_ , ranging from
human abnormal behavior detection, power monitoring, healthcare and fraud
detection in finance (see Table 1). As the anomalies in _2d-gesture_ , _Power
demand_ , and _ECG_ are mainly sequence outliers, we apply the _DeepFIB-s_
model on these datasets. In contrast, the _Credit Card_ dataset only contains
point outliers, and hence we use _DeepFIB-p_ model on it.

To make a fair comparison with existing models, we use the standard evaluation
metrics on the corresponding datasets. For _2d-gesture_ , _Power demand_ and
_Credit Card_ , we use precision, recall, and F1-score following (Shen & Kwok,
2020). For _ECG_ datasets, we use the AUROC (area under the ROC curve), AUPRC
(area under the precision-recall curve) and F1-score, following (Shen et al.,
2021). To detect anomalies, we use the maximum anomaly score in each sub-
models over the validation dataset to set the threshold.

More details on experimental settings, additional experimental results and
discussions (e.g., hyperparameter analysis) are presented in the supplementary
materials.

###  4.1 Q1: Comparison with state-of-the-art methods

Table 2: Comparison of anomaly detection performance (as %), on _2d-gesture_
and _Power demand_ datasets. The best results are in bold and the second best
results are underlined.

Methods | 2d-gesture | Power demand
---|---|---
precision | recall | F1-score | precision | recall | F1-score
DAGMM | 25.66 | _80.47_ | 38.91 | 34.37 | 41.72 | 37.69
EncDec-AD | 24.88 | 100.0 | 39.85 | 13.98 | 54.20 | 22.22
LSTM-VAE | 36.62 | 67.72 | 47.54 | 8.00 | 56.66 | 14.03
MADGAN | 29.41 | 76.4 | 42.47 | 13.20 | 60.57 | 27.67
AnoGAN | 57.85 | 46.50 | 51.55 | 20.28 | 44.41 | 28.85
BeatGAN | 55.11 | 45.33 | 49.74 | 8.04 | 76.58 | 14.56
OmniAnomaly | 27.70 | 79.67 | 41.11 | 8.55 | _78.73_ | 15.42
MSCRED | _61.26_ | 59.11 | 60.17 | _55.80_ | 34.32 | 42.50
THOC | 54.78 | 75.00 | _63.31_ | 61.50 | 36.34 | _45.68_
DeepFIB | 93.90 $\pm$ 0.35 | 60.77 $\pm$ 0.24 | 73.79 $\pm$ 0.19 | 52.21 $\pm$ 0.31 | 99.99 $\pm$ 0.01 | 68.60 $\pm$ 0.15

  * •

\- The results of other baselines in the table are extracted from (Shen &
Kwok, 2020)

2d-gesture and Power demand: The results in Table 2 show that the proposed
_DeepFIB-s_ achieves 16.55$\%$ and 50.18$\%$ F1-score improvements on
_2d-gesture_ and _Power demand_ , respectively, compared with the second best
methods.

For _2d-gesture_ , the available training data is limited and the temporal
relations contained in the data are complex (body jitter), making it difficult
to obtain a discriminative representation in AD models. DAGMM (Zong et al.,
2018) shows low performance since it does not consider the temporal
information of the time-series data at all. As for the AD solutions based on
generative models (EncDecAD (Malhotra et al., 2016), LSTM-VAE (Park et al.,
2018), MAD-GAN (Li et al., 2019), AnoGAN (Schlegl et al., 2017), BeatGAN (Zhou
et al., 2019), OmniAnomaly (Su et al., 2019)), they usually require a large
amount of training data, limiting their performance in data-scarce scenario.
Compared to the above methods, the encoder-decoder architecture MSCRED (Zhang
et al., 2019) is relatively easier to train and its AD performance is
considerably higher. Moreover, the recent THOC (Shen & Kwok, 2020) work
further improves AD performance by fusing the multi-scale temporal information
to capture the complex temporal dynamics.

The proposed _DeepFIB-s_ model outperforms all the above baseline methods
since the proposed self-imputation technique allows the model to learn more
robust temporal relations from much more self-generated training samples.
Notably, we also observe that the _precision_ of the _DeepFIB-s_ dominates the
other baselines. We attribute it to the anomaly localization algorithm that
can locate the anomaly’s precise start and end positions, significantly
reducing the false positive rate.

Figure 7: For Power demand dataset, (a) shows two cycles of normal data
(0-1200 frame) wherein each cycle contains 5 peaks. (b) shows two cycles with
anomaly with missing peaks highlighted using red arrows. The waveform of
original data (light blue) is overlaid on the prediction result (light red).
Lower color bars show the ground truth (GT) label and our detection result
(Result).

For _Power demand_ , the data contains many _contextual anomaly 111Contextual
anomalies are observations or sequences that deviate from the expected
patterns within the time series however if taken in isolation they are within
the range of values expected for that signal (Cook et al., 2020)._
subsequences (see Fig. 7). It is quite challenging for existing AD approaches
to learn such context information by extracting temporal features from the
entire time series as a whole. In contrast, the proposed sequence-wise masking
strategy facilitates learning different kinds of temporal patterns, which is
much more effective in detecting such contextual anomalies. As shown in Table
2, the _recall_ of our _DeepFIB-s_ model almost reaches 100$\%$, indicating
all anomalies have been detected. The _precision_ is not the best, and we
argue that some of the false positives are in fact resulted from the poorly
labeled test set (see our supplementary material).

Table 3: Comparison of anomaly detection performance (as %), on ECG datasets.

Metrics | Methods | ECG | Average
---|---|---|---
A | B | C | D | E | F
AUROC | RAE | 64.95 | 75.24 | 68.27 | 60.71 | 77.92 | 44.68 | 65.29
RRN | 69.50 | 72.07 | 68.49 | 47.05 | 78.81 | 47.87 | 63.97
BeatGAN | 66.51 | 73.14 | 58.69 | 59.33 | 82.98 | 44.19 | 64.14
RAE-ensemble | 68.26 | 77.63 | 70.55 | 64.64 | 83.14 | 39.66 | 67.31
RAMED | _73.58_ | _78.82_ | _78.79_ | _69.44_ | _83.36_ | _55.64_ | _73.27_
DeepFIB | 87.60 $\pm$0.85 | 84.40 $\pm$1.23 | 94.05 $\pm$0.73 | 72.55 $\pm$0.54 | 84.81 $\pm$0.62 | 63.23 $\pm$0.12 | 81.11
AUPRC | RAE | 51.84 | 40.32 | 31.23 | 15.54 | 24.17 | 7.76 | 28.48
RRN | 54.90 | 43.13 | 33.49 | 11.63 | 37.68 | 7.93 | 31.46
BeatGAN | 52.50 | 44.94 | 19.01 | 14.84 | 34.46 | 7.66 | 28.90
RAE-ensemble | 56.23 | 54.21 | _49.90_ | _18.47_ | 38.48 | 7.25 | _37.42_
RAMED | 56.23 | 54.23 | 34.63 | 17.78 | 45.78 | _10.59_ | 36.54
DeepFIB | 85.18$\pm$0.63 | 75.48 $\pm$0.56 | 73.47 $\pm$0.67 | 23.14 $\pm$0.45 | _38.27 $\pm$ 0.72_ | 13.16 $\pm$ 0.23 | 51.45
F1 | RAE | 52.51 | 49.03 | 32.79 | 25.43 | 33.63 | 15.47 | 34.81
RRN | 56.08 | 43.48 | 38.30 | 20.64 | 44.37 | 15.47 | 36.39
BeatGAN | 51.93 | 45.18 | 27.99 | 23.67 | 47.02 | 16.68 | 35.41
RAE-ensemble | _56.42_ | _52.40_ | _58.68_ | 27.75 | 44.98 | 15.47 | _42.62_
RAMED | 54.27 | 51.03 | 34.45 | _30.87_ | _52.23_ | _20.63_ | 40.58
DeepFIB | 80.90 $\pm$0.63 | 78.06 $\pm$0.82 | 78.37 $\pm$0.13 | 44.71 $\pm$0.19 | 58.00 $\pm$0.26 | 34.08 $\pm$0.73 | 62.35

  * •

\- The results of other baselines in the table are referred from (Shen et al.,
2021)

ECG(A-F): Compared with (A),(B),(C) datasets, (D),(E),(F) are clearly noisy,
which affect the performance of the anomaly detectors significantly.
Nevertheless, Table 3 shows that _DeepFIB-s_ achieves an average 46.3$\%$
F1-score improvement among all datasets and an impressive 65.2$\%$ improvement
for ECG(F) dataset. There are mainly two reasons: (1) the data is scarce (See
Table 1). Existing AD methods are unable to learn robust temporal relations
under such circumstances. In contrast, the self-imputation training strategy
together with the bidirectional forecasting mechanism used in our _DeepFIB-s_
model can well address this issue; (2) the proposed DTW anomaly score is more
effective in detecting the anomaly sequence than the previous point-wise
residual scoring (see Section 4.2.1 ). Notably, the AUPRC of _DeepFIB_ in
_ECG(E)_ is slightly lower than RAMED (Shen et al., 2021), and we attribute to
the fact that some unlabeled sub-sequences are too similar to labeled
anomalies in the raw data.

Credit Card: Due to the nature of this application, this dataset is stochastic
and the temporal relation is not significant. Therefore, as shown in Table 4,
traditional AD solutions without modeling the underlying temporal dependency
achieve fair performance, e.g., OCSVM (Ma & Perkins, 2003), ISO Forest (Liu et
al., 2008). Besides, the AR (Rousseeuw & Leroy, 1987) with a small window size
(e.g., 3, 5) can also identify the local change point without considering
longer temporal relations. However, the large _recall_ and small _precision_
values show its high false positive rates. The prediction-based method, LSTM-
RNN (Bontemps et al., 2016) tries to learn a robust temporal relation from the
data, which is infeasible for this dataset. In contrast, the reconstruction-
based method, RAE (recurrent auto-encoder) (Malhotra et al., 2016) performs
better since it can estimate the outliers based on the local contextual
information. The proposed _DeepFIB-p_ model outperforms all baseline methods,
because it can better extract local correlations with the proposed self-
imputation strategy. At the same time, compared to our results on other
datasets, the relative $26.3\%$ improvement over the second best solution
(_AR_) is less impressive and the F1-score with our _DeepFIB-p_ model is still
less than $25\%$. We attribute it to both the dataset complexity and the lack
of temporal relations in this dataset.

Table 4: Comparison of anomaly detection performance (as %), on Credit Card dataset. Methods | Credit Card
---|---
precision | recall | F1-score
AR | 11.30 | 65.20 | _19.20_
ISO Forest | 9.80 | 56.90 | 16.80
OCSVM | 1.70 | _62.00_ | 18.30
LSTM-RNN | 0.40 | 11.00 | 0.70
RAE | 16.90 | 21.52 | 18.89
DeepFIB-p | _16.52 $\pm$ 0.31_ | 46.57 $\pm$ 0.41 | 24.25$\pm$ 0.37
_RAE ∗_ | 13.93 $\pm$ 0.12 | 53.36 $\pm$ 0.21 | 22.07$\pm$ 0.36
_DeepFIB-p †_ | 16.55 $\pm$ 0.22 | 21.08 $\pm$ 0.12 | 18.50$\pm$ 0.21

###  4.2 Q2: Ablation study

In this section, we first evaluate the impact of various components in our
_DeepFIB-s_ and _DeepFIB-p_ models. Next, we replace the SCINet with other
sequence models to evaluate its impact.

####  4.2.1 Component analysis

DeepFIB-p: To demonstrate the impact of the proposed mask-and-impute mechanism
in point outlier detection. We add two baseline methods: (1) _DeepFIB-p †_,
wherein we remove the self-imputation strategy; (2) _RAE ∗_, we implement the
same mask-and-impute strategy and apply it to the baseline method _RAE_. In
Table 4, the performance improvement and degradation of the corresponding
variants compared to _DeepFIB-p_ and _RAE_ clearly demonstrate the
effectiveness of the proposed self-imputation strategy for point outlier
detection.

Figure 8: Ablation study of DeepFIB-s for five datasets (F1-score $\%$). _2nd
best_ denotes the previous SOTA results for each dataset.

DeepFIB-s: To investigate the impact of different modules of _DeepFIB-s_ , we
compare two variants of the _DeepFIB-s_ on five datasets. The details of the
variants are described as below: For _w/o. localization_ , we remove the
anomaly localization algorithm from our _DeepFIB-s_ model. The _w/o.
localization $\&$ DTW_ further removes the DTW scoring mechanism, and the
anomalies are determined based on point-wise residual errors. As shown in Fig.
8, all these components are essential for achieving high anomaly detection
accuracy. At the same time, the proposed self-imputation training strategy is
still the main contributor to the performance of our _DeepFIB-s_ model, as the
results of _w/o. localization $\&$ DTW_ are still much better than those of
the 2nd best solution. Besides, the performance gain of the DTW anomaly
scoring indicates that the point-wise outlier estimation is not suitable for
evaluating sequence-wise anomalies.

####  4.2.2 Impact of SCINet

In our _DeepFIB_ framework, we revise SCINet for time series imputation. To
show its impact, we replace it with other sequence models in _DeepFIB-s_. As
we can see in Table 5, compared with TCN (Bai et al., 2018) and LSTM
(Hochreiter & Schmidhuber, 1997), using SCINet indeed brings significant
improvements on F1-sc, whiorech clearly shows its strong temporal relation
extraction capability and the effectiveness of the revised architecture for TS
imputation. At the same time, compared to the previous SOTA methods (_2nd
best_) for the corresponding dataset, with the same mask-and-impute strategy,
we can still achieve remarkable performance without using SCINet, indicating
the effectiveness of the proposed self-imputation concept itself.

Table 5: The comparison of different sequence models. _2nd best_ denotes the previous SOTA methods in each datasets ( THOC in _2d-gesture_ and RAE-ensemble in _ECG(A)_ ). Methods | ECG(A) | 2d-gesture
---|---|---
SCINet | 80.90 $\pm$ 0.63 | 73.79 $\pm$ 0.19
TCN | 69.86 $\pm$ 0.22 | 69.55 $\pm$ 0.28
LSTM | 64.16 $\pm$ 0.21 | 66.83 $\pm$ 0.55
_2nd best_ | 56.42 | 63.31

##  5 Conclusion

In this paper, we propose a novel self-imputation framework _DeepFIB_ for time
series anomaly detection. Considering the two types of common anomalies in TS
data, we implement two mask-and-impute models biased towards them, which
facilitate extracting more robust temporal relations than existing AD
solutions. Moreover, for sequence-wise anomalies, we propose a novel anomaly
localization algorithm that dramatically improves AD detection accuracy.
Experiments on various real-world TS datasets demonstrate that DeepFIB
outperforms state-of-the-art AD approaches by a large margin, achieving up to
more than $65\%$ relative improvement in F1-score.

## References

  * Angiulli & Pizzuti (2002) F. Angiulli and C. Pizzuti.  Fast outlier detection in high dimensional spaces.  In _PKDD_ , 2002.
  * Bai et al. (2018) S. Bai, J.Z. Kolter, and V. Koltun.  An empirical evaluation of generic convolutional and recurrent networks for sequence modeling.  _ArXiv_ , abs/1803.01271, 2018.
  * Bl’azquez-Garc’ia et al. (2021) A. Bl’azquez-Garc’ia, U. Conde, A.and Mori, and J.A. Lozano.  A review on outlier/anomaly detection in time series data.  _ACM Computing Surveys (CSUR)_ , 54:1 – 33, 2021.
  * Bontemps et al. (2016) L. Bontemps, V. L. Cao, J. McDermott, and N. A. Le-Khac.  Collective anomaly detection based on long short-term memory recurrent neural networks.  _International conference on future data and security engineering_ , abs/1703.09752, 2016.
  * Cao et al. (2018) W. Cao, D. Wang, J. Li, H. Zhou, L. Li, and Y. Li.  Brits: Bidirectional recurrent imputation for time series.  _NeurIPS_ , abs/1805.10572, 2018.
  * Chen et al. (2020) T. Chen, S. Kornblith, M. Norouzi, and G. Hinton.  A simple framework for contrastive learning of visual representations.  _PMLR_ , abs/2002.05709, 2020.
  * Chen et al. (2021) Z. Chen, D. Chen, X. Zhang, Z. Yuan, and X. Cheng.  Learning graph structures with transformer for multivariate time series anomaly detection in iot.  _IEEE Internet of Things Journal_ , abs/2104.03466, 2021.
  * Cook et al. (2020) A. A. Cook, G. Mısırlı, and Z. Fan.  Anomaly detection for iot time-series data: A survey.  _IEEE Internet of Things Journal_ , 7:6481–6494, 2020.
  * Deldari et al. (2021) S. Deldari, D. V. Smith, H. Xue, and F. D. Salim.  Time series change point detection with self-supervised contrastive predictive coding.  _Proceedings of the Web Conference 2021_ , 2021.
  * Deng & Hooi (2021) A. Deng and B. Hooi.  Graph neural network-based anomaly detection in multivariate time series.  _ArXiv_ , abs/2106.06947, 2021.
  * Devlin et al. (2019) J. Devlin, Ming-Wei Chang, Kenton Lee, and Kristina Toutanova.  Bert: Pre-training of deep bidirectional transformers for language understanding.  In _NAACL_ , 2019.
  * Falck et al. (2020) F. Falck, S.K. Sarkar, S. Roy, and S.L. Hyland.  Contrastive representation learning for electroencephalogram classification.  2020\.
  * Fan et al. (2020) H. Fan, F. Zhang, and Y. Gao.  Self-supervised time series representation learning by inter-intra relational reasoning.  _arXiv_ , abs/2011.13548, 2020.
  * Fang & Wang (2020) C. Fang and C. Wang.  Time Series Data Imputation: A Survey on Deep Learning Approaches.  _arXiv_ , 2020.
  * Feng & Tian (2021) C. Feng and P. Tian.  Time series anomaly detection for cyber-physical systems via neural system identification and bayesian filtering.  _Proceedings of the 27th ACM SIGKDD Conference on Knowledge Discovery & Data Mining_, 2021.
  * Gupta et al. (2014) M. Gupta, J. Gao, C.C. Aggarwal, and J. Han.  Outlier Detection for Temporal Data: A Survey.  _IEEE Transactions on Knowledge and Data Engineering_ , 26:2250–2267, 2014.
  * Hochreiter & Schmidhuber (1997) S. Hochreiter and J. Schmidhuber.  Long short-term memory.  _Neural Computation_ , 9:1735–1780, 1997.
  * Holt (2004) C. Holt.  Forecasting seasonals and trends by exponentially weighted moving averages.  _International Journal of Forecasting_ , 20:5–10, 2004\.
  * Hu et al. (2020) W. Hu, B. Liu, J. Gomes, M. Zitnik, P. Liang, V.S. Pande, and J. Leskovec.  Strategies for pre-training graph neural networks.  _ICLR_ , 2020.
  * Hundman et al. (2018) K. Hundman, V. Constantinou, Christopher Laporte, Ian Colwell, and T. Söderström.  Detecting spacecraft anomalies using lstms and nonparametric dynamic thresholding.  _Proceedings of the 24th ACM SIGKDD International Conference on Knowledge Discovery & Data Mining_, 2018.
  * Keogh et al. (2005) E. Keogh, J. Lin, and A. Fu.  Hot sax: Efficiently finding the most unusual time series subsequence.  In _ICDM_ , pp. 226–233, 2005.
  * Kieu et al. (2019) T. Kieu, B. Yang, C. Guo, and C.S. Jensen.  Outlier detection for time series with recurrent autoencoder ensembles.  In _IJCAI_ , 2019.
  * Li et al. (2019) D. Li, D. Chen, L. Shi, B. Jin, J. Goh, and S. Ng.  Mad-gan: Multivariate anomaly detection for time series data with generative adversarial networks.  In _ICANN_ , 2019.
  * Liu et al. (2008) F. Liu, K. Ting, and Z. Zhou.  Isolation forest.  _2008 Eighth IEEE International Conference on Data Mining_ , pp. 413–422, 2008.
  * Liu et al. (2021) M. Liu, Q. Zeng, A.and Lai, and Q Xu.  Time series is a special sequence: Forecasting with sample convolution and interaction.  _ArXiv_ , abs/2106.09305, 2021.
  * Liu et al. (2019) Y. Liu, R. Yu, S. Zheng, E. Zhan, and Y. Yue.  Naomi: Non-autoregressive multiresolution sequence imputation.  In _NeurIPS_ , 2019.
  * Luo et al. (2019) Y. Luo, Y. Zhang, X. Cai, and X. Yuan.  E²gan: End-to-end generative adversarial network for multivariate time series imputation.  In _IJCAI_ , 2019.
  * Ma & Perkins (2003) J. Ma and S. Perkins.  Time-series novelty detection using one-class support vector machines.  _Proceedings of the International Joint Conference on Neural Networks, 2003._ , 3:1741–1745 vol.3, 2003.
  * Malhotra et al. (2016) P. Malhotra, A. Ramakrishnan, G. Anand, L. Vig, P. Agarwal, and G.M. Shroff.  Lstm-based encoder-decoder for multi-sensor anomaly detection.  _ArXiv_ , abs/1607.00148, 2016.
  * Park et al. (2018) D. Park, Y. Hoshi, and C.C. Kemp.  A multimodal anomaly detector for robot-assisted feeding using an lstm-based variational autoencoder.  _IEEE Robotics and Automation Letters_ , 3:1544–1551, 2018\.
  * Pathak et al. (2016) D. Pathak, P. Krähenbühl, J. Donahue, T. Darrell, and A.A. Efros.  Context encoders: Feature learning by inpainting.  In _CVPR_ , 2016.
  * Rousseeuw & Leroy (1987) P. Rousseeuw and A. Leroy.  Robust regression and outlier detection.  1987\.
  * Ruff et al. (2021) L. Ruff, J. Kauffmann, R.A. Vandermeulen, G. Montavon, W. Samek, M. Kloft, T.G. Dietterich, and K. Muller.  A Unifying Review of Deep and Shallow Anomaly Detection.  _Proceedings of the IEEE_ , 109(5):756–795, 2021\.  ISSN 0018-9219.
  * Saeed et al. (2021) A. Saeed, F.D. Salim, T. Ozcelebi, and J.J. Lukkien.  Federated self-supervised learning of multisensor representations for embedded intelligence.  _IEEE Internet of Things Journal_ , 8:1030–1040, 2021.
  * Sakoe & Chiba (1978) H. Sakoe and S. Chiba.  Dynamic programming algorithm optimization for spoken word recognition.  _IEEE Transactions on Acoustics, Speech, and Signal Processing_ , 26:159–165, 1978.
  * Schlegl et al. (2017) T. Schlegl, P. Seeböck, Waldstein, S.M., U.M. Schmidt-Erfurth, and G. Langs.  Unsupervised anomaly detection with generative adversarial networks to guide marker discovery.  In _IPMI_ , 2017.
  * Shen & Kwok (2020) Z. Shen, L.and Li and J.T. Kwok.  Timeseries anomaly detection using temporal hierarchical one-class network.  In _NeurIPS_ , 2020.
  * Shen et al. (2021) Z. Shen, L.and Yu, Q. Ma, and J.T. Kwok.  Time series anomaly detection with multiresolution ensemble decoding.  In _AAAI_ , 2021.
  * Shyu et al. (2003) M. Shyu, Shu-Ching Chen, Kanoksri Sarinnapakorn, and LiWu Chang.  A novel anomaly detection scheme based on principal component classifier.  2003\.
  * Su et al. (2019) Y. Su, Y. Zhao, C. Niu, R. Liu, W. Sun, and D. Pei.  Robust anomaly detection for multivariate time series through stochastic recurrent neural network.  _Proceedings of the 25th ACM SIGKDD International Conference on Knowledge Discovery & Data Mining_, 2019.
  * Yoo et al. (2021) Y. Yoo, U. Kim, and J. Kim.  Recurrent reconstructive network for sequential anomaly detection.  _IEEE Transactions on Cybernetics_ , 51:1704–1715, 2021\.
  * Yu et al. (2016) Q. Yu, L. Jibin, and L. Jiang.  An improved arima-based traffic anomaly detection algorithm for wireless sensor networks.  _International Journal of Distributed Sensor Networks_ , 12, 2016\.
  * Zhang et al. (2019) C. Zhang, D. Song, Y. Chen, X. Feng, C. Lumezanu, W. Cheng, J. Ni, B. Zong, H. Chen, and N.a Chawla.  A deep neural network for unsupervised anomaly detection and diagnosis in multivariate time series data.  _AAAI_ , abs/1811.08055, 2019.
  * Zhang et al. (2016) R. Zhang, P. Isola, and A.A. Efros.  Colorful image colorization.  In _ECCV_ , 2016.
  * Zheng et al. (2018) Y. Zheng, X. Zhou, W. Sheng, Y. Xue, and S. Chen.  Generative adversarial network based telecom fraud detection at the receiving bank.  _Neural networks : the official journal of the International Neural Network Society_ , 102:78–86, 2018.
  * Zhou et al. (2019) B. Zhou, S. Liu, B. Hooi, X. Cheng, and J. Ye.  Beatgan: Anomalous rhythm detection using adversarially generated time series.  In _IJCAI_ , 2019.
  * Zong et al. (2018) B. Zong, Q. Song, M. Min, W. Cheng, C. Lumezanu, D. Cho, and H. Chen.  Deep autoencoding gaussian mixture model for unsupervised anomaly detection.  In _ICLR_ , 2018.