ローカルレベルとローカルトレンドモデル

🎓 レベル：標準　|　重要度：A（必須）

📎 前提：カルマンフィルタ　|　関連：ETSモデルと状態空間表現（ETS⟺ARIMA）・予測の評価指標と時系列CV（ホールドアウト）

要点（BLUF）

ローカルレベルモデル＝水準がランダムウォークで漂い、観測はそれにノイズが乗る最小の構造時系列。ローカル線形トレンドモデル＝水準に加え傾き $\nu_t$ も確率的に動く（傾きもランダムウォーク）。
信号対雑音比 $Q/H$ （状態ノイズ÷観測ノイズ）が状態の滑らかさを決めます。 $Q/H$ 大＝水準が機敏に動く、小＝滑らか。UnobservedComponents がこの分散をデータから復元します。
ローカルレベルは ARIMA(0,1,1)、ローカル線形トレンドは ARIMA(0,2,2) と等価（ETSモデルと状態空間表現の ETS⟺ARIMA と接続）。トレンド系列ではローカル線形トレンドが傾きを外挿し、平坦予測のローカルレベルを大きく上回ります。

1. ローカルレベルモデル（水準のランダムウォーク）

最小の構造時系列。隠れた水準 $\mu_t$ がランダムウォーク、観測 $y_t$ はそれに測定ノイズが乗ったもの（状態空間モデルの枠組みの最小例）。

\begin{aligned} y_t &= \mu_t + \varepsilon_t, & \varepsilon_t &\sim N(0, H)\\ \mu_t &= \mu_{t-1} + \eta_t, & \eta_t &\sim N(0, Q) \end{aligned}

信号対雑音比 $q=Q/H$ が挙動を決めます。

$q$ 大（状態ノイズが優勢）：水準が素早く動くので、推定はほぼ生の観測に追従（カルマンゲイン $K\to1$ ）。
$q$ 小（観測ノイズが優勢）：水準は滑らか、推定は観測を強く均す（ $K$ 小、カルマンフィルタのゲイン実験）。
$q\to0$ ：水準が動かない＝定数平均。 $q\to\infty$ ：観測ノイズ無視＝純ランダムウォーク。

ローカルレベルの $h$ 期先点予測は最終水準で平坦（ランダムウォークの最良予測は現在値）。だから点予測は素朴予測に近いのですが、状態空間なので較正された予測区間が出るのが値打ちです。

コード①：ローカルレベルの分散を復元し、予測区間つきで予測

真の $Q=0.5,\ H=4.0$ を仕込んだローカルレベル系列に UnobservedComponents(level="local level") を当て、最尤で $Q,H$ を復元します。フィルタ水準が真の水準と相関するか、get_forecast().conf_int() の予測区間がホールドアウトを捕らえ、ホライズンとともに広がるかを確かめます。

import warnings
warnings.simplefilter("ignore")
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import japanize_matplotlib
from statsmodels.tsa.statespace.structural import UnobservedComponents

# 真のローカルレベル：水準=ランダムウォーク, 観測=水準+ノイズ
np.random.seed(2)
n = 200
Q_true, H_true = 0.5, 4.0          # 状態ノイズ Q, 観測ノイズ H
alpha = np.zeros(n); alpha[0] = 30.0
for t in range(1, n):
    alpha[t] = alpha[t-1] + np.random.normal(0, np.sqrt(Q_true))
y = alpha + np.random.normal(0, np.sqrt(H_true), n)
y = pd.Series(y, index=pd.date_range("2010-01-01", periods=n, freq="MS"))

h = 24
train, test = y.iloc[:-h], y.iloc[-h:]

