본문 바로가기

Gaussian Process

Adaptive Gaussian Process Change Point Detection(ADAGA): 일변량 시계열 데이터 이상치 탐지 방법론


Time seris data에서 change point를 탐지(detect)하는 것은 현실 세계에서 많은 분야의 근본적인 업무입니다. Change Point(CP)란, 데이터를 표현하는데 사용된 모델의 hyper-parameter가 갑자기 변하는 지점(point)을 의미합니다.

본 글은 논문 Adaptive Gaussian Process Change Point Detection(2022)을 바탕으로 gaussian process(가우시안 프로세스)와 statistical hypothesis testing(가설 검정)을 결합한 ADAGA를 통해 time series data에서 어떻게 CP를 탐지하는지에 대해서 설명하고자 합니다.

 

ADAGA의 경우, 크게 두 개의 파트로 나눌 수 있으며 자세한 설명은 <Section 2>에서 하겠습니다.

< " Key point: Adaptive Gaussian Process Change Point Detection(ADAGA) "

: GP-based CP detection algorithm, based on statistical testing > 

[1st Part] Gaussian Process-based solution

  • temporal prcoess의 mean(평균)과 covariance structure(공분산 구조)가 함께 변화합니다.
    • 즉, 평균 또는 공분산 구조가 상수이어야 한다는 가정에 의존하는 것 없이 데이터로부터 GP hyper-parameter를 추론합니다.

[2nd Part] statistical hypothesis testing basd on a generalized likelihood ratio test

  • 데이터로부터 추론된 GP hyper-parameter를 기반으로 현재 sliding window 상에 CP가 존재하는지 여부에 대해 가설검정을 수행합니다.

0. Introduction

<Section 0>에서는 'CP detection in time series'에 대한 기존의 방법론들의 간략한 소개 및 한계점을 토대로 ADAGA의 등장 배경에 대해 설명하고자 합니다.

 

 Bayesian message-passing alogrithm(BOCPD) + Gaussian Process(GP)

CP detection의 기본적인 solution은 Adams&MacKay(2007)가 제안한 방법론인 BOCPD입니디다. BOCPD는 subsequent CPs 사이의 시간인 run-length에 의존하는 방법론이며 CP들끼리 independent identically distribution(i.i.d) 샘플이라는 가정을 합니다. 하지만, 대부분 time series data는 서로 temporal correlation이 존재하기 때문에 i.i.d 가정이 성립하지 않습니다. 이에 대해서 Saatci et al.(2010)은 BOCPD와 Gaussian Process(GP)의 결합을 제안합니다. GP의 경우, 모델을 직접으로 time domain에서 signal domain으로 매핑시키기 때문에 flexible하고 powerful한 모델입니다. 

 

▶ Statistical Hypothesis Testing + Gaussian Process(GP)

[BOCPD+GP의 단점]

GP와 결합한 BOCPD-based detection은 flexible하지만 prior hyper-parameter에 극심하게 의존한다는 단점이 존재합니다. 이는 성능 저하로 이어집니다.

 

이에 대하여 Han et al(2019)는 error threshold를 제공하는 statistical hypothesis testing(가설검정)이 CP detection에 있어서 더 robust함을 주장합니다. Test based approach는 deep learning methods 또는 Gaussian processes(GP)와 결합하여 사용될 수 있습니다.

GP와 결합한 test based approach의 경우, GP의  mean function(평균 함수) 혹은 convariance function(공분산 함수)가 변하는 지점 혹은 시점을 가설검정을 토대로 CP 여부 판단합니다. 기존의 generalized likelihood ratio test 기반인 가설검정과 GP와 결합한 방법론은 평균 함수 혹은 공분산 함수 둘 중의 하나를 상수로 고정하고 가설검정을 통해 CP를 탐지합니다.

  • Keshavarz et al.(2018) 고정된 covariance function(→ 상수로 설정)인 GP의 mean function 변화를 기반으로 가설검정을 통해 single CP를 탐지합니다.
  • Han et al.(2019)는 정된 mean function (→ 상수로 설정) 인 GP의 covariance function 변화를 기반으로 가설검정을 통해 single CP를 탐지합니다.

