モデル選択と残差診断

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

📎 前提：ARMA・ARIMAモデル　|　関連：予測の評価指標と時系列CV・訓練・検証・テストと交差検証（機械学習）

要点（BLUF）

次数選択： $(p,q)$ をグリッド探索し、AIC または BIC が最小のモデルを選ぶ。AIC $=-2\ln\hat L+2k$ 、BIC $=-2\ln\hat L+k\ln n$ 。BIC の方が罰則が重く倹約的（小さいモデルを選びがち）。
残差診断：良いモデルの残差は白色ノイズ（自己相関なし）。Ljung-Box 検定（帰無＝白色）の p 値が大きく、残差 ACF がほぼ 0 なら合格。
当てはめ不足だと残差に自己相関が残り、Ljung-Box が棄却（p 値小）。標本内 AIC が良くても標本外で当たるとは限らない——必ずバックテスト（予測の評価指標と時系列CV）と併用します。

1. 次数選択：AIC と BIC

ACF/PACF の指紋（AR・MAモデル）で見当はつきますが、ARMA は両方とも尾を引いて切れないので（ARMA・ARIMAモデル）、最終的には情報量規準で $(p,q)$ を選びます。情報量規準は「当てはまり（尤度）」と「複雑さ（パラメータ数 $k$ ）」のトレードオフを1つの数字にします。

\mathrm{AIC}=-2\ln\hat L+2k,\qquad \mathrm{BIC}=-2\ln\hat L+k\ln n

ここで $\hat L$ は最大尤度、 $k$ は推定パラメータ数、 $n$ は標本数。第1項 $-2\ln\hat L$ は当てはまるほど小さく、第2項が複雑さへの罰則です。違いは罰則の重み：AIC は $2k$ 、BIC は $k\ln n$ 。 $\ln n>2$ （つまり $n>e^2\approx 7.4$ ）なら BIC の罰則が AIC より重いので、現実のデータでは BIC の方が複雑なモデルを嫌い、倹約的なモデルを選ぶ傾向があります。AIC は予測誤差最小を狙う設計でやや過剰適合寄り、BIC は真のモデルを一致推定する設計で倹約寄り、と覚えます。

コード①：真の次数を AIC/BIC で選べるか

真の次数 ARMA(2,1)（ $\phi_1=0.6,\phi_2=-0.3,\theta=0.5$ ）を仕込み、 $(p,q)$ を $0\sim3$ でグリッド探索して AIC・BIC を計算し、それぞれの最小を選びます。

import numpy as np
import warnings
warnings.simplefilter("ignore")
from statsmodels.tsa.arima.model import ARIMA
from statsmodels.tsa.arima_process import arma_generate_sample

# 真の ARMA(2,1): φ1=0.6, φ2=-0.3, θ=0.5
np.random.seed(31)
y = arma_generate_sample(np.array([1, -0.6, 0.3]), np.array([1, 0.5]),
                         nsample=2000, scale=1.0)

rows = []
for p in range(4):
    for q in range(4):
        try:
            r = ARIMA(y, order=(p, 0, q), trend="c").fit()
            rows.append((p, q, r.aic, r.bic))
        except Exception:
            pass

aic_best = min(rows, key=lambda x: x[2])
bic_best = min(rows, key=lambda x: x[3])
print(f"{'p':>2}{'q':>2}{'AIC':>12}{'BIC':>12}")
for p, q, aic, bic in rows:
    mark = ""
    if (p, q) == (aic_best[0], aic_best[1]): mark += " <-AIC最小"
    if (p, q) == (bic_best[0], bic_best[1]): mark += " <-BIC最小"
    print(f"{p:>2}{q:>2}{aic:>12.1f}{bic:>12.1f}{mark}")
print(f"真の次数=(2,1)  AIC選択=({aic_best[0]},{aic_best[1]})  BIC選択=({bic_best[0]},{bic_best[1]})")

出力：

 p q         AIC         BIC
 0 0      7317.5      7328.7
 0 1      5908.1      5924.9
 0 2      5682.5      5704.9
 0 3      5684.5      5712.5
 1 0      6360.1      6376.9
 1 1      5734.6      5757.0
 1 2      5684.5      5712.5
 1 3      5680.8      5714.4
 2 0      5825.9      5848.3
 2 1      5668.4      5696.4 <-AIC最小 <-BIC最小
 2 2      5670.1      5703.7
 2 3      5671.7      5710.9
 3 0      5708.9      5736.9
 3 1      5670.1      5703.7
 3 2      5672.0      5711.2
 3 3      5673.7      5718.5