res = UnobservedComponents(train, level="local level").fit(disp=False)
p = dict(zip(res.param_names, res.params))
print("分散パラメータの復元（真値 -> 推定値）:")
print(f"  Q 状態（水準）: {Q_true:.3f} -> {p['sigma2.level']:.3f}")
print(f"  H 観測       : {H_true:.3f} -> {p['sigma2.irregular']:.3f}")
print(f"  推定 信号対雑音比 Q/H = {p['sigma2.level']/p['sigma2.irregular']:.3f}（真値 {Q_true/H_true:.3f}）")

level_filt = res.filtered_state[0]
corr = np.corrcoef(level_filt, alpha[:-h])[0, 1]
print(f"フィルタ水準と真の水準の相関 = {corr:.3f}")

fc = res.get_forecast(steps=h)
mean = fc.predicted_mean; ci = fc.conf_int(alpha=0.05)
ci_lo, ci_hi = ci.iloc[:, 0].values, ci.iloc[:, 1].values
covered = ((test.values >= ci_lo) & (test.values <= ci_hi)).sum()
print(f"95%予測区間カバレッジ = {covered}/{h}")
print(f"区間幅 1期先={ci_hi[0]-ci_lo[0]:.2f}  24期先={ci_hi[-1]-ci_lo[-1]:.2f}（先ほど広がる）")

fig, ax = plt.subplots(figsize=(9, 4.5))
ax.plot(train.index, train.values, color="0.6", lw=1, label="訓練データ")
ax.plot(test.index, test.values, color="k", lw=1.2, label="実測（検証）")
ax.plot(mean.index, mean.values, color="C1", lw=2, label="点予測（最終水準で平坦）")
ax.fill_between(mean.index, ci_lo, ci_hi, color="C1", alpha=0.25, label="95%予測区間")
ax.axvline(train.index[-1], ls=":", color="k")
ax.set_xlabel("年月"); ax.set_ylabel("値"); ax.legend()
ax.set_title("ローカルレベルモデル：分散復元と予測区間つき予測")
plt.tight_layout(); plt.show()

出力：

分散パラメータの復元（真値 -> 推定値）:
  Q 状態（水準）: 0.500 -> 0.465
  H 観測       : 4.000 -> 3.886
  推定 信号対雑音比 Q/H = 0.120（真値 0.125）
フィルタ水準と真の水準の相関 = 0.870
95%予測区間カバレッジ = 24/24
区間幅 1期先=9.18  24期先=15.76（先ほど広がる）

出力の意味：最尤推定は $Q=0.465$ （真値 $0.5$ ）、 $H=3.886$ （真値 $4.0$ ）と、2つの分散をきれいに復元。信号対雑音比も $0.120$ で真値 $0.125$ にほぼ一致——「観測のどれだけが本物の水準変化で、どれだけが測定ノイズか」をデータから当てられたわけです。フィルタ水準は真の水準と相関 $0.870$ 。予測は最終水準で平坦（橙の水平線）で、95%区間が $24/24$ を捕らえ、幅は $9.18\to15.76$ とホライズンとともに広がります（和分系列の予測分散累積、ARMA・ARIMAモデル）。点予測は素朴予測（最終値）とほぼ同じ——ローカルレベルの値打ちは「点」ではなく、分散構造の復元と較正された区間にあります。

2. ローカル線形トレンドモデル（傾きも動く）

水準に加え、傾き $\nu_t$ も確率的に動く（傾きもランダムウォーク）モデル。トレンドの向きや勾配がゆっくり変わってよくなります。

\begin{aligned} y_t &= \mu_t + \varepsilon_t, & \varepsilon_t &\sim N(0, H)\\ \mu_t &= \mu_{t-1} + \nu_{t-1} + \xi_t, & \xi_t &\sim N(0, Q_{\text{level}})\\ \nu_t &= \nu_{t-1} + \zeta_t, & \zeta_t &\sim N(0, Q_{\text{trend}}) \end{aligned}

状態は $\alpha_t=(\mu_t,\nu_t)'$ の2次元。 $Q_{\text{trend}}\to0$ なら傾きは一定（決定論的トレンド）、大きいほど勾配が機敏に変わります。ローカルレベルと違い、 $h$ 期先予測は最終の水準から最終の傾きで直線的に外挿されます——だから本物のトレンドがある系列で強い。

