Motivation
단변량 시계열을 모델링하는 방법을 서베이하고 있다. 가장 쉽게 접근할 수 있는 방법이 ARIMA, Holt-winter's method이고, 이 중 holt 방법을 이전 포스트에서 살펴보았다. 여기에 통계적 가정을 추가하여 prediction interval까지 얻을 수 있는 모형으로 ETS가 있으며 이 또한 이 포스트에서 살펴보았다.
하지만 세상에 시계열의 종류는 정말 다양하고.. 그만큼 복잡한 친구들도 많으니.. 내가 회사에서 다루고 있는 시계열 데이터가 간단한 시계열은 아니었고 이에 맞는 방법을 찾기 시작했다.
Why not Holt's method?
모델링하려는 데이터가 air passengers와 같으면 얼마나 좋을까.. 이 데이터는 1) 명확한 주기성이 보이고 2) variance가 커지는 것이 명확히 보이며 3) seasonal period가 작다. 이런 류의 데이터를 모델링하기에는 arima나 holt 방법이 적절할 수 있다. 그러나 현실의 데이터는 생각보다 더러운 경우가 많다.
- noise가 많이 껴있어 주기성이 흐릿하게 보일 수 있다. 이 경우에는 smoothing을 통해 주기성을 명확하게 만들어줘야 한다.
- seasonal period가 매우 클 수 있다. 한달 주기의 시계열은 seasonal period가 12인 반면에 분 주기의 시계열은 seasonal period가 60*24이다 (single seasonality만 고려)
가지고 있는 시계열이 두번째 경우와 같다면 시계열을 upsampling하여 주기성을 줄일 수 있다. 예를 들어, 1분 단위로 관측된 어떤 값을 30분 단위로 평균을 내서 마치 30분마다 관측된 시계열처럼 만드는 것이다. 말은 간단해보이지만 실제로 활용하려면 여러 고민되는 지점이 있다.
- upsampling은 시계열의 발생 빈도를 줄이는 것과 같다. 즉, 자주 발생하는 시계열을 덜 자주 발생하는 시계열로 만드는 것이다. 만약 시계열을 통해 시의성 있는 의사결정을 하고자 한다면, upsampling을 통해 만든 시계열은 활용도가 떨어질 수 있다.
- 목적에 따라 다르겠지만 시계열을 통해 anomaly point을 탐지하고자 한다면, anomaly point가 upsampling을 통해 뭉개질 수 있다.
이런 여러 이유 때문에 시계열의 주기성을 그대로 보존하거나, 크게 줄이고 싶지 않을 수 있다. 문제는 arima나 holt's 방법은 주기성이 높은 시계열을 모델링하는데 적절하지 않다는 것이다.
R의 forecasting 라이브러리 공동 저자인 Hyndman의 글에 따르면 arima나 holt's method는 주기성이 긴 시계열에 적합하지 않다고 한다.
- arima: 이론적으로 주기성이 긴 데이터를 모델링할 수 있다. 그러나 주기성이 긴 데이터의 정상성을 충족하기 위해 차분을 하는 것이 그렇게 추천되는 일은 아니다.
- holt's method: 각 주기성의 initial state을 모두 추정해야한다. 만약, 주기성이 96이라면 추정해야할 파라미터가 96개가 되는 것이다. 이러한 상황에서는 likelihood가 수렴하더라도 표준편차가 매우 크거나 공분산 행렬의 역행렬이 불안정하는 등의 적절하지 못한 결과를 얻게될 수 있다.
이런 motivation을 가지고, 오늘 리뷰하려는 방법론 (BATS, TBATS)가 개발되었다.
What is BATS?
BATS: Box-Cox transform, Arma error, Trend, and Seasonal components의 약자이다.
state space model의 measurement equation과 state equation을 살펴보자.
$$ y^{(\omega)}_t = \begin{cases} \dfrac{y^{\omega}_t - 1}{\omega} & \omega \neq 0 \\ \log y_t & \omega = 0 \end{cases} $$
$$ \begin{align} \ell_t &= \ell_{t-1} + \phi b_{t-1} + \alpha d_t \\ b_t &= (1-\phi)b + \phi b_{t-1} + \beta d_t \\ s^{(i)}_t &= s^{(i)}_{t-m_i} + \gamma_i d_t \\ d_t &= \sum^p_{i=1} \psi_i d_{t-i} + \sum^q_{i=1} \theta_i \epsilon_{t-i} + \epsilon_t \end{align} $$
notation은 아래와 같다.
- $\omega$: box-cox parameter
- $m_1, \cdots, m_T$: $T$개의 seasonal periods
- $\ell_t$: $t$ 시점의 local level
- $b$: long-run trend
- $b_t$: $t$ 시점의 short-run trend
- $s^{(i)}_t$: $t$ 시점의 $i$번째 seasonal component
- $d_t$: $t$ 시점의 ARMA($p,q$) 과정
- $\epsilon_t$: $N(0, \sigma^2)$ 과정
- $\alpha, \beta, \phi, \gamma_i$: smoothing parameter
BATS가 풀고자 했던 문제
비선형 state space model이 가지는 문제
exponential smoothing에 근간을 두는 비선형 state space model은 그렇게 좋은 성능을 내지 못했다. 불안정한 추정치를 보였고 일정 forecasting horizon을 넘으면 infinite variance가 나오기도 했다. 이를 해결하기 위해 BATS는 Box-Cox 변환을 이용하여 비선형 데이터까지 모델링할 수 있는 여지를 남겼다.
measurement equation의 error가 uncorrelated 하다는 가정
exponential smoothing의 state space model은 기본적으로 에러항이 uncorrelated하다고 가정한다 ($\epsilon \sim NID(0, \sigma^2)$ 하지만, 현실 데이터에서는 이렇지 않은 경우가 더 많다. 따라서 어느정도 correlation을 허용하는 것이 현실적인 가정이다. 이를 해결하기 위해 BATS는 에러항이 ARMA($p,q$)을 따른다고 가정한다.
BATS가 여전히 풀지 못하는 문제
seasonal periods가 커질수록 파라미터 수가 늘어나는 문제
기본적으로 BATS의 계절성 부분은 holt's method와 동일하기 때문에 이 문제에서 자유로울 수는 없다. 각 계절성 부분의 seasonal period를 $m_i$라고 한다면, 계절성 부분에서만 총 $\sum^T_{i=1} m_i$개의 initial state estimate가 필요하다.
non-integer seasonal period을 가지는 시계열을 모델링할 수 없는 문제
마찬가지로 BATS의 계절성 부분은 holt's method와 동일하기 때문에 이 문제에서 자유로울 수 없다.
BATS가 가지는 이런 문제들을 해결하기 위해 TBATS가 개발되었다.
What is TBATS?
Trigonometric BATS. BATS의 특징을 그대로 가져가면서 seasonal component만 Fourier series을 이용해서 표현하는 방법이다.
TBATS의 state equation은 아래와 같다.
$$ y^{(\omega)}_t = \begin{cases} \dfrac{y^{\omega}_t - 1}{\omega} & \omega \neq 0 \\ \log y_t & \omega = 0 \end{cases} $$
$$ \begin{align} \ell_t &= \ell_{t-1} + \phi b_{t-1} + \alpha d_t \\ b_t &= (1-\phi)b + \phi b_{t-1} + \beta d_t \\ s^{(i)}_t &= \sum^{k_i}_{j=1} s^{(i)}_{j,t} \\ d_t &= \sum^p_{i=1} \psi_i d_{t-i} + \sum^q_{i=1} \theta_i \epsilon_{t-i} + \epsilon_t \end{align} $$
여기서 seasonal component $s^{(i)}_t$을 구성하는 $s^{(i)}_{j,t}$는 아래와 같이 fourier series에 기반하여 정의된다.
$$ \begin{align} s^{(i)}_t &= \sum^{k_i}_{j=1} s^{(i)}_{j,t} \\ s^{(i)}_{j,t} &= s^{(i)}_{j,t-1} cos \lambda^{(i)}_j + s^{*(i)}_{j,t-1} sin \lambda^{(i)}_j + \gamma^{(i)}_1 d_t \\ s^{*(i)}_{j,t} &= s^{(i)}_{j,t-1} sin \lambda^{(i)}_j + s^{*(i)}_{j,t-1} cos \lambda^{(i)}_j + \gamma^{(i)}_2 d_t \end{align}$$
notation은 아래와 같의 정의된다.
- $\gamma^{(i)}_1, \gamma^{(i)}_2$: smoothing parameters
- $\lambda^{(i)}_j = 2 \pi j /m_i$
- $s^{(i)}_{j,t}$: stochastic level of the ith seasonal component
- $s^{*(i)}_{j,t}$: stochastic growth in the level of the ith seasonal component
- $ k_i = \begin{cases} m_i / 2 & m_i \text{ is even} \\ (m_i - 1)/2 & m_i \text{ is odd} \end{cases} $: number of harmonics required for the ith seasonal component
위 식을 종합하여 TBATS의 measurement equation을 도출하면 아래와 같다.
$$ \begin{align} y^{\omega}_t = \ell_{t-1} + \phi b_{t-1} + \sum^T_{i=1} s^{(i)}_{t-1} + d_t \end{align}$$
TBATS가 풀고자 했던 문제
추정해야하는 파라미터 수를 줄였다.
BATS 모형까지는 seasonal periods가 커지면 그만큼의 initial states을 추정해야했기 때문에 추정해야할 파라미터 수가 많았다. TBATS는 trigonometric equation을 통해 $2(k_1 + k_2 + \cdots + k_T)$개 만큼의 initial seasonal values의 추정을 필요로 한다. 이는 결국 ith seasonal component의 필요한 harmonics 수인데.. 논문에서는 이거를 추정함으로써 파라미터수가 줄어든다고 했지만 아직 이에 대해 정확히 이해를 하지는 못했다.
non-integer seasonal periods을 가지는 시계열도 모델링할 수 있다.
trigonometric equation으로 계절성 부분을 모델링하기 때문에 non-integer인 시계열도 모델링할 수 있다.
Using tbats model in python
sktime 라이브러리에서 tbats의 api를 제공한다. 정확히는, sktime도 tbats을 직접 구현한 것은 아니고 다른 레포에서 구현한 tbats을 api로 가져와 쓰는 형식이다. sktime 공식 문서에 있는 예시를 그대로 가져와봤다.
from sktime.datasets import load_airline
from sktime.forecasting.tbats import TBATS
y = load_airline()
forecaster = TBATS(
use_box_cox=False,
use_trend=False,
use_damped_trend=False,
sp=12,
use_arma_errors=False,
n_jobs=1)
forecaster.fit(y)
y_pred = forecaster.predict(fh=[1,2,3])
print(y_pred)
```
1961-01 442.565236
1961-02 431.995779
1961-03 466.500825
```box-cox, trend, arma error에 대한 옵션을 지정하는 argument가 따로 있다. 사실, 나는 tbats가 적합한 파라미터를 보고 싶었는데.. get_fitted_params 메써드를 사용해도 결과가 나오지 않아서 당황했다. 코드를 살펴본 결과.. sktime의 메써드에서 적합한 파라미터를 보기보다는 가져온 구현체의 메써드를 이용해서 적합된 파라미터를 살펴봐야하는 것 같았다.
print(forecaster._forecaster.summary())
```
Use Box-Cox: False
Use trend: False
Use damped trend: False
Seasonal periods: [12.]
Seasonal harmonics [5]
ARMA errors (p, q): (0, 0)
Smoothing (Alpha): 1.186745
Seasonal Parameters (Gamma): [-2.90794894e-07 1.66634461e-08]
AR coefficients []
MA coefficients []
Seed vector [135.68553355 -45.90885595 19.00442987 -4.14834791 3.62942987
1.98716013 3.89859135 15.49990498 -8.88445902 -6.69310179
-5.81142743]
AIC 1554.897669
```내가 잘못 본건지.. 아직 sktime의 api가 완성되지 않은 부분이 있는 것인지는.. 모르겠지만 우선은 이렇게 사용해야할 것 같다.
tbats는 state space model이므로 likelihood을 이용해서 모형을 적합하고 따라서 prediction interval을 구할 수 있다. sktime에서는 predict_interval 메써드를 통해 이를 수행할 수 있다.
forecaster.predict_interval(fh=[1,2,3])
interval 까지 제공을 하니 anomaly detection 등에 유용하게 쓰일 수 있을 것 같다~
단점
시간이 많이 걸린다.. 꽤나 많이 걸린다.. holt-winter's method에 비해 확실히 많이 걸린다. 이 medium 글에는, 시간을 단축학 위한 방법으로 box-cox, arma error의 옵션을 False로 두는 것을 제안한다. 음.. 우선 코드적인 부분에서 시간을 줄일 수 있는지 생각해보고, 그래도 불가능하다면 모형의 가정을 하나씩 푸는 것도 괜찮을 것 같다.
Reference
'ML&DL > Times Series' 카테고리의 다른 글
| [ML / Time series] State space models for exponential smoothing (0) | 2023.03.07 |
|---|---|
| [ML / Time series] Exponential Smoothing (0) | 2023.01.02 |
댓글