Wikipediaアクセス数の推移を予測してみた

目次
1.はじめに
　1.1 検証内容
　1.2 実装環境

2.データ前処理

3.学習モデルの作成
　3.1 SARIMA
　3.1.1 パラメーターの決定
　3.1.2 モデルの構築と学習・予測
3.2 LSTM
　3.2.1 説明変数と目的変数の生成
　3.2.2 モデルの構築
　3.2.3 モデルの学習
　3.2.4 モデルの予測

4.比較
　4.1 データの編集
　 4.2 グラフで可視化
　4.3 各モデルの精度を数値化

5.まとめ

1.はじめに

1.1 検証内容

元アメリカンフットボール選手のペイトン・マニングのウィキペディアのページへのアクセス数の推移を時系列解析で予測してみました。
今回は、SARIMAモデルとLSTMモデルを使って、データを予測してみました。

1.2 実装環境

・Google Colaboratory
・Python

2.データ前処理

DataFrameでcsvを読み込んだ後、indexを日ごとにして時系列データを作成したいので、dsの列をindexに指定します。

# アメリカンフットボールプレーヤのPayton ManningのWikiアクセス数
df = pd.read_csv('http://logopt.com/data/peyton_manning.csv')
df = pd.DataFrame(df)
df = df.set_index("ds")

3.学習モデルの作成

3.1　SARIMAモデル

SARIMAモデルで、アクセス数の時系列解析を行う。

3.1.1 パラメーターの決定

まずはSARIMAモデルの学習に使う最適なパラメーターを決定する。
今回は、BIC(ベイズ情報量基準)によって最適なパラメーターを決定していく。
コードは以下のコードで、最適なパラメーターを求める。

def selectparameters(DATA,s):
    p=d=q=range(0,2)
    pdq = list(itertools.product(p,d,q))
    seasonal_pdq =[(x[0],x[1],x[2],s)  for x in pdq]
    parameters =[]
    BICs = np.array([])
    for param in pdq:
      for param_seasonal in seasonal_pdq:
        try:
          mod = sm.tsa.statespace.SARIMAX(DATA,order=param, seasonal_order=param_seasonal)
          results = mod.fit()
          parameters.append([param,param_seasonal,results.bic])
          BICs = np.append(BICs,results.bic)
        except:
          continue

    return parameters[np.argmin(BICs)]

3.1.2　モデルの構築と学習・予測

続いて、決定したパラメーターを用いてモデルを構築し、そのモデルで学習・予測を行っていく。

SARIMA_model = sm.tsa.statespace.SARIMAX(train, order = (1,1,1), seasonal_order = (0,0,1,12)).fit()
pred_SARIMA = SARIMA_model.predict("2015-10-23","2016-01-20")


fig = plt.figure(figsize = (20,8))
plt.title("Payton Manning wiki access count")
plt.xlabel("time")
plt.ylabel("access count")
plt.plot(pred_SARIMA)
plt.show()

3.2 LSTM

次は、LSTMで時系列解析を行っていく。

3.2.1 説明変数と目的変数の生成

LSTMで時系列解析を行っていくために、データから説明変数目と的変数を生成する。
LSTMはニューラルネットワークを使った手法のため、入力層に説明変数を入力し、出力層で出力されるデータと目的変数との誤差に着目し、その誤差が小さくなるように学習していく手法のため、説明変数と目的変数を生成する必要がある。

# 目的変数と説明変数の作成
def make_dataset(low_data, maxlen):

    data, target = [], []

    for i in range(len(low_data)-maxlen):
        data.append(low_data[i:i + maxlen])
        target.append(low_data[i + maxlen])

    re_data = np.array(data).reshape(len(data), maxlen, 1)
    re_target = np.array(target).reshape(len(data), 1)

    return re_data, re_target

# RNNへの入力データ数
window_size = 12

# 入力データと教師データへの分割
X, y = make_dataset(df["y"], window_size)

データとテストデータに分割し、訓練データをさらに、訓練データと検証データに分割する。

