【ベンチマーク】複数の周期が重なるデータ

最終更新 2026-03-06 読了時間 6 分

まとめ

電力消費量や Web トラフィックのように、曜日の波（周期 7）と年間の波（周期 365）が同時に乗っている日次データを扱います。Holt-Winters や SARIMAX の基本形は seasonal_periods を 1 つしか持てないので、週次に合わせると年次を取りこぼし、年次に合わせると 365 個の季節係数で推定が不安定になります。このページでは振幅 10 の週次と振幅 25 の年次を重ねた 730 日の合成データに 6 つのモデルを当て、末尾 28 日の MASE で比べます。年次の波を sin(2πkt/365.25) と cos(2πkt/365.25) の 3 次までのフーリエ項として SARIMAX の exog に渡す方法が、週次と年次を同時に捉えて唯一 MASE を大きく下げます。残差の ACF を見ると、週次だけのモデルには年次成分の名残が長周期の相関として残ることが確認できます。

電力消費量・気温・Webトラフィックなど、「曜日の波」と「年間の波」が同時に存在するデータは多い。ここでは2つの周期が重なる合成データを生成し、季節パラメーターが1つしか設定できないモデルの限界と、フーリエ特徴量による回避策を検証する。

1. データの特徴 #

このシナリオの難しさは、seasonal_periods を 1 つしか指定できないモデル（Holt-Winters や SARIMAX の基本形）が大半であることです。週次パターン（sp=7）を指定すれば年次パターンを見逃し、年次パターン（sp=365）を指定すれば季節係数が 365 個になってパラメータ数が爆発し、2 年分 730 点では 1 係数あたり 2 点しか使えずフィッティングが不安定になります。回避策は、長い周期の方を少数の sin/cos 項で近似して外生変数に渡すことで、3 次のフーリエ項なら 6 列で年次の滑らかな波を表せます。

週次 7 日 + 年次 365 日の波を重ねる短い周期 7 は seasonal_order で扱う長い周期 365 は sin/cos 3 次 = 6 列で近似 SARIMAX の exog に渡す 28 日先を予測して MASE で比較残差 ACF に周期の名残がないか確認

何が起きているか（ステップで見る） #

コードを読む前に、この手法が図の上で何をしているかを追ってみましょう。「再生」を押すか、スライダーを動かすと 1 手順ずつ進みます。

見るポイント：水準 → 傾き → 季節と要素を足すほど、未来への伸ばし方が豊かになる。

2. 合成データの生成 #

次のコードは、水準 200 に周期 7 日・振幅 10 の波、周期 365.25 日・振幅 25 の波、1 日 0.02 のトレンド、標準偏差 5 のノイズを足した 730 日の系列を作ります。上段に全期間、下段に直近 60 日を拡大して描き、週次の波が見えるようにします。

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

rng = np.random.default_rng(42)
n = 730  # 2年分
dates = pd.date_range("2022-01-01", periods=n, freq="D")
t = np.arange(n)

# 週次周期 + 年次周期 + 緩やかなトレンド
weekly = 10 * np.sin(2 * np.pi * t / 7)
yearly = 25 * np.sin(2 * np.pi * t / 365.25)
trend = 0.02 * t
noise = rng.normal(0, 5, n)
series = pd.Series(200 + trend + weekly + yearly + noise, index=dates)

h = 28
train_raw = series.iloc[:-h]
test_raw = series.iloc[-h:]

fig, axes = plt.subplots(2, 1, figsize=(12, 6), gridspec_kw={"height_ratios": [2, 1]})
axes[0].plot(series.index, series.values, color="#94a3b8", linewidth=0.5)
axes[0].axvline(series.index[-h], color="#ef4444", linestyle="--")
axes[0].set_title("週次 + 年次の複合季節データ（赤線より右がテスト）")
axes[0].grid(alpha=0.2)

# 直近2か月を拡大
recent = series.iloc[-60:]
axes[1].plot(recent.index, recent.values, color="#2563eb", linewidth=1)
axes[1].set_title("直近2か月の拡大（週次の波が見える）")
axes[1].grid(alpha=0.2)

fig.tight_layout()
plt.show()

コードの解説番号をクリックすると、コード内の該当行がハイライトされます。

