Categories: Python時系列分析

SARIMAXモデルでPVを予測してみる

北海道電力の推定した時にSARIMAXモデルを使って予測したのですが、SARIMAXモデルでARMAモデルで予測したPVも予測できるんじゃないの?と思ったので試してみました。

■季節調整済みARIMAモデルで電力使用状況を推定してみる
http://jbclub.xii.jp/?p=695

■季節調整済みARIMAモデルで電力使用状況を推定してみる part2 予測表示
http://jbclub.xii.jp/?p=718

%pylab
import pandas as pd
import pandas_datareader.data as pdr
import numpy as np
import matplotlib.pyplot as plt
import statsmodels.api as sm
import seaborn as sns

df = pd.read_csv("C:\\Users\\jbc20170101-20170228.csv",skiprows=1, names=['Date', 'PV']).dropna()
tmd = pd.to_datetime(df['Date'])
df.index = tmd
del df["Date"]
df.plot()

plt.close()

arma_mod0107 = sm.tsa.ARMA(df['PV'], order=(1,7))
arma_res0107 = arma_mod0107.fit(trend='nc',disp=False) 

sarimax_mod17 = sm.tsa.statespace.SARIMAX(df['PV'], trend='n', order=(1,0,7))
sarimax_res17 = sarimax_mod17.fit(trend='nc',disp=False) 

前回、ARMA(1,7)モデルでグラフにしてみたのでSARIMAX(1,0,7)で試してみました。

print (arma_res0107.aic,arma_res0107.bic)
print (sarimax_res17.aic,sarimax_res17.bic)

481.15822431206027 498.719417779
480.311313145 497.872506612

ARMA(1,7)とSARIMAX(1,0,7)のAICとBICを比較してみると、値が異なる結果となりました。同じ値になると予想していたのですが。「trend」のパラメータは「n」「c」「t」「ct」があってそれぞれ試して見たのですが同じ値にはなりませんでした。

sarimax_mod17ct = sm.tsa.statespace.SARIMAX(df['PV'], trend='ct', order=(1,0,7))
sarimax_res17ct = sarimax_mod17ct.fit(trend='nc',disp=False) 
print (sarimax_res17ct.aic,sarimax_res17ct.bic)

472.491043427 493.954724331

今回、「ct」が一番よい結果となったので「ct」でグラフにしてみます。

予測の値を取得します。

# Dynamic predictions
predict_dy = sarimax_res17ct.get_prediction(start='2017-01-01', end='2017/3/15',dynamic='2017/2/21')
predict_dy_ci = predict_dy.conf_int()

グラフにします。

# Graph
fig, ax = plt.subplots()
ax.set(title='JBClub.xii.jp PV', xlabel='Date', ylabel='PV')

# Plot data points
df.ix['2017/1/1':, 'PV'].plot(ax=ax, style='b', label='Value')

# Plot predictions
predict_dy.predicted_mean.ix['2017/2/21':].plot(ax=ax, style='g', label='Dynamic forecast')
ci = predict_dy_ci.ix['2017/2/21':]
ax.fill_between(ci.index, ci.ix[:,0], ci.ix[:,1], color='g', alpha=0.1)

legend = ax.legend(loc='lower right')

トレンドを含め予測された結果となりました。

ところで、パラメータがよくわからないまま使っているところがありました。サンプルを見てそのまま使ってみたのですが、例えば

 arma_mod0107.fit(trend='nc',disp=False) 

の「trend=’nc’」とか「disp=False」とかなんぞや?って感じなんですが、設定しないとエラーになるし、よくわかりませんでした。パラメータについてはもう少し調べてみたいと思います。

※参考
■statsmodels.tsa.arima_model.ARMA
http://statsmodels.sourceforge.net/devel/generated/statsmodels.tsa.arima_model.ARMA.html
■statsmodels.tsa.statespace.sarimax.SARIMAX
http://www.statsmodels.org/dev/generated/statsmodels.tsa.statespace.sarimax.SARIMAX.html

関連する投稿:

kaz

View Comments

Share
Published by
kaz

Recent Posts

よく使うショートカットキー

土曜日の日経新聞のNIKKEI…

3年 ago

気になるETFのメモ

ETFを色々調べていたので銘柄…

4年 ago

SPDRゴールド・シェアETF(GLD)を買いました

無事に米国株取引口座を開設でき…

4年 ago

インフレになったらどうなるのか

最近、インフレになったらどうな…

4年 ago

複利の話とiDeCoのすすめ

転職した会社には確定拠出年金が…

4年 ago