# データの分割
X_train, X_test, y_train, y_test = temporal_train_test_split(X, y, test_size=TIME_SPAN)
# validation データの作成
X_train, X_val, y_train, y_val = temporal_train_test_split(X_train, y_train, test_size=TIME_SPAN)

3.2.2 LSTMモデルの構築

LSTMモデルを構築していく。
今回は、損失関数にはMSE(平均2乗誤差)を使用し、最適化関数にはadamを使用する。

# ネットワークの構築
model = Sequential() # Sequentialモデル

model.add(LSTM(50, batch_input_shape=(None, window_size, 1))) # LSTM 50層
model.add(Dense(1)) # 出力次元数は1

model.compile(loss='mean_squared_error', optimizer=Adam() , metrics = ['accuracy'])

3.2.3 モデルの学習

# 学習用パラメータ
batch_size = 20
n_epoch = 150

# 学習
hist = model.fit(X_train, y_train,
                 epochs=n_epoch,
                 validation_data=(X_val, y_val),
                 verbose=0,
                 batch_size=batch_size)

# 損失値(Loss)の遷移のプロット
plt.plot(hist.history['loss'],label="train set")
plt.plot(hist.history['val_loss'],label="test set")
plt.title('model loss')
plt.xlabel('epoch')
plt.ylabel('loss')
plt.legend()
plt.show()

3.2.4 モデルの予測

# 予測
pred_lstm = model.predict(X_test)
plt.plot(pred_lstm)
plt.show()

4.比較

最後に各モデルでの予測結果と実際の値をグラフにして、可視化する。
また、各モデルの予測の精度を数値化する。

4.1 データの編集

まずは、予測期間の実際の値とSARIMAモデル、LSTMモデルでの予測の値をひとつのデータにまとめる。

compare_pred = test.reset_index()[["ds","y"]]
compare_pred = compare_pred.rename(columns={"y":"test_data"})
compare_pred["SARIMA"] = pred_SARIMA.to_numpy()
compare_pred["LSTM"] = pred_lstm

グラフで可視化
化すると以下のようになりました。

# 全グラフ表示
fig = plt.figure(figsize=(20, 8))

plt.title("Payton Manning Wiki Access Count 2007-12  2016-1")
plt.xlabel("time")
plt.ylabel("access count")
plt.plot(compare_pred["ds"], compare_pred["test_data"], label="test data")
plt.plot(compare_pred["ds"], pred_SARIMA, label="SARIMA") # SARIMA
plt.plot(compare_pred["ds"], pred_lstm, label="LSTM") # LSTM
plt.legend()
plt.show()

4.3 各モデルの精度を数値化

決定係数R²と2乗平均平方根誤差の値でモデルの精度を算出してみました。

print("SARIMA")
print("R2スコア",r2_score(compare_pred["test_data"], compare_pred["SARIMA"]))
print("RMSE", np.sqrt(mean_squared_error(compare_pred["test_data"], compare_pred["SARIMA"])),"\n")

print("LSTM")
print("R2スコア",r2_score(compare_pred["test_data"], compare_pred["LSTM"]))
print("RMSE", np.sqrt(mean_squared_error(compare_pred["test_data"], compare_pred["LSTM"])),"\n")

SARIMA
R2スコア 0.23510242560205252
RMSE 0.6846035144169879

LSTM
R2スコア 0.3777788155994898
RMSE 0.6174614934374342

5 まとめ

今回は、初めてSARIMAモデルとLSTMモデルで、ウィキペディアのアクセス数の時系列解析を行っていきました。
グラフで可視化するだけだと、どちらのモデルが精度が高いのか少しわかりづらいですね。。。
ただ、数値化すると、R2スコアの値もRMSEの値もどちらもLSTMモデルのほうが精度が高かったです。

今回の2つのモデルで解析してみると、LSTMモデルのほうが精度が高かったですが、どちらもR2スコアが少し低めなので、精度としては微妙な感じがしますね。
SARIMAモデルもLSTMモデルも、パラメーターや層の数を調節したら、より高い精度のものが作れるのかなとも思います。
これからもっと勉強して、より様々な分析を出来るように頑張っていきたいと思います。