コード②：トレンド系列で傾きを外挿（ローカルレベル・素朴と比較）

水準と傾きの両方が確率的に動く系列に、ローカル線形トレンドを当てて 24ヶ月先を予測します。比較としてローカルレベル（傾きを持たない）と素朴予測（最終値）も並べ、トレンドの外挿が効くことを RMSE で示します。予測区間つき・時間順ホールドアウトで評価します。

import warnings
warnings.simplefilter("ignore")
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import japanize_matplotlib
from statsmodels.tsa.statespace.structural import UnobservedComponents

# 真のローカル線形トレンド：水準＋ゆっくり変わる傾き
np.random.seed(5)
n = 160
Q_level, Q_trend, H_true = 0.02, 0.008, 3.0
mu = np.zeros(n); nu = np.zeros(n)
mu[0], nu[0] = 20.0, 0.4
for t in range(1, n):
    mu[t] = mu[t-1] + nu[t-1] + np.random.normal(0, np.sqrt(Q_level))
    nu[t] = nu[t-1] + np.random.normal(0, np.sqrt(Q_trend))
y = mu + np.random.normal(0, np.sqrt(H_true), n)
y = pd.Series(y, index=pd.date_range("2012-01-01", periods=n, freq="MS"))

h = 24
train, test = y.iloc[:-h], y.iloc[-h:]

# ローカル線形トレンド と 比較用に ローカルレベル も当てる
res_trend = UnobservedComponents(train, level="local linear trend").fit(disp=False)
res_level = UnobservedComponents(train, level="local level").fit(disp=False)

fc = res_trend.get_forecast(steps=h)
mean = fc.predicted_mean; ci = fc.conf_int(alpha=0.05)
ci_lo, ci_hi = ci.iloc[:, 0].values, ci.iloc[:, 1].values

rmse_trend = np.sqrt(np.mean((mean.values - test.values)**2))
fc_level = res_level.get_forecast(steps=h).predicted_mean.values
rmse_level = np.sqrt(np.mean((fc_level - test.values)**2))
naive = np.full(h, train.values[-1])
rmse_naive = np.sqrt(np.mean((naive - test.values)**2))
covered = ((test.values >= ci_lo) & (test.values <= ci_hi)).sum()
slope_est = res_trend.filtered_state[1][-1]   # 終端の推定傾き
print(f"終端の推定傾き nu = {slope_est:.3f}（真の終端 {nu[-h-1]:.3f}）")
print("ホールドアウト RMSE:")
print(f"  ローカル線形トレンド = {rmse_trend:.3f}")
print(f"  ローカルレベル       = {rmse_level:.3f}（トレンドを外挿できず平坦）")
print(f"  素朴(最終値)         = {rmse_naive:.3f}")
print(f"95%予測区間カバレッジ（トレンド）= {covered}/{h}")

fig, ax = plt.subplots(figsize=(9, 4.5))
ax.plot(train.index, train.values, color="0.6", lw=1, label="訓練データ")
ax.plot(test.index, test.values, color="k", lw=1.2, label="実測（検証）")
ax.plot(mean.index, mean.values, color="C1", lw=2, label="点予測（傾きを外挿）")
ax.fill_between(mean.index, ci_lo, ci_hi, color="C1", alpha=0.25, label="95%予測区間")
ax.axvline(train.index[-1], ls=":", color="k")
ax.set_xlabel("年月"); ax.set_ylabel("値"); ax.legend()
ax.set_title("ローカル線形トレンドモデル：傾きを外挿した予測")
plt.tight_layout(); plt.show()

出力：

終端の推定傾き nu = 1.623（真の終端 1.809）
ホールドアウト RMSE:
  ローカル線形トレンド = 1.853
  ローカルレベル       = 24.068（トレンドを外挿できず平坦）
  素朴(最終値)         = 24.163