11-12 行目振幅の違う 2 つの正弦波を重ねる
週次は振幅 10、年次は振幅 25 なので、年次の方が 2.5 倍大きな波です。週次だけを捉えるモデルは、この振幅 25 の波を残差として抱えることになり、28 日のテスト期間では年次の波の局面（上昇中か下降中か）に応じて予測が一方向にずれます。周期を 365 ではなく 365.25 にしているのは、うるう年を含む実データに合わせた習慣です。
13-14 行目トレンドとノイズは控えめに
1 日 0.02 のトレンドは 2 年で 14.6 の上昇で、年次の振幅 25 より小さく、季節性の比較を邪魔しない程度です。ノイズの標準偏差 5 は週次の振幅 10 の半分で、週次の波がノイズに埋もれず、かつ完全には滑らかでない現実的な見え方になります。
17-19 行目テストは 28 日 = 週次 4 周期
h = 28 は週次周期 7 のちょうど 4 倍で、テスト期間に曜日の波が 4 回そろいます。一方で 28 日は年次周期の 8% にすぎないので、テスト期間中の年次成分はほぼ単調な上昇か下降に見えます。この「週次は周期として、年次はトレンドのように見える」構図が、モデルごとの外れ方を決めます。
21-25 行目全期間では年次の波しか見えない
730 日を 12 インチ幅に描くと 1 日あたり 0.016 インチなので、周期 7 の波は潰れて見えず、振幅 25 の年次の波だけが 2 山の形で現れます。height_ratios=[2, 1] で上段を大きくし、年次の形を主役にしています。
28-30 行目直近 60 日を拡大すると週次が見える
series.iloc[-60:] で最後の 2 か月だけを切り出すと、7 日ごとの山谷が 8〜9 回並び、週次の波がはっきり見えます。全体図と拡大図を重ねることで「このデータには 2 つの周期がある」ことを図だけで納得できる作りです。

複数周期が重なるデータ

3. 比較するパイプライン #

モデル	季節性の扱い	備考
季節ナイーブ (sp=7)	週次のみ	年次は無視
SES	なし	両方無視
Holt-Winters 加法 (sp=7)	週次のみ	年次は残差に漏出
ARIMA(2,1,1)	なし	両方無視
SARIMAX(1,0,1)(1,1,1,7)	週次のみ	年次は捉えられない
SARIMAX + フーリエ(年次)	週次 + 年次	sin/cos項を外生変数に追加

4. 実験と結果 #

次のコードは年次周期を 3 次のフーリエ項（6 列）に変換する make_fourier を定義し、6 つのモデルで 28 日先を予測して MASE を横棒グラフにします。季節ナイーブ基準の m=7 で MASE を計算し、1 を下回った棒だけ緑で塗ります。

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
import warnings
from statsmodels.tsa.holtwinters import ExponentialSmoothing, SimpleExpSmoothing
from statsmodels.tsa.statespace.sarimax import SARIMAX
from statsmodels.tsa.arima.model import ARIMA

warnings.filterwarnings("ignore")

def compute_mase(actual, predicted, train, m=7):
    e = np.abs(actual - predicted)
    scale = np.abs(train.values[m:] - train.values[:-m]).mean()
    return e.mean() / scale if scale > 0 else np.nan

def make_fourier(index, period, order):
    """sin/cos のフーリエ特徴量を生成。"""
    t = np.arange(len(index))
    cols = {}
    for k in range(1, order + 1):
        cols[f"sin_{period}_{k}"] = np.sin(2 * np.pi * k * t / period)
        cols[f"cos_{period}_{k}"] = np.cos(2 * np.pi * k * t / period)
    return pd.DataFrame(cols, index=index)

# フーリエ特徴量（年次周期をsin/cosで表現）
fourier_all = make_fourier(series.index, period=365.25, order=3)
fourier_train = fourier_all.iloc[:-h]
fourier_test = fourier_all.iloc[-h:]

results = {}

# 1. 季節ナイーブ
last7 = train_raw.values[-7:]
results["季節ナイーブ(7)"] = np.tile(last7, 4)[:h]

# 2. SES
results["SES"] = SimpleExpSmoothing(
    train_raw, initialization_method="estimated"
).fit().forecast(h).values

# 3. HW加法 (sp=7)
results["HW加法(7)"] = ExponentialSmoothing(
    train_raw, trend="add", seasonal="add", seasonal_periods=7,
    initialization_method="estimated",
).fit().forecast(h).values

# 4. ARIMA(2,1,1)
results["ARIMA"] = ARIMA(train_raw, order=(2, 1, 1)).fit().forecast(steps=h).values

# 5. SARIMAX(1,0,1)(1,1,1,7)
results["SARIMAX(7)"] = SARIMAX(
    train_raw, order=(1, 0, 1), seasonal_order=(1, 1, 1, 7),
    enforce_stationarity=False, enforce_invertibility=False,
).fit(disp=False).forecast(steps=h).values