真の次数=(2,1)  AIC選択=(2,1)  BIC選択=(2,1)

出力の意味： $n=2000$ と十分なデータがあるので、AIC・BIC とも真の $(2,1)$ を一発で当てました（最小 AIC $5668.4$ ・最小 BIC $5696.4$ ）。罰則の違いは「より大きいモデルへの嫌い方」に表れます——たとえば過剰な $(3,1)$ は AIC で $5670.1$ （ $(2,1)$ との差わずか $1.7$ ）なのに、BIC では $5703.7$ （差 $7.3$ ）。AIC は大きいモデルにわずかしか罰しないので過剰適合しやすく、BIC ははっきり罰して倹約的。データが少ないと AIC が無駄に大きい次数を選ぶことがあるため、両方見て、最後は標本外（バックテスト）で確かめます。

2. 残差診断：残差は白色ノイズか

次数を選んだら、それで十分かを残差で点検します。モデルが系列の依存をすべて吸収できていれば、残差 $\hat\varepsilon_t=y_t-\hat y_t$ には自己相関が残らず白色ノイズになるはずです。逆に残差にまだ自己相関があれば、取り残した構造がある＝当てはめ不足のサイン。これを定量化するのが Ljung-Box 検定です。

Q_{\mathrm{LB}}=n(n+2)\sum_{k=1}^{m}\frac{\hat\rho_k^2}{\,n-k\,}\ \overset{H_0}{\sim}\ \chi^2_{\,m-r}

$\hat\rho_k$ は残差のラグ $k$ 自己相関、 $m$ は検定するラグ数、 $r$ は推定した ARMA パラメータ数（自由度の調整）。帰無仮説 $H_0$ ＝残差は白色（自己相関なし）。残差の自己相関 $\hat\rho_k$ がどれも小さければ $Q_{\mathrm{LB}}$ は小さく p 値は大きい（棄却せず＝合格）。どこかに自己相関が残れば $Q_{\mathrm{LB}}$ が大きく p 値が小さい（棄却＝当てはめ不足）。

コード②：正しいモデル vs 当てはめ不足モデル

コード①と同じ ARMA(2,1) 系列に、正しい ARMA(2,1) とわざと小さくした AR(1) を当て、Ljung-Box の p 値と残差 ACF を対比します。

import numpy as np
import warnings
warnings.simplefilter("ignore")
from statsmodels.tsa.arima.model import ARIMA
from statsmodels.tsa.arima_process import arma_generate_sample
from statsmodels.stats.diagnostic import acorr_ljungbox
from statsmodels.tsa.stattools import acf

# 真の ARMA(2,1)（コード①と同じ生成）
np.random.seed(31)
y = arma_generate_sample(np.array([1, -0.6, 0.3]), np.array([1, 0.5]),
                         nsample=2000, scale=1.0)

# 正しいモデル ARMA(2,1) と 過小モデル AR(1)
res_ok = ARIMA(y, order=(2, 0, 1), trend="c").fit()
res_bad = ARIMA(y, order=(1, 0, 0), trend="c").fit()

lb_ok = acorr_ljungbox(res_ok.resid, lags=[10], return_df=True)
lb_bad = acorr_ljungbox(res_bad.resid, lags=[10], return_df=True)
p_ok = lb_ok["lb_pvalue"].iloc[0]
p_bad = lb_bad["lb_pvalue"].iloc[0]
print("Ljung-Box検定（帰無＝残差は白色ノイズ, lag=10）")
print(f"  正しい ARMA(2,1): p値={p_ok:.3f}  -> {'白色(棄却せず)' if p_ok>0.05 else '自己相関残る(棄却)'}")
print(f"  過小   AR(1)    : p値={p_bad:.3f}  -> {'白色(棄却せず)' if p_bad>0.05 else '自己相関残る(棄却)'}")

# 残差ACF（ラグ1,2）で当てはめ不足を可視化
a_ok = acf(res_ok.resid, nlags=5, fft=True)
a_bad = acf(res_bad.resid, nlags=5, fft=True)
print(f"残差ACF lag1: 正しい={a_ok[1]:.3f}  過小={a_bad[1]:.3f}")
print(f"残差ACF lag2: 正しい={a_ok[2]:.3f}  過小={a_bad[2]:.3f}")

出力：

Ljung-Box検定（帰無＝残差は白色ノイズ, lag=10）
  正しい ARMA(2,1): p値=0.956  -> 白色(棄却せず)
  過小   AR(1)    : p値=0.000  -> 自己相関残る(棄却)
