えんじにあのじゆうちょう

勉強したことを中心にアウトプットしていきます。

時系列モデル: ARモデルを実装する

はじめに

最近、なぜか自分の中で時系列モデルがブームになっています。
今回は、時系列データのモデリングの1つであるARモデルについて、勉強したことをまとめてみます。

解説

定義

ARモデルは過去の値から現在の値を求める自己回帰モデルです。
数式としては以下のように表されます。

 y_t = \sum_k \theta_ky_{t-k}

要は線形回帰の各 xがそのまま過去データになったものと解釈できそうです。
つまり、扱う値が違うだけで最小二乗法を解けば求められそうですね。

実装

それではいつものごとく実装してみようと思います。
とはいえ、最小二乗法は以前以下のエントリでやっていたりするので、今回は実装せずに線形回帰のアルゴリズムとしてはscikit-learnのLinearRegressionを使うことにします。

t.marufeuille.dev

まず、このあたりをインポートしておいてください。

import numpy as np
from sklearn.linear_model import LinearRegression
import matplotlib.pyplot as plt

データの生成

まずは回帰に使うデータを作成します。
2つくらい作ってみます。


1. ホワイトノイズ + 異常値1つ

np.random.seed(5)
series1 = np.random.normal(loc=0, scale=1, size=(N,))
series1[66] = 7

f:id:marufeuillex:20191206223714p:plain

2. 基本ホワイトノイズで、あるタイミングを堺に平均が大きく変化 + 異常値1つ

np.random.seed(4)
series_a = np.random.normal(loc=0, scale=1, size=(N//2,))
series_b = np.random.normal(loc=5, scale=1, size=(N//2,))
series_b[6] = 0
series2 = np.hstack([series_a, series_b])

f:id:marufeuillex:20191206223725p:plain

モデル作成の前に

ARモデルは過去何個前までの時点のデータを入れるかにより、どれだけ周期性を捉えられるかが決まります。
ここでは一旦2個前までのデータを使うAR(2)モデルとして構築します。

そのために、データを変形する関数を定義します。

def convert(series, lag=2):
    data = series[:, np.newaxis]
    for j in range(lag):
        data = np.hstack([data, np.roll(series,j)[:, np.newaxis]])
    return data[lag:]

これにより、元々は y_{t}, y_{t+1}, y_{t+2}, y_{t+3}, ...と並んでいる一連のデータを (y_{t}, y_{t-1}, y_{t-2}), (y_{t+1}, y_{t}, y_{t-1}), (y_{t+2}, y_{t+1}, y_{t}), ...と言った形にまとめ直しています。

そして、

data1 = convert(series1)
data2 = convert(series2)

とします。

モデルの作成 ~ 予測

あとは簡単で、モデルをフィットさせるだけです。

lag = 2

# 線形回帰にfitさせる
lr = LinearRegression()
lr.fit(data1[:, 1:], data1[:, 0])

# 直近のデータを使って未来を予測&予測結果の先頭(予測できない部分)に0を詰め、開始番号を合わせる
preds =  lr.predict(data1[:, 1:])
for i in range(10):
    preds = list(preds) + list(lr.predict([preds[-2:]]))
preds = np.hstack([[0 for j in range(lag)],preds])

plt.plot(np.arange(N), series1)
plt.plot(np.arange(N+10), preds, color="orange")
plt.legend()

f:id:marufeuillex:20191206231437p:plain
結果を見ると、ホワイトノイズは期待値=0のただのランダムなので、うまく予測できていないことがわかります。
(ある意味では期待値0付近をうまく予測できていると言えなくもない気がしますが)

series2も同様に予測できます。コードは省きますが、結果は以下のようになります。

f:id:marufeuillex:20191206231449p:plain

とてもきれいにfittingできています。ただ、未来予測はギザギザしているだけで、これでいいのかはちょっと判断がつかないですね。

まとめ

結局いつものことですが、ARモデルでもデータの形や特徴によってはうまくモデリングできなさそうということがわかりました。
始める前は身構えていましたが、基礎から積み上げていくと簡単に理解できるなぁ、としみじみ感じました。