▶ ADAGA: Statistical Hypothesis Testing + Gaussian Process(GP)

[기존에 존재하는 방법론의 Statistical Hypothesis Testing+GP의 단점]

Mean function(평균 함수) 혹은 covariance function(공분산 함수)을 상수로 고정시키는 경우, CP detection algorithm이 찾을 수 있는 비정상(anomalies)의 개수가 감소할 수 있습니다. 

 

논문에서 제안하는 ADAGA의 모델의 경우 GP-baesd CP detection 알고리즘이며 기존의 방법론과는 다르게 mean function 혹은 covariance function을 상수로 고정시키지 않습니다. 즉, 평균 함수와 공분산 함수 둘 중의 하나가 상수라는 가정을 두지 않고 데이터로부터 추론된 GP hyper-parameter를 기반으로 CP를 탐지하는 가설검정을 수행합니다. 

 

▶ Summary

요약하면, CP detection algorithm의 기본적인 방법론인 "GP based BOCPD dection algorithm"은 prior hyper-parameter에 크게 의존하여 'GP based statistical hypothesis testing'이 제안되었습니다. 기존에 존재하는 'GP based statistical hypothesis tesing' 방법론은 mean function 혹은 covariance function을 상수로 고정시켰습니다. 이는 탐지하는 비정상 개수가 감소할 수 있다는 문제점이 존재하여 평균 함수와 공분산 함수 둘 중의 하나가 상수라는 가정을 두지 않는 ADAGA 모델을 제안하였습니다.


1. Background: GP regression

<Section 1>에서는 Gaussian Process Regression의 간단한 소개를 하고자 합니다. 자세한 내용은 다음 아래의 블로그를 참조하시기를 바랍니다.

▶ [Notation]

  • Given Dataset: $\{x_i,y_i\}^n_{i=1}, \, x_i \in \mathbb{R}^{p}$
  • $\mathbf{x}=(x_1,x_2,\cdot,x_n)^T, \, \mathbf{y}=(y_1,y_2,\cdots,y_n)^T$
  • $\mathbf{y}=f(\mathbf{x})+\epsilon_n, \, \epsilon_n \sim N_n(0, \sigma^2)$
    • 함수 $f: \mathbb{R}^p \to \mathbb{R}$

▶ [Prior Distribution: Gaussian Process or Multivariate Gaussian Distribution]

Gaussian Process는 함수에 대한 분포입니다.