# 6. SARIMAX + フーリエ（年次）
model_f = SARIMAX(
    train_raw, exog=fourier_train, order=(1, 0, 1), seasonal_order=(1, 1, 1, 7),
    enforce_stationarity=False, enforce_invertibility=False,
).fit(disp=False)
results["SARIMAX+フーリエ"] = model_f.forecast(steps=h, exog=fourier_test).values

# MASE 計算
mase_scores = {name: compute_mase(test_raw.values, pred, train_raw, m=7)
               for name, pred in results.items()}

# 棒グラフ
names = list(mase_scores.keys())
scores = [mase_scores[n] for n in names]
colors = ["#10b981" if s < 1 else "#ef4444" for s in scores]

fig, ax = plt.subplots(figsize=(10, 5))
bars = ax.barh(names, scores, color=colors, edgecolor="white")
ax.axvline(1.0, color="#1e40af", linestyle="--", linewidth=1, label="MASE = 1（季節ナイーブ水準）")
ax.set_xlabel("MASE（低いほど良い）")
ax.set_title("複数季節データ: モデル別 MASE 比較")
ax.legend()
for bar, s in zip(bars, scores):
    ax.text(bar.get_width() + 0.05, bar.get_y() + bar.get_height() / 2,
            f"{s:.2f}", va="center", fontsize=9)
ax.grid(alpha=0.2, axis="x")
fig.tight_layout()
plt.show()

コードの解説番号をクリックすると、コード内の該当行がハイライトされます。

8-11 行目MASE の基準は周期 7 の季節ナイーブ
m=7 なので分母は「7 日前の値をそのまま使った予測」の学習期間内 MAE です。日次データでは前週同曜日が最も素朴で強い基準になるため、これを 1.0 として各モデルを相対評価します。年次を基準にしないのは、365 日前の値を使う予測が 2 年分のデータでは 1 年分しか計算できず、基準として不安定だからです。
13-20 行目長い周期を sin/cos の少数の列で近似する
周期 365.25 の波を sin(2πkt/P) と cos(2πkt/P) の k = 1〜3 の組、計 6 列で表します。sin と cos を対で入れるのは、位相（山がいつ来るか）を係数の比で自由に表せるようにするためで、どちらか一方だけでは山の位置が固定されます。季節係数 365 個を推定する代わりに 6 個の回帰係数で済むのが、フーリエ近似の最大の利点です。
23-25 行目フーリエ項は学習とテストを同じ時間軸で切る
make_fourier を全 730 日分で一度に作ってから iloc で学習と予測に分けるのは、t の起点を揃えるためです。学習用とテスト用を別々に np.arange から作ると、テスト側の t が 0 から始まって年次の位相がずれ、予測が丸ごと間違った季節の値になります。
48-51 行目SARIMAX(1,0,1)(1,1,1,7) は週次だけを持つ
seasonal_order=(1, 1, 1, 7) は周期 7 の季節差分と季節 AR・MA を各 1 次入れた指定で、曜日の波は捉えますが年次を表す項はありません。enforce_stationarity=False と enforce_invertibility=False は、係数の制約を外して最適化が止まらないようにする設定で、合成データの比較では収束を優先しています。
54-58 行目exog にフーリエ項を渡すと年次が回帰項として入る
同じ SARIMAX に exog=fourier_train を足すだけで、6 列のフーリエ項の係数が ARMA 部分と同時に推定されます。予測時には exog=fourier_test で未来 28 日分の sin/cos を渡す必要があり、これを忘れると forecast が例外になります。週次は seasonal_order、年次は exog という分担で 2 つの周期を同時に扱えるのがこのモデルだけです。
30-31 行目季節ナイーブは直近 7 日を 4 回繰り返す
np.tile(last7, 4)[:h] で学習期間最後の 1 週間を 28 日に伸ばします。前週と同じ曜日の値を使うので週次は完璧に再現しますが、年次の上昇・下降には追随せず、テスト期間中に年次成分が動いた分だけ一方向に外れます。

モデル別MASE比較

読み方のポイント #

SARIMAX + フーリエだけが週次と年次の両方を捉え、MASE が最も低い緑の棒になります。SARIMAX(7) と HW加法(7) は週次は捉えるものの、振幅 25 の年次成分が残差に漏れるため MASE が 1 を切るかどうかの境目にとどまり、テスト 28 日間の年次の局面によっては季節ナイーブに負けることもあります。SES と ARIMA(2,1,1) は季節性を持たないので曜日の波に毎日 ±10 ずれ、MASE が 1 を大きく超えて赤くなります。棒グラフの順序は実行環境の statsmodels のバージョンで多少入れ替わりますが、「フーリエ付きが最良、季節性なしが最悪」という構図は変わりません。