95%予測区間カバレッジ（トレンド）= 24/24

出力の意味：ローカル線形トレンドは終端の傾きを $\nu=1.623$ （真値 $1.809$ ）と推定し、それを延長して予測。ホールドアウト RMSE は $1.853$ と、ローカルレベル（ $24.068$ ）・素朴（ $24.163$ ）を桁違いに上回ります。理由は明快で、系列は上昇トレンドを持つのに、ローカルレベルと素朴は「最終水準で平坦」に予測するため、24ヶ月で系列がぐんぐん離れて RMSE が爆発する。トレンドモデルだけが傾きを状態として持ち、外挿できる。区間も $24/24$ で較正済み。モデルに正しいトレンド構造を入れているかが、これほど予測精度を分けます（モデル選択はモデル選択と残差診断）。

3. ARIMA との等価性

構造時系列の主役は ARIMA と裏でつながっています（ETSモデルと状態空間表現の ETS⟺ARIMA と同じ話）。

ローカルレベル ⟺ ARIMA(0,1,1)：1階差分すると MA(1) になる。ETS(A,N,N)＝SES とも等価。
ローカル線形トレンド ⟺ ARIMA(0,2,2)：2階差分で MA(2)。ETS(A,A,N)＝Holt とも等価。

同じ予測を出すのに、ARIMA は差分と MA で「自己相関」として、構造時系列は水準・傾きの「状態」として表します。後者の利点は、状態（水準・傾き）を取り出して分解・解釈でき、欠測や不規則間隔にカルマンで素直に対応できること（平滑化と欠測補間）。前者の利点は定常化理論の明快さと外生変数の扱いやすさ。同じ家族を別の言葉で見ているわけです。

4. 数式の直観

$Q/H$ は「どれだけ過去を信じるか」のつまみ：状態ノイズ $Q$ が大きいほど「世界は本当に動いている」とみなし観測に追従、 $H$ が大きいほど「観測は当てにならない」と均す。最尤推定はこのつまみをデータに語らせて決めます。
傾きを状態に持つ＝外挿できる：ローカルレベルは「今の高さ」しか覚えていないので予測は平坦。ローカル線形トレンドは「今の高さ＋今の勾配」を覚えるので直線で伸ばせる。何を状態に入れるかが、何を外挿できるかを決める。
構造時系列＝分解しながら予測：水準・傾き（さらに季節・周期）を別々の状態にすれば、ETS/STL（STL分解）と同じく「分解した成分」を見ながら予測区間つきで予測できます。

⚠️ よくある誤解

「ローカルレベルが素朴に勝てないのは失敗」ではない：ローカルレベル＝ARIMA(0,1,1) の最適予測は本質的に「賢い素朴予測」。点予測で素朴を大きく超えないのは正常で、値打ちは分散の復元と較正された区間にあります。
「トレンドモデルは常に上位」ではない：本物のトレンドが無い系列にローカル線形トレンドを当てると、傾きを過剰に外挿して遠未来で暴れます。トレンドの有無はホールドアウトで確認を（予測の評価指標と時系列CV）。
「分散パラメータは必ず復元できる」ではない：ローカル線形トレンドは水準ノイズと傾きノイズが識別しにくく、片方が0に潰れることがあります（滑らかトレンド化）。短い系列では特に不安定。AIC や予測精度で妥当性を点検します。
「予測区間は外挿でも信頼できる」ではない：区間は「モデルが正しく、ノイズが正規」が前提。トレンドの外挿は構造変化に弱く、遠未来ほど区間は名目より楽観的になりがちです。
「ローカルレベルとローカル線形トレンドだけ」ではない：季節成分や循環を加えた**基本構造モデル（BSM）**へ自然に拡張できます（seasonal= 等）。第7章のベイズ構造時系列でさらに柔軟化します。