$$f(x) \sim GP(m(x), k(x,x'))$$

  • $m(x)$: mean function(평균함수)
  • $k(x,x') = E[(f(x)-m(x))(f(x')-m(x'))]$: covariance function
    • GPR의 경우, 커널 함수 $k(\cdot, \cdot)$에 의해 함수 $f$의 형태(shape)가 달라집니다([Figure 1 참조).

[Figure 1] Prior Model

보편적으로 평균 함수 $m(x)=0$으로 설정합니다. 만약, $f$가 평균함수 $m(x)$가 $0$인 GP를 따르는 경우, $n$차원의 함수 $f(x) = \begin{bmatrix} f(x_1) \\ f(x_2) \\ \vdots \\ f(x_n) \end{bmatrix} \in \mathbb{R}^n$은 $n$차원의 다변량 정규분포를 따릅니다.

$$\begin{bmatrix} f(x_1) \\ f(x_2) \\ \vdots \\ f(x_n) \end{bmatrix}  \sim N_n(0, \begin{bmatrix} k(x_1, x_1) + \sigma^2 & k(x_1,x_2) & \cdots & k(x_1,x_n) \\ k(x_2,x_1) & k(x_2,x_2) +  \sigma^2 & \cdots & k(x_2,x_n) \\  \vdots& \dots & \vdots  \\ k(x_n, x_1) & k(x_n, x_2) & \cdots & k(x_n,x_n) +  \sigma^2 \end{bmatrix})$$

 

▶ [Update: Posterior]

만약, 새로운 데이터 $x_*$이 입력된 경우, 관측된 데이터 $f=(f(x_1),f(x_2),\cdots, f(x_n))$와 새로운 데이터 $f_* = f(x_*)$의 결합정규분포를 정의할 수 있습니다.

$$\begin{bmatrix} f(x_1) \\ f(x_2) \\ \vdots \\ f(x_n) \\ f(x_*) \end{bmatrix}  \sim N(0, \begin{bmatrix} k(x_1, x_1) + \sigma^2 & k(x_1,x_2) & \cdots & k(x_1,x_n) & k(x_1,x_*) \\ k(x_2,x_1) & k(x_2,x_2) + \sigma^2 & \cdots & k(x_2,x_n)  & k(x_2, x_*) \\  \vdots& \dots & \vdots  \\ k(x_n, x_1) & k(x_n, x_2) & \cdots & k(x_n,x_n)  + \sigma^2 & k(x_n, x_*)\\ k(x_*,x_1) & k(x_*,x_2) & \cdots & k(x_*, x_n) & k(x_*, x_*) \end{bmatrix}) \tag{1}$$

식 $(1)$을 더 간단하게 다음과 같이 표현 가능합니다($f(\mathbf{x}) = (f(x_1),f(x_2),\cdots,f(x_n))^T$).

$$\begin{bmatrix} f(\mathbf{x}) \\ f(x_*) \end{bmatrix}  \sim N(0, \begin{bmatrix} K(\mathbf{x},\mathbf{x}) + \sigma^2I_n& K(\mathbf{x},x_*) \\ K(x_*,\mathbf{x}) & k(x_*, x_*) \end{bmatrix}) \tag{2}$$

우리의 주 목적은 기존의 관측치 $\mathbf{x}$에 더해 새로운 관측치 $x_*$이 주어졌을 때 기존의 함수(사전 분포)  $f$를 새로운 함수(사후 분포) $f_*$로 업데이트 하는 것입니다.

$(2)$ 조건부 분포로 변환하면 다음과 같습니다(사후 분포).

$$f_*|x_*,\mathbf{x},f \sim N(\mu_*, \sigma^2_*)$$

  • $\mu_* = K(x_*,\mathbf{x})(K(\mathbf{x},\mathbf{x})+\sigma^2I_n)^{-1}f$
  • $\sigma_* = k(x_*,x_*) - K(x_*,\mathbf{x})(K(\mathbf{x},\mathbf{x})+\sigma^2I_n)^{-1}K(\mathbf{x},x_*)$

[FIgure 2] Gaussian Process Regression

[Figure 2]빨간색 선(red line)이 사후 평균 $\mu_*$, 하늘색 음영(blue area)이 사후 표준편차 $\sigma_*$입니다.

GPR의 가장 큰 특징은 불확실성(uncertainty, $\sigma_*$)을 다룰 수 있다는 점입니다. 관측된 데이터 ×가 존재하는 부분은 상대적으로 불확실성이 적은 반면, 관측된 데이터 ×가 존재하지 않은 부분은 불확실성이 큰 것을 확인할 수 있습니다.

 

GPR의 경우 불확실성을 다룰 수 있다는 큰 특징이 존재하지만 커널 계산 비용에 있어서 큰 비용이 든다는 단점이 존재합니다. 이에 대한 해결방안으로 'inducing points(Snelson & Ghahramani, 2006; Titsias, 2009; Hensman et al., 2013)''feature approximations'가 approximation 방법론이 제안되었습니다.


2. Adaptive Streaming GP Regression

<Section 2>에서는 GP-based change point detection 방법론인 ADAGA에 대해 본격적으로 설명하고자 합니다. 

2.1 CPs and Kernel Function

▶ CP is strictly related to the kernel function 

GP based CP detection에 있어서 CP는 커널 함수(kernel function)와 밀접하게 연관되어 있다는 중요한 성질이 존재합니다.

<Section 1>에서 보았듯이 [Figure 1]을 통해 커널 함수에 따라 사전 분포(prior distribution)로부터 샘플링된 함수 $f$의 형태(shape)가 달라지는 것을 확인할 수 있습니다. 이는 GP based algorithm에 의해 탐지(detect)되는 CP의 type은 어떤 prior model(혹은 어떤 커널 함수)를 사용하는지에 따라 결정된다는 것을 의미합니다.

 

2.2 ADAGA's partitioning strategy and Hypothesis

< Section 2.2>에서는 window 및 subwindow에 대해 살펴보고자 합니다.

■ Windows and Subwindows

ADAGA는 streaming data에서 작동하도록 디자인된 모델입니다. 즉, 무한개의 데이터가 존재하는 상황이라고 생각하시면 될 거 같습니다.

※ Streaming data는 연속적으로 생성되고 전송되는 데이터를 나타내는 용어로 일정한 흐름으로 제공되며 일반적으로 실시간으로 처리됩니다. 

※ Streaming data(시간의 순서가 존재): $\forall (x_h, y_h), (x_k, y_k) \in D, h < k$ → $x_h < x_k$

 

▶ Notation

  • sequence of datapoints $(x_i, y_i), \, x_i \in \mathbb{R}^d, \, y_i \in \mathbb{R}$: infinite amount of time으로부터 샘플링된 객체
  • Make a parition of dataset $D$ = $\{x_i,y_i\}^M_{i=1}$: 전체 데이터셋의 partition dataset 생성
    • 샘플링할 데이터의 수 ${}^{\exists}M, \, M \in \{1,2, \cdots, n\}$는 임의로 설정합니다.

▶ Assumption

샘플링하는 과정에 있어서 각 데이터 포인트들의 위치(location)에 대해서는 어떠한 가정도 하지 않습니다. 즉, 데이터 포인트간의 공간적 배열에 대한 가정을 하지 않고 무작위로 샘플링을 하겠다는 것을 의미합니다.

 

▶ Definition of Window

[Figure 3] Window

ADAGA는 계속해서 새로운 데이터 점들을 샘플링함으로써 데이터셋 $D$의 partition dataset을 업데이트합니다. 우리는 이 partition dataset에 속하는 부분집합(subset) $\mathbf{\color{lightblue}{W_1, W_2, \cdots, W_p}}$를 window라고 명시합니다. ADAGA에 의해 생성된 window $W_i, \, i \in \{1,2, \cdots, p\}$에는 샘플링된 순서에 따라 인접한(adjacent) 데이터 포인트가 포함되어 있습니다. Window를 수식으로 정의하면 다음과 같습니다.

$$\forall{i} \in \{1,2,\cdots,p\}, {}^{\exists} m, M \in \{1,2,\cdots,n\}, s.t.\\ W_i = \{(x_w, y_w)\in D | m \leq w \leq M\}$$

 

Partition dataset은 계속해서 증가하는 데이터셋으로 현재 window $W_i$에 데이터들이 계속해서 추가되어집니다. 즉, window size $|W|$는 고정된 크기가 아닙니다($|W_1| \ne |W_2| \ne \cdots |W_p|$).

Question: 그럼 언제까지 현재 window $W_i$에 데이터셋들( batch)이 추가되는가?

Answer: 현재 window $W_i$에서 CP가 포함된 것을 발견할 때까지입니다.

이와 같이 window안에 CP가 포함된 경우 window는 spoiled되었다고 말합니다.

※ [요약] The partition dataset is created dynamically. During Streaming, the batches are appended to the same window, until it is found to contain a CP

 

▶ Window가 spoiled된 경우

현재 window $W_i$가 spoiled된 경우,  가장 최근에 추가된 데이터를 포함한 새로운 window $W_{i+1}$를 생성합니다. 

 

Question: 그렇다면 새로운 batch가 추가될 때마다 현재 window $W_i$가 spoiled되었는지 어떻게 판단하는가?

Answer: window의 부분집합인 subwindow $S$를 통해서 판단합니다

 

▶ Criterion: Subwindow $S$

Subwindow $S$는 가장 최근에 추가된 데이터들인 window의 마지막 부분을 의미합니다.

※ 시간의 순서가 존재하는 streaming data의 성질로 window 내에 존재하는 점들과 subwindow에 존재하는 점들은 모두 함수 도메인 안에서 서로 인접해 있습니다.

[Figure 4] Subwindow

Window의 크기 $|W|$는 고정되어있지 않은 반면, subwindow 크기 $|S|$는 고정되어 있다는 특징이 존재하며 hyper-parameter입니다. 우리는 데이터의 분포에 대한 사전적인 지식을 바탕으로 $|S|$를 설정해야 하며 만약 잘못된 $|S|$를 사용할 경우, CP detection에 큰 영향을 끼칩니다. 각 streaming step마다(새로운 batch가 추가될 때마다) subwindow의 크기는 다음 조건을 만족해야 합니다.

$$|W_i| \geq |S| \tag{3}$$

[식 $(3)$의 효과]

식 $(3)$은 ADAGA가 너무 작은 window 내에서 작동하지 않도록 해주며 overfitting 방지 효과도 존재합니다.

 

만약, 전체 window 상에서의 global GP의 regression보다 subwindow 상에서의 local GP의 regression이 향상된다면 window는 spoiled되었다고 판단합니다. 즉, 전체 window 상에서 GP를 훈련했을 때(global GP)와 오직 subwindow 상에서 GP를 훈련했을 때(local GP)의 효과를 비교해 보는 것입니다.

[Figure 4] Hypothesis

▶ Null Hypothesis

    "The null hypothesis assumes that the function values in the subwindow
come from the same observational model as the rest of the window"

 

귀무가설의 경우, subwindow 상에서의 함수 값이 나머지 window 상에서의 관측 모형(observational model)과 동일한 관측 모형에서 나온다고 가정합니다. 즉, 현재 subwindow 상에는 CP가 없기 때문에 window는 spoiled 되어있지 않다고 주장하는 가설입니다.


"Definition 1"

  • $\phi_{H_0}$: whole currrent window로부터 학습한 GP hyperparameters
  • $\sigma^2_{H_0} $: whole currrent window로부터 학습한 noise variance
  • $H_0$하에서, subwindow 상에서 관측치 $y_s$의 marginal likelihood는 다음과 같습니다.

$$p(y_s|H_0)=p(y_s|x_s, \phi_{H_0},\sigma^2_{H_0}) \tag{4}$$


▶ Alternative Hypothesis

"The alternative hypothesis assumes that the function values in the subwindow
come from the superposition of two GP experts:"

대립가설은 귀무가설과는 반대로 현재 window가 spoiled되었다고 주장하며 2개의 GP 모델(GP expert)을 결합하여 사용합니다("a combination of a mixure-of-experts").

 

$\color{darkred}{\mathbf{M_c}}$: potentially spoiled, which has the same parameters as the rest of the window

$M_c: \, p(y_s|x_s, \phi_{H_0}, \sigma^2_{H_0}) \tag{5}$

$\color{darkred}{\mathbf{M_n}}$: potentially much better whose parameters are inferred from the subwindow only

$M_n: \, p(y_s|x_s, \phi_{new}, \sigma^2_{new}) \tag{6}$

첫번째 GP expert $M_c$는 subwindow 상에서와 나머지 window 상에서 같은 parameter를 가지며 subwindow 구간은 잠재적으로 spoiled되어 있다고 가정합니다. [Figure 4-b]를 통해 subwindow 구간에서 전체 window로부터 학습한 GP hyperparameter $\phi_{H_0}, \sigma^2_{H_0}$으로부터 구해진 posterior distribution($M_c$의 사후 분포)의 평균 함수(red line) $\color{red}{u_*}$는 실제 데이터 곡선(green line)으로부터 벗어나는 것을 확인할 수 있습니다.

 

두번째 GP expert $M_n$는 오직 subwindow로부터 추론된 parameter $\phi_{new}, \sigma^2_{new}$가 subwindow 상에서 데이터를 훨씬 더 잘 설명한다고 주장합니다. [Figure 4-c]를 통해 subwindow 구간에서 subwindow 상에서 학습한 GP hyperpameter $ \phi_{new}, \sigma^2_{new} $로부터 구해진 posterior distribution의 평균 함수(red line) $\color{red}{u_*}$는 실제 데이터 곡선(green line)과 정확하게 일치하는 것을 확인할 수 있습니다.

 

두 개의 GP model $M_c$와 $M_n$의 likelihood 곱(식 $(5)$와 식 $(6)$의 곱)을 통해 $H_1$하에서 subwindow 상에서 관측치 $y_s$의 marginal likelihood를 구할 수 있습니다.

$$p(y_s|H_1) = \dfrac{p(y_s|x_s, \phi_{H_0}, \sigma^2_{H_0})p(y_s|x_s,\phi_{new}\sigma^2_{new})}{Z_1} \tag{7}$$

  • $Z_1$: normalization term

$H_1$하에서는 결국 전체 데이터셋을 두 개의 GP expert의 곱을 통해 표현합니다.

※ [참고] 두번째 GP expert $M_n$의 경우는 subwindow 구간에 존재하는 데이터로부터 추론되기 때문에 데이터의 first part(초기 부분)에서는 오직 첫번째 GP $M_c$로만 모델링합니다.

 

▶ Likelihood Ratio Test

우리는 현재 subwindow 상에 존재하는 객체 $y_s$에 대해서 두 개의 서로 다른 liklihood model $p(y_s|H_0)$와 $p(y_s|H_1)$가 존재합니다. 따라서, 우리는 likelihood ratio test로 가설검정을 할 수 있습니다.


"Definition 2"

Given null and alternative hypothesis, the likelihood ratio is

$$R = 2 ln \dfrac{p(y_s|H_1)}{p(y_s|H_0)}=2 ln \dfrac{1}{Z_1} - 2ln Z_{new} - y^T_sV^{-1}_{new}y_s  \tag{8}$$

  • $Z_{new}$: he normalization constant for the GP expert whose hyper-parameters are learned from the subwindow only
  • $V_{new}$: its covariance matrix

임의의 threshold $\tau$가 존재하는 경우,

  • $R \geq \tau$이면 귀무가설 $H_0$를 기각합니다. 즉, 현재 window는 spoiled되었다고 판단합니다.
  • $R < \tau$이면 귀무가설 $H_0$을 채택합니다. 즉, 현재 window는 spoiled되지 않았다고 판단합니다.

만약 현재 window가 spoiled되었다고 판단된 경우 subwindow $|S|$의 마지막 데이터들을 제외한 나머지 데이터 $(x_{window}, y_{window})$를 저장하고 새로운 window를 생성합니다.

  • $x_{window}$ = Last $|S|$ points in $x_{window}$
  • $y_{window}$ = Last $|S|$ points in $y_{window}$

▶ Algorithm of ADAGA

전체적인 알고리즘은 아래와 같습니다.


본 글은 GP based statistical hypothesis test 방법론 중 ADAGA에 대해 설명한 글입니다.ADAGA가 window와 subwindow를 생성하는 방법부터 시작해서 subwindow 구간에 존재하는 $y_s$에 대해 2개의 서로 다른 likelihood $p(y_s|H_0)$와 $p(y_s|H_1)$를 정의함으로써 likelihood ratio test를 수행하였습니다. 임의의 treshold $\tau$를 설정하는 부분은 따로 다루지 않았는데 궁금하신 분은 논문을 참고하시면 될 거 같습니다.