明日の気温を予測するとき、1年前の気温と昨日の気温のどちらをより重視するべきでしょうか?ほとんどの人は「昨日の気温の方が参考になる」と答えるでしょう。指数平滑法(Exponential Smoothing)は、まさにこの直感を数学的に定式化した手法です。最近のデータほど大きな重みを与え、古いデータの影響を指数的に減衰させます。
指数平滑法は1950年代にロバート・ブラウンによって提案された古典的な手法ですが、その発展形であるHolt-Winters法は、現代でも需要予測や在庫管理の現場で広く使われています。指数平滑法を理解することで、以下のような場面で活用できます。
- 需要予測・在庫管理 — 小売業やサプライチェーンにおけるリアルタイムの需要予測と在庫最適化
- ビジネスKPIの予測 — 売上高、ユーザー数、コストなどのトレンドと季節性を考慮した将来予測
- SARIMAモデルとの比較 — 状態空間モデルとの関係を通じて、現代的な時系列分析の理解を深める
本記事の内容
- 移動平均の限界と指数平滑法の動機
- 単純指数平滑法(SES)の再帰的更新式と重み構造
- Holtの二重指数平滑法(トレンド追従)
- Holt-Winters法(加法モデル・乗法モデル)
- パラメータの最適化(尤度最大化)
- Pythonでのスクラッチ実装とstatsmodelsの使い方
- SARIMAとの比較
前提知識
この記事を読む前に、以下の記事を読んでおくと理解が深まります。
移動平均の限界
時系列の平滑化において最も素朴な方法は単純移動平均(SMA)です。過去 $m$ 期のデータの平均を取ります。
$$ \hat{x}_t = \frac{1}{m} \sum_{i=0}^{m-1} x_{t-i} = \frac{x_t + x_{t-1} + \cdots + x_{t-m+1}}{m} $$
しかし、単純移動平均には以下の問題があります。
- 全てのデータに等しい重み — 昨日の値も $m$ 日前の値も同じ重み $1/m$ で扱われます。しかし、直感的には新しいデータほど予測に重要なはずです。
- 窓の長さの選択が難しい — $m$ が小さすぎるとノイズに敏感、大きすぎると変化への追従が遅くなります。
- 過去 $m$ 期分のデータを保持する必要がある — メモリの観点で非効率です。
これらの問題を解決するのが指数平滑法です。過去のデータに対して指数的に減衰する重みを与えることで、最近のデータを重視しつつ、過去の情報も活かす仕組みを実現します。しかも、再帰的な更新式を使うことで、データを1つずつ処理するだけで計算が完結します。
単純指数平滑法(SES)
直感的な理解
単純指数平滑法の基本的なアイデアは非常にシンプルです。「予測値を更新するとき、新しい観測値と前回の予測値の加重平均をとる」というものです。
$$ \begin{equation} \hat{x}_{t+1} = \alpha x_t + (1 – \alpha) \hat{x}_t \end{equation} $$
ここで $\alpha \in (0, 1)$ は平滑化パラメータです。この式は次のように解釈できます。
- $\alpha$ が1に近い場合: 最新の観測値 $x_t$ を強く信じる → 変化に素早く追従するが、ノイズにも敏感
- $\alpha$ が0に近い場合: 前回の予測値 $\hat{x}_t$ を強く信じる → 安定した予測だが、変化への追従が遅い
重み構造の導出
SESの式を再帰的に展開すると、過去のデータに対する重み構造が明らかになります。
まず、$\hat{x}_t = \alpha x_{t-1} + (1-\alpha)\hat{x}_{t-1}$ を元の式に代入すると、
$$ \hat{x}_{t+1} = \alpha x_t + (1-\alpha)[\alpha x_{t-1} + (1-\alpha)\hat{x}_{t-1}] $$
右辺を整理すると、
$$ \hat{x}_{t+1} = \alpha x_t + \alpha(1-\alpha) x_{t-1} + (1-\alpha)^2 \hat{x}_{t-1} $$
さらに $\hat{x}_{t-1}$ を展開して代入を繰り返すと、
$$ \begin{equation} \hat{x}_{t+1} = \alpha \sum_{j=0}^{t-1} (1-\alpha)^j x_{t-j} + (1-\alpha)^t \hat{x}_1 \end{equation} $$
$t$ が十分大きいとき $(1-\alpha)^t \approx 0$ なので、初期値 $\hat{x}_1$ の影響は無視できます。
この展開から重要なことがわかります。時刻 $t-j$ のデータに対する重みは $\alpha(1-\alpha)^j$ です。$0 < 1-\alpha < 1$ なので、$j$ が大きくなるほど重みは指数的に減衰します。これが「指数平滑法」という名前の由来です。
重みの正規化
全ての重みの和が1になることを確認しましょう。$t \to \infty$ の極限で、
$$ \alpha \sum_{j=0}^{\infty} (1-\alpha)^j = \alpha \cdot \frac{1}{1-(1-\alpha)} = \alpha \cdot \frac{1}{\alpha} = 1 $$
等比級数の公式を使いました。つまり、指数平滑法は正規化された重み付き平均であり、重みが指数的に減衰するという性質を持っています。
誤差修正形式
SESの更新式は、別の形に書き換えることもできます。$\hat{x}_{t+1} = \alpha x_t + (1-\alpha)\hat{x}_t$ の右辺に $\hat{x}_t$ を加えて引くと、
$$ \hat{x}_{t+1} = \hat{x}_t + \alpha(x_t – \hat{x}_t) $$
ここで $e_t = x_t – \hat{x}_t$ は時刻 $t$ における予測誤差です。したがって、
$$ \begin{equation} \hat{x}_{t+1} = \hat{x}_t + \alpha e_t \end{equation} $$
この形式は直感的に理解しやすいです。「新しい予測 = 前の予測 + 誤差の $\alpha$ 倍だけ修正」ということです。予測が外れた分を、$\alpha$ の割合でフィードバックして修正するのです。
単純指数平滑法はレベル(水準)の変化を捉えることはできますが、トレンドがある場合には常に後追いになってしまいます。次に、トレンドへの追従能力を追加したHoltの手法を見ていきましょう。
Holtの二重指数平滑法(Holt’s Linear Trend Method)
トレンドへの追従問題
単純指数平滑法は、トレンド(上昇傾向や下降傾向)がある時系列に対しては系統的にバイアスが生じます。上昇トレンドの場合、予測値は常に実際の値よりも低くなります。
この問題を解決するために、チャールズ・ホルト(C.C. Holt)は、レベル成分に加えてトレンド成分を別途追跡するアイデアを提案しました。
定式化
Holtの手法では、2つの状態変数を更新します。
レベル更新式: $$ \begin{equation} \ell_t = \alpha x_t + (1-\alpha)(\ell_{t-1} + b_{t-1}) \end{equation} $$
トレンド更新式: $$ \begin{equation} b_t = \beta(\ell_t – \ell_{t-1}) + (1-\beta) b_{t-1} \end{equation} $$
予測式: $$ \begin{equation} \hat{x}_{t+h} = \ell_t + h \cdot b_t \end{equation} $$
ここで、$\ell_t$ は時刻 $t$ におけるレベル(水準)、$b_t$ は時刻 $t$ におけるトレンド(傾き)、$\alpha \in (0,1)$ はレベルの平滑化パラメータ、$\beta \in (0,1)$ はトレンドの平滑化パラメータです。
各式の直感的解釈
レベル更新式(式4)の右辺は2項からなります。第1項 $\alpha x_t$ は最新の観測値、第2項 $(1-\alpha)(\ell_{t-1} + b_{t-1})$ はトレンドを考慮した前回の予測値です。つまり、トレンドを加味した予測と最新の観測の加重平均をとっているのです。
トレンド更新式(式5)では、$\ell_t – \ell_{t-1}$ がレベルの最新の変化量(最新の「傾き」の推定値)です。これと前回のトレンド推定値 $b_{t-1}$ の加重平均をとることで、トレンドを滑らかに更新しています。
予測式(式6)は、現在のレベルから $h$ ステップ先を線形外挿するものです。$b_t > 0$ なら上昇トレンド、$b_t < 0$ なら下降トレンドを予測します。
減衰トレンド(Damped Trend)
Holtの線形トレンドは、長期予測では過度に楽観的(または悲観的)になりがちです。実際には、トレンドは時間とともに減衰することが多いです。Gardner & McKenzie(1985)はトレンドに減衰係数 $\phi \in (0, 1)$ を導入しました。
$$ \hat{x}_{t+h} = \ell_t + (\phi + \phi^2 + \cdots + \phi^h) b_t = \ell_t + \frac{\phi(1 – \phi^h)}{1 – \phi} b_t $$
$\phi < 1$ のとき、$h \to \infty$ で予測値は $\ell_t + \phi/(1-\phi) \cdot b_t$ に収束します。つまり、長期的にはトレンドが飽和し、予測値が一定値に近づきます。
Holtの手法でトレンドを追従できるようになりましたが、季節性のある時系列にはまだ対応できていません。次に、季節成分を追加したHolt-Winters法を見ていきましょう。
Holt-Winters法
加法モデルと乗法モデル
季節変動の振幅が時間とともに一定である場合(例えば、売上の季節変動が常に100万円程度)は加法モデル(additive model)が適切です。一方、季節変動の振幅が全体のレベルに比例して増大する場合(例えば、売上が倍になると季節変動も倍になる)は乗法モデル(multiplicative model)が適切です。
加法Holt-Winters法
加法モデルでは、観測値 $x_t$ がレベル $\ell_t$、トレンド $b_t$、季節成分 $s_t$、誤差 $e_t$ の和で表されます。
$$ x_t = \ell_t + b_t + s_t + e_t $$
更新式は以下の3つです。ここで $s$ は季節周期(月次データなら $s = 12$)です。
レベル更新: $$ \begin{equation} \ell_t = \alpha(x_t – s_{t-s}) + (1-\alpha)(\ell_{t-1} + b_{t-1}) \end{equation} $$
トレンド更新: $$ \begin{equation} b_t = \beta(\ell_t – \ell_{t-1}) + (1-\beta)b_{t-1} \end{equation} $$
季節更新: $$ \begin{equation} s_t = \gamma(x_t – \ell_t) + (1-\gamma)s_{t-s} \end{equation} $$
予測式: $$ \begin{equation} \hat{x}_{t+h} = \ell_t + h \cdot b_t + s_{t+h-s} \end{equation} $$
ここで $\gamma \in (0, 1)$ は季節成分の平滑化パラメータです。
レベル更新式(式7)では、$x_t – s_{t-s}$ で季節成分を除去してからレベルを推定しています。前の季節周期の同じ位相の季節成分 $s_{t-s}$ を差し引くことで、「季節性を取り除いた真のレベル」を推定しているのです。
季節更新式(式9)では、$x_t – \ell_t$ が最新の季節成分の推定値です。これと1周期前の季節成分 $s_{t-s}$ の加重平均をとることで、季節パターンを緩やかに更新しています。
乗法Holt-Winters法
乗法モデルでは、観測値が各成分の積で表されます。
$$ x_t = \ell_t \cdot b_t \cdot s_t \cdot e_t $$
ただし、実用的にはトレンドは加法的に扱い、季節成分のみ乗法的に扱う形が多く使われます。
レベル更新: $$ \begin{equation} \ell_t = \alpha \frac{x_t}{s_{t-s}} + (1-\alpha)(\ell_{t-1} + b_{t-1}) \end{equation} $$
トレンド更新: $$ \begin{equation} b_t = \beta(\ell_t – \ell_{t-1}) + (1-\beta)b_{t-1} \end{equation} $$
季節更新: $$ \begin{equation} s_t = \gamma \frac{x_t}{\ell_t} + (1-\gamma)s_{t-s} \end{equation} $$
予測式: $$ \begin{equation} \hat{x}_{t+h} = (\ell_t + h \cdot b_t) \cdot s_{t+h-s} \end{equation} $$
乗法モデルでは、季節成分 $s_t$ は比率(例えば、1月は年間平均の0.8倍、7月は1.3倍)として解釈されます。加法モデルの季節成分が「差」であるのに対し、乗法モデルの季節成分は「比」であるという違いがあります。
理論的な定式化を理解したところで、まずはスクラッチ実装でアルゴリズムの挙動を確認しましょう。
Pythonでのスクラッチ実装
単純指数平滑法のスクラッチ実装
import numpy as np
import matplotlib.pyplot as plt
def simple_exponential_smoothing(data, alpha, initial=None):
"""単純指数平滑法のスクラッチ実装"""
n = len(data)
smoothed = np.zeros(n)
# 初期値: 最初の観測値を使用
smoothed[0] = initial if initial is not None else data[0]
for t in range(1, n):
smoothed[t] = alpha * data[t] + (1 - alpha) * smoothed[t - 1]
return smoothed
# サンプルデータ生成(レベルシフトのあるノイズ)
np.random.seed(42)
n = 200
level = np.concatenate([np.ones(80) * 10, np.ones(60) * 15, np.ones(60) * 12])
data = level + np.random.normal(0, 1.5, n)
# 異なるαでの平滑化
alphas = [0.1, 0.3, 0.7, 0.95]
fig, axes = plt.subplots(2, 2, figsize=(12, 8))
axes = axes.flatten()
for idx, alpha in enumerate(alphas):
smoothed = simple_exponential_smoothing(data, alpha)
axes[idx].plot(data, alpha=0.4, linewidth=0.8, label='Observed')
axes[idx].plot(smoothed, linewidth=2, color='red', label=f'SES (α={alpha})')
axes[idx].set_title(f'α = {alpha}')
axes[idx].legend()
axes[idx].set_ylim(5, 20)
plt.suptitle('Simple Exponential Smoothing with Different α', fontsize=13)
plt.tight_layout()
plt.show()
4つのパネルから、$\alpha$ の値がSESの挙動に与える影響が明確にわかります。
- $\alpha = 0.1$(左上): 平滑化が非常に強く、ノイズがよく除去されていますが、レベルシフトへの追従が遅いです。レベルが15に変わった後も、しばらくの間は低い値を予測し続けています。
- $\alpha = 0.3$(右上): バランスの取れた平滑化です。レベルシフトへの追従がやや速くなっています。
- $\alpha = 0.7$(左下): 最新のデータを重視するため、レベルシフトに素早く追従しますが、ノイズも拾いやすくなっています。
- $\alpha = 0.95$(右下): ほぼ最新の観測値そのものに近い予測となり、平滑化の効果がほとんどありません。
加法Holt-Winters法のスクラッチ実装
import numpy as np
import matplotlib.pyplot as plt
def holt_winters_additive(data, s, alpha, beta, gamma, n_forecast=0):
"""加法Holt-Winters法のスクラッチ実装"""
n = len(data)
# 初期化
level = np.zeros(n)
trend = np.zeros(n)
season = np.zeros(n + n_forecast)
fitted = np.zeros(n + n_forecast)
# 初期値: 最初の季節周期から推定
level[0] = np.mean(data[:s])
trend[0] = (np.mean(data[s:2*s]) - np.mean(data[:s])) / s
# 季節成分の初期値
for i in range(s):
season[i] = data[i] - level[0]
# 更新
for t in range(1, n):
# レベル更新
level[t] = alpha * (data[t] - season[t - s]) + (1 - alpha) * (level[t-1] + trend[t-1])
# トレンド更新
trend[t] = beta * (level[t] - level[t-1]) + (1 - beta) * trend[t-1]
# 季節更新
season[t] = gamma * (data[t] - level[t]) + (1 - gamma) * season[t - s]
# フィッティング値
fitted[t] = level[t] + trend[t] + season[t]
# 予測
for h in range(1, n_forecast + 1):
fitted[n - 1 + h] = level[n-1] + h * trend[n-1] + season[n - 1 + h - s]
return fitted, level, trend, season[:n]
# AirPassengersデータ風の合成データ
np.random.seed(42)
n = 144 # 12年分の月次データ
s = 12
time = np.arange(n)
trend_true = 100 + 2.5 * time
season_true = 30 * np.sin(2 * np.pi * time / 12) + 15 * np.cos(4 * np.pi * time / 12)
data = trend_true + season_true + np.random.normal(0, 10, n)
# フィッティングと予測
n_forecast = 24
fitted, level, trend, season = holt_winters_additive(
data, s=12, alpha=0.3, beta=0.05, gamma=0.3, n_forecast=n_forecast
)
fig, axes = plt.subplots(3, 1, figsize=(12, 10))
# 全体のフィット
time_all = np.arange(n + n_forecast)
axes[0].plot(time, data, 'o', markersize=2, alpha=0.5, label='Observed')
axes[0].plot(time_all[:n], fitted[:n], linewidth=1.5, color='red', label='Fitted')
axes[0].plot(time_all[n:], fitted[n:], linewidth=2, color='green',
linestyle='--', label='Forecast')
axes[0].axvline(x=n-1, color='gray', linestyle=':', alpha=0.5)
axes[0].set_title('Additive Holt-Winters (Scratch Implementation)')
axes[0].legend()
# レベルとトレンド
axes[1].plot(time, level, label='Level', linewidth=1.5)
axes[1].plot(time, trend_true[:n], '--', label='True Trend', alpha=0.5)
axes[1].set_title('Level Component')
axes[1].legend()
# 季節成分
axes[2].plot(time, season, label='Estimated Season', linewidth=1.5)
axes[2].plot(time, season_true, '--', label='True Season', alpha=0.5)
axes[2].set_title('Seasonal Component')
axes[2].legend()
plt.tight_layout()
plt.show()
スクラッチ実装の結果から、Holt-Winters法の各コンポーネントの挙動を確認できます。
- 全体のフィット(上段): 合成データの季節パターンとトレンドをよく捉えており、24ステップ先の予測も合理的です。
- レベル成分(中段): 推定されたレベルが真のトレンドにうまく追従しています。$\beta = 0.05$ と小さく設定したため、トレンドが安定的に推定されています。
- 季節成分(下段): 推定された季節パターンが真の季節成分とよく一致しています。$\gamma = 0.3$ で季節パターンを適度に更新しています。
次に、statsmodelsのExponentialSmoothingクラスを使って、実データでの分析を行いましょう。
statsmodelsによる実装
AirPassengersデータでの分析
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from statsmodels.datasets import get_rdataset
from statsmodels.tsa.holtwinters import ExponentialSmoothing
# データ準備
data = get_rdataset("AirPassengers").data
data.columns = ['time', 'passengers']
passengers = data['passengers'].values
time_index = pd.date_range('1949-01', periods=len(passengers), freq='MS')
ts = pd.Series(passengers, index=time_index)
# 訓練/テスト分割
train = ts[:'1958-12']
test = ts['1959-01':]
# 3つのモデルを比較
models = {}
# 1. 加法モデル
hw_add = ExponentialSmoothing(train, trend='add', seasonal='add',
seasonal_periods=12).fit()
models['Additive'] = hw_add
# 2. 乗法モデル
hw_mul = ExponentialSmoothing(train, trend='add', seasonal='mul',
seasonal_periods=12).fit()
models['Multiplicative'] = hw_mul
# 3. 減衰トレンド + 乗法季節
hw_damped = ExponentialSmoothing(train, trend='add', seasonal='mul',
seasonal_periods=12, damped_trend=True).fit()
models['Damped Mul'] = hw_damped
fig, axes = plt.subplots(3, 1, figsize=(12, 12))
colors = {'Additive': 'red', 'Multiplicative': 'blue', 'Damped Mul': 'green'}
for idx, (name, model) in enumerate(models.items()):
forecast = model.forecast(len(test))
axes[idx].plot(train.index, train, label='Train', linewidth=1.5, color='gray')
axes[idx].plot(test.index, test, label='Test', linewidth=1.5, color='black')
axes[idx].plot(forecast.index, forecast, label=f'{name} Forecast',
linewidth=2, color=colors[name], linestyle='--')
axes[idx].set_title(f'Holt-Winters {name} Model')
axes[idx].legend()
# 精度指標
mae = np.mean(np.abs(test.values - forecast.values))
mape = np.mean(np.abs((test.values - forecast.values) / test.values)) * 100
axes[idx].text(0.02, 0.95, f'MAE={mae:.1f}, MAPE={mape:.1f}%',
transform=axes[idx].transAxes, fontsize=10,
verticalalignment='top',
bbox=dict(boxstyle='round', facecolor='wheat', alpha=0.5))
plt.tight_layout()
plt.show()
# パラメータの表示
for name, model in models.items():
print(f"\n{name} Model Parameters:")
print(f" α (level) = {model.params['smoothing_level']:.4f}")
print(f" β (trend) = {model.params['smoothing_trend']:.4f}")
print(f" γ (seasonal) = {model.params['smoothing_seasonal']:.4f}")
if 'damping_trend' in model.params:
print(f" φ (damping) = {model.params['damping_trend']:.4f}")
print(f" AIC = {model.aic:.2f}")
3つのモデルの比較から、重要な知見が得られます。
- 加法モデル: 季節変動の振幅が一定と仮定するため、AirPassengersのような分散が増大するデータでは、後半の予測精度が低下しがちです。
- 乗法モデル: 季節変動の振幅がレベルに比例するため、AirPassengersデータの特性に合致しています。MAPEが最も小さくなることが期待されます。
- 減衰トレンド+乗法季節: 長期予測ではトレンドが徐々に飽和するため、過度に楽観的な予測を避けられます。2年先程度の予測であれば、通常のトレンドモデルとの差は小さいですが、さらに長期の予測では減衰トレンドが有利になります。
成分分解の可視化
Holt-Winters法がデータをどのように分解しているかを可視化します。
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from statsmodels.datasets import get_rdataset
from statsmodels.tsa.holtwinters import ExponentialSmoothing
# データ準備
data = get_rdataset("AirPassengers").data
data.columns = ['time', 'passengers']
passengers = data['passengers'].values
time_index = pd.date_range('1949-01', periods=len(passengers), freq='MS')
ts = pd.Series(passengers, index=time_index)
# 乗法モデルのフィッティング
model = ExponentialSmoothing(ts, trend='add', seasonal='mul',
seasonal_periods=12).fit()
fig, axes = plt.subplots(4, 1, figsize=(12, 10))
# 元データとフィッティング
axes[0].plot(ts.index, ts, linewidth=1.5, label='Observed')
axes[0].plot(ts.index, model.fittedvalues, linewidth=1.5,
color='red', alpha=0.7, label='Fitted')
axes[0].set_title('Observed vs Fitted')
axes[0].legend()
# レベル成分
axes[1].plot(ts.index, model.level, linewidth=1.5, color='darkorange')
axes[1].set_title('Level Component')
# トレンド成分
axes[2].plot(ts.index, model.trend, linewidth=1.5, color='green')
axes[2].set_title('Trend Component (Slope)')
# 季節成分
axes[3].plot(ts.index, model.season, linewidth=1.5, color='purple')
axes[3].axhline(y=1.0, color='gray', linestyle='--', alpha=0.5)
axes[3].set_title('Seasonal Component (Multiplicative)')
plt.tight_layout()
plt.show()
Holt-Winters法による成分分解の結果から、データの構造が明確に読み取れます。
- レベル成分: 全体的な水準が1949年の約100から1960年の約400まで一貫して上昇しています。
- トレンド成分: 傾き(1か月あたりの増加量)が時間とともにわずかに変動しながらも、概ね安定しています。
- 季節成分: 乗法モデルなので1.0を中心に振動しています。夏季(7〜8月)に1.2程度のピーク、冬季に0.85程度の谷を持つパターンが安定的に繰り返されています。
SARIMAとの比較
指数平滑法とSARIMAモデルは、しばしば競合する手法として比較されます。両者の違いを実験的に確認しましょう。
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from statsmodels.datasets import get_rdataset
from statsmodels.tsa.holtwinters import ExponentialSmoothing
from statsmodels.tsa.statespace.sarimax import SARIMAX
import warnings
warnings.filterwarnings('ignore')
# データ準備
data = get_rdataset("AirPassengers").data
data.columns = ['time', 'passengers']
passengers = data['passengers'].values
time_index = pd.date_range('1949-01', periods=len(passengers), freq='MS')
ts = pd.Series(passengers, index=time_index)
# 訓練/テスト分割
train = ts[:'1958-12']
test = ts['1959-01':]
n_test = len(test)
# Holt-Winters(乗法季節)
hw_model = ExponentialSmoothing(train, trend='add', seasonal='mul',
seasonal_periods=12).fit()
hw_forecast = hw_model.forecast(n_test)
# SARIMA(対数変換)
log_train = np.log(train)
sarima_model = SARIMAX(log_train, order=(0, 1, 1),
seasonal_order=(0, 1, 1, 12)).fit(disp=False)
sarima_forecast_log = sarima_model.forecast(n_test)
sarima_forecast = np.exp(sarima_forecast_log)
fig, ax = plt.subplots(figsize=(12, 5))
ax.plot(train.index, train, label='Train', linewidth=1.5, color='gray')
ax.plot(test.index, test, label='Test', linewidth=2, color='black')
ax.plot(test.index, hw_forecast, label='Holt-Winters (Mul)',
linewidth=2, linestyle='--', color='blue')
ax.plot(test.index, sarima_forecast, label='SARIMA(0,1,1)(0,1,1)₁₂',
linewidth=2, linestyle='--', color='red')
ax.set_title('Holt-Winters vs SARIMA — AirPassengers Forecast')
ax.legend()
plt.tight_layout()
plt.show()
# 精度比較
for name, forecast_vals in [('Holt-Winters', hw_forecast.values),
('SARIMA', sarima_forecast.values)]:
mae = np.mean(np.abs(test.values - forecast_vals))
mape = np.mean(np.abs((test.values - forecast_vals) / test.values)) * 100
rmse = np.sqrt(np.mean((test.values - forecast_vals) ** 2))
print(f"{name:15s}: MAE={mae:.1f}, MAPE={mape:.1f}%, RMSE={rmse:.1f}")
比較結果から、両モデルの特性が読み取れます。
- 予測精度: AirPassengersデータではHolt-Winters(乗法)とSARIMAの予測精度は比較的近いことが多いです。データの性質によってどちらが優れるかは変わります。
- 解釈性: Holt-Winters法はレベル・トレンド・季節成分の3つに分解されるため、ビジネスユーザーにとって解釈しやすいです。SARIMAのAR/MA係数は統計的な知識がないと解釈しにくいです。
- 計算速度: 指数平滑法はSARIMAよりも一般に計算が速く、リアルタイム処理に適しています。
- 柔軟性: SARIMAは外生変数を追加したSARIMAXへの拡張が容易です。指数平滑法にも状態空間モデルとしての定式化があり、理論的には等価な拡張が可能です。
まとめ
本記事では、指数平滑法からHolt-Winters法までの理論と実装を詳しく解説しました。
- 単純指数平滑法: 最近のデータほど大きな重みを与え、重みが指数的に減衰する加重平均。パラメータ $\alpha$ が唯一の調整項
- Holtの二重指数平滑法: レベルとトレンドの2つの状態変数を追跡し、線形トレンドに追従できる
- Holt-Winters法: 季節成分を追加し、加法モデルまたは乗法モデルでトレンド+季節性のある時系列を予測する
- 乗法モデル vs 加法モデル: 季節変動の振幅がレベルに比例する場合は乗法モデル、一定の場合は加法モデルが適切
指数平滑法は状態空間モデルの特殊ケースとして位置づけることもでき、ARIMAモデルとの理論的なつながりが深いです。
次のステップとして、以下の記事も参考にしてください。