Pythonで時系列解析をする(ホルトウィンターズ)

時系列データはどこの現場にも存在します。
ただし、説明変数が存在しないことも多いです。

その場合、1変量から構造を取り出す手法が存在します。

基本的なパターンはトレンド + 周期性 + ノイズに分割することです。今回利用するホルトウィンターズもそのの方法の一つです。

基本的な使い方

from statsmodels.tsa.holtwinters import ExponentialSmoothing as hw

mdl = hw(
   data,
   trend = 'additive',
   seasonal = 'additive',
   seasonal_periods = 7
)

res = mdl.fit()

res.summary()

トレンドと周期性を使うかとその形、周期の頻度を渡すことで各種パラメータを探してくれます。また、summaryを使うことで各種、値や係数を確認できます。

試してみる

テストデータの作成

以下の特徴をもつデータをつくって試してみます

  • 7日間の曜日の周期性
  • 長い頻度でくるcos波
  • 長期のシグモイドの上昇トレンド
import math
import random
import pandas as pd


trend1 = [ 1/(1 + math.exp(-1*(i - N/2) * 0.1)) for i in  range(N)]
trend2 = [ math.cos(i * 0.25) for i in range(N)]

season_base = [ random.random() for i in range(7)]
season = [season_base[i%7] for i in range(N)]
noise = [random.random() for i in range(100)]

df = pd.DataFrame(
    {
        "trend1":trend1,
        "trend2":trend2,
        "season":season,
        "noise":noise
    }
)


#上記の要素を足し合わせる

ts = df.apply(lambda x : sum(x),axis=1)
ts.plot()

各要素
f:id:esu-ko:20200719113328p:plain 足し合わせ
f:id:esu-ko:20200719111312p:plain

大きいcos波と細かい7日周期にシグモイドがのって徐々に上昇する形になりました。

フィッティングしてみる

mdl = hw(
   ts,
   trend = 'additive',
   seasonal = 'additive',
   seasonal_periods = 7
)

res = mdl.fit()

res.summary()

みるべきところは、smoothing_**のところで、α、β、γが1に近ければ直近の値の影響をうけるモデルになっています。

今回は、αとβが0.4,γが0.2なのでトレンドはやや短期の影響をγはあまり短期の影響を受けていないようです。

各値を実際の値と比較していきます

#トレンド
plt.plot(res.level + res.slope)
plt.plot(df['trend1'] + df['trend2'] )

#季節性
plt.plot(res.season)
plt.plot(df['season'])

#全体
plt.plot(res.fittedvalues)
plt.plot(df['trend1'] + df['trend2'] + df['season'])
ts.plot()

f:id:esu-ko:20200719112431p:plain

トレンドの組み合わせはざっくりした形は捉えているようですが、ノイズを広いギザギザした形が多くなっています。加えて数値が上振れしています。

f:id:esu-ko:20200719112605p:plain そのためか、季節性をみてみると実際の値より低く推定されています。形としては後半は良さそうですが、前半はブレ幅が大きくなっています。

f:id:esu-ko:20200719112956p:plain 全体でみてみると、オレンジのノイズ以外を載せたモデルより、緑のノイズを載せた形にちかくなっています。
基本的な形は抑えれているようですが、今回のようにノイズが季節性と同じくらいの値のデータの場合はそちらに引っ張られるようです。