トレンドの影響

import japanize_matplotlib
import pandas as pd
import matplotlib.pyplot as plt
import numpy as np
import seaborn as sns

サンプルデータを作成

date_list = pd.date_range("2021-01-01", periods=720, freq="D")
value_list = [
    10
    + np.cos(np.pi * i / 28.0) * (i % 3 > 0)
    + np.cos(np.pi * i / 14.0) * (i % 5 > 0)
    + np.cos(np.pi * i / 7.0)
    + (i / 10) ** 1.1 / 20
    for i, di in enumerate(date_list)
]

df = pd.DataFrame(
    {
        "日付": date_list,
        "観測値": value_list,
    }
)

df.head(10)
日付観測値
02021-01-0111.000000
12021-01-0212.873581
22021-01-0312.507900
32021-01-0411.017651
42021-01-0511.320187
52021-01-0610.246560
62021-01-079.350058
72021-01-089.740880
82021-01-099.539117
92021-01-108.987155
plt.figure(figsize=(10, 5))
sns.lineplot(x=df["日付"], y=df["観測値"])

png

XGBoostで時系列データを予測してみる

df["曜日"] = df["日付"].dt.weekday
df["年初からの日数%14"] = df["日付"].dt.dayofyear % 14
df["年初からの日数%28"] = df["日付"].dt.dayofyear % 28


def get_trend(timeseries, deg=3, trainN=0):
    """時系列データのトレンドの線を作成する

    Args:
        timeseries(pd.Series) : 時系列データ。
        deg(int) : 多項式の次数
        trainN(int): 多項式の係数を推定するために使用するデータ数

    Returns:
        pd.Series: トレンドに相当する時系列データ。
    """
    if trainN == 0:
        trainN = len(timeseries)

    x = list(range(len(timeseries)))
    y = timeseries.values
    coef = np.polyfit(x[:trainN], y[:trainN], deg)
    trend = np.poly1d(coef)(x)
    return pd.Series(data=trend, index=timeseries.index)


trainN = 500
df["トレンド"] = get_trend(df["観測値"], trainN=trainN, deg=2)

plt.figure(figsize=(10, 5))
sns.lineplot(x=df["日付"], y=df["観測値"])
sns.lineplot(x=df["日付"], y=df["トレンド"])
<AxesSubplot:xlabel='日付', ylabel='観測値'>

png

X = df[["曜日", "年初からの日数%14", "年初からの日数%28"]]
y = df["観測値"]

trainX, trainy = X[:trainN], y[:trainN]
testX, testy = X[trainN:], y[trainN:]
trend_train, trend_test = df["トレンド"][:trainN], df["トレンド"][trainN:]

トレンドを考慮せずに予測する

import xgboost as xgb
from sklearn.metrics import mean_squared_error

regressor = xgb.XGBRegressor(max_depth=5).fit(trainX, trainy)
prediction = regressor.predict(testX)

plt.figure(figsize=(10, 5))
sns.lineplot(x=df["日付"][trainN:], y=prediction)
sns.lineplot(x=df["日付"][trainN:], y=testy)

plt.legend(["モデルの出力", "正解"], bbox_to_anchor=(0.0, 0.78, 0.28, 0.102))
print(f"誤差(MSE) = {mean_squared_error(testy, prediction)}")
誤差(MSE) = 2.815118389938834

png

トレンドを考慮して予測する

先にトレンドに相当する分を観測値から取り除いて、トレンドがない状態の数値を予測します。 その後、トレンドに相当する分を足して最終的な予測値としています。

regressor = xgb.XGBRegressor(max_depth=5).fit(trainX, trainy - trend_train)
prediction = regressor.predict(testX)
prediction = [pred_i + trend_i for pred_i, trend_i in zip(prediction, trend_test)]

plt.figure(figsize=(10, 5))
sns.lineplot(x=df["日付"][trainN:], y=prediction)
sns.lineplot(x=df["日付"][trainN:], y=testy)

plt.legend(["モデルの出力+トレンド", "正解"], bbox_to_anchor=(0.0, 0.78, 0.28, 0.102))
print(f"誤差(MSE) = {mean_squared_error(testy, prediction)}")
誤差(MSE) = 0.46014173311011325

png