残差ACF lag1: 正しい=0.002  過小=0.300
残差ACF lag2: 正しい=-0.002  過小=-0.324

出力の意味：正しい ARMA(2,1) の残差は Ljung-Box の p 値 $0.956$ ——帰無（白色）を棄却せず合格。残差 ACF もラグ1で $0.002$ 、ラグ2で $-0.002$ とほぼ 0 で、依存を取り尽くしています。一方、わざと小さくした AR(1) は p 値 $0.000$ ——棄却＝残差に自己相関が残る。残差 ACF はラグ1で $0.300$ 、ラグ2で $-0.324$ と明らかに構造が残っています（MA 項と2つ目の AR 項を取りこぼした分）。「残差が白色か」がモデル妥当性の核心で、Ljung-Box が棄却されたら次数を上げる（または別の構造を疑う）合図です。

3. コード③：残差 ACF を見比べる

数値だけでなく、残差 ACF を棒グラフで見ると「白色か／構造が残るか」が一目で分かります。

import numpy as np
import warnings
warnings.simplefilter("ignore")
import matplotlib.pyplot as plt
import japanize_matplotlib  # 日本語ラベル用
from statsmodels.tsa.arima.model import ARIMA
from statsmodels.tsa.arima_process import arma_generate_sample
from statsmodels.graphics.tsaplots import plot_acf

# 真の ARMA(2,1)（コード①と同じ生成）
np.random.seed(31)
y = arma_generate_sample(np.array([1, -0.6, 0.3]), np.array([1, 0.5]),
                         nsample=2000, scale=1.0)
res_ok = ARIMA(y, order=(2, 0, 1), trend="c").fit()
res_bad = ARIMA(y, order=(1, 0, 0), trend="c").fit()

fig, ax = plt.subplots(1, 2, figsize=(11, 4))
plot_acf(res_ok.resid, lags=15, ax=ax[0]);  ax[0].set_title("正しいARMA(2,1)の残差ACF（白色＝帯の中）")
plot_acf(res_bad.resid, lags=15, ax=ax[1]); ax[1].set_title("過小AR(1)の残差ACF（自己相関が残る）")
plt.tight_layout(); plt.show()

グラフの意味：左（正しいモデル）は、すべての棒が信頼帯（青い影）の中＝残差に自己相関なし＝白色。右（過小モデル）は、ラグ1・2の棒が帯を大きく超えており、取り残した依存がはっきり見えます。Ljung-Box の p 値（数値）と残差 ACF（図）はセットで使います——数値で合否、図でどのラグに問題があるかを読みます。

4. モデル選択の作法（まとめ）

同定：ACF/PACF（定常性と自己相関）と差分・季節差分（ランダムウォークと単位根・季節ARIMA（SARIMA））で候補の $(p,d,q)(P,D,Q)_s$ を絞る。
選択：AIC/BIC のグリッドで次数を選ぶ（AIC は予測寄り・BIC は倹約寄り）。
診断：残差の Ljung-Box（白色か）・残差 ACF・正規性を確認。棄却なら次数を見直す。
検証：標本外でバックテスト（ウォークフォワード・素朴予測比較、予測の評価指標と時系列CV）。標本内の AIC が良くても標本外で勝つとは限らない。

⚠️ よくある誤解

「AIC 最小が最良モデル」ではない：AIC は標本内の指標。過剰適合しやすく、標本外で劣ることがあります。必ずバックテストで標本外性能を確認します（予測の評価指標と時系列CV）。
「AIC と BIC は同じ結論」ではない：罰則が違うので、特に小標本や複雑なモデルでBIC はより小さい次数を選びます。食い違ったら、目的（予測精度なら AIC 寄り・解釈や倹約なら BIC 寄り）と標本外性能で決めます。
「グリッドで最小を選べば客観的」ではない：多数の候補を比べると、偶然 AIC が小さいモデルを拾う多重比較の罠があります。候補を絞る・標本外で再検証することが大事。
「Ljung-Box が通れば完璧」ではない：自己相関が無いだけで、**分散不均一（ボラティリティ・クラスタリング＝ARCH 効果、第6章）**や非正規性は別途点検が要ります。残差の二乗の自己相関も見ます。
「自由度調整は無視してよい」ではない：Ljung-Box の自由度は $m-r$ （ $r$ ＝推定 ARMA パラメータ数）。ラグ数 $m$ は推定パラメータ数より十分大きく取ります。