5. 誤差パターン分析 #

MASE では「どの周期を取りこぼしたか」が分かりません。次のコードは SARIMAX + フーリエ、SARIMAX(7)、SES の 3 つについて、左にテスト 28 日の実測と予測、右に残差の ACF（ラグ 14 まで）を並べ、残差に週次や年次の名残が残っているかを確かめます。

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
from statsmodels.graphics.tsaplots import plot_acf

targets = [
    ("SARIMAX+フーリエ", "#2563eb"),
    ("SARIMAX(7)", "#f97316"),
    ("SES", "#ef4444"),
]

fig, axes = plt.subplots(len(targets), 2, figsize=(13, 3.5 * len(targets)))
test_idx = test_raw.index

for i, (name, color) in enumerate(targets):
    pred = results[name]
    residuals = test_raw.values - pred
    mase = mase_scores[name]

    ax_l = axes[i, 0]
    ax_l.plot(test_idx, test_raw.values, "o-", color="#94a3b8", markersize=3, label="実測")
    ax_l.plot(test_idx, pred, "s--", color=color, markersize=3, label="予測")
    ax_l.set_title(f"{name}  (MASE={mase:.2f})")
    ax_l.legend(fontsize=8)
    ax_l.grid(alpha=0.2)

    ax_r = axes[i, 1]
    plot_acf(residuals, ax=ax_r, lags=14, color=color, title="")
    ax_r.set_title("残差ACF")

fig.suptitle("複数季節データにおける誤差パターン比較", fontsize=13, y=1.02)
fig.tight_layout()
plt.show()

コードの解説番号をクリックすると、コード内の該当行がハイライトされます。

3-7 行目両方・週次だけ・なしの 3 段階を選ぶ
2 周期を捉えたモデル、1 周期だけのモデル、季節性なしのモデルを 1 本ずつ選んでいます。この 3 段階を縦に並べると、残差の ACF が「ほぼ白色雑音 → 長周期の名残 → ラグ 7 のスパイク」と段階的に悪くなる様子が 1 枚で読めます。
14 行目残差 = 実測 − 予測の符号で外れの向きを読む
残差が正ならモデルは低く見積もり、負なら高く見積もっています。週次だけのモデルでは、テスト期間中に年次成分が上昇局面なら残差が全体に正へ、下降局面なら負へ偏り、左の図で予測の破線が実測の点列の片側に寄ることで見えます。
25 行目lags=14 で週次 2 周期分の相関を見る
plot_acf(residuals, lags=14) はラグ 1〜14 の自己相関を棒で描き、信頼区間を帯で示します。ラグ 7 と 14 に帯を超える棒があれば週次を取りこぼしており、ラグ 1 から数日にわたって正の相関が緩やかに続けば年次成分（28 日ではほぼ直線に見える）が残っています。テストが 28 点なので、ラグ 14 より先は推定が粗くなり信頼できません。

誤差パターン分析

SARIMAX + フーリエの残差はほぼ白色雑音になり、ACF の 14 本の棒は信頼区間の帯の中に収まります。SARIMAX(7) の残差には年次成分の名残が「ラグ 1〜5 で正の相関が緩やかに減る」形の長周期の自己相関として現れ、左の図では予測の破線が実測の点列から一定方向にずれています。SES は曜日の波を全く持たないので、残差に周期 7 の波がそのまま残り、ACF のラグ 7 と 14 に顕著なスパイクが立ちます。

6. よくある失敗パターン #

最初の失敗は seasonal_periods=365 を直接指定することで、Holt-Winters なら 365 個の季節係数を 2 年分 730 点から推定することになり、1 係数あたり 2 点では収束しないか過学習します。年次のような長い周期は 3〜5 次のフーリエ項（6〜10 列）で近似し、exog に渡すのが現実的です。二つ目は季節成分を 1 つ指定して満足することで、seasonal_periods=7 で MASE が 0.9 を切っても、残差 ACF に長周期の相関が残っていれば年次を取りこぼしています。plot_acf で残差を確認し、帯を超える棒の周期から足りない成分を特定します。三つ目はフーリエ項の次数を上げすぎることで、10 次（20 列）にすると年次の波に細かいギザギザまで当てはめ、2 年分のデータでは学習期間だけに合う形になります。

ウォークフォワード検証 — 時系列の正しい評価手法
MASE — 季節ナイーブ基準の精度指標
ARIMA モデル — 自己回帰と移動平均
ホルト・ウィンター法 — 季節性を含む指数平滑法