変数選択の規準

数理統計　第9章 §9.3　★★ 変数を増やすと $S_e$ は必ず減る。予測誤差の不偏推定として $C_p$・AIC・交差検証を導く

なぜ $S_e$ だけで選べないか

説明変数の候補が多いとき、どれを使うかを決めたい。$S_e$ や $R^2$ はこの目的に使えない。

$X$ に列を1本足したモデルの残差平方和を $S_e'$ とする。これは、もとのモデルの係数に新しい列の係数0を加えたものも候補に含む最小化の結果である。したがって $S_e'\le S_e$ が必ず成り立つ。

つまり $S_e$ は変数を足すほど下がり、$S_e$ 最小を規準にすると常に全部入りのモデルが選ばれる。しかし不要な変数を入れると推定量の分散が増え、将来のデータに対する予測は悪くなる。そこで規準は「手元のデータへの当てはまり」ではなく、「同じ計画で取り直したデータへの予測誤差」を推定するように作る。

Mallows の $C_p$

真の平均を $\bm{\mu}=E[\bm{Y}]$ とする。これは候補モデルの列空間に入っているとは限らない。候補モデルの射影行列を $H$（階数 $p$）とする。

同じ計画で独立に取り直したデータを $\bm{Y}^*$ とする（分布は $\bm{Y}$ と同じ）。その予測誤差の期待値を求める。

証明

予測の誤差を3項に分ける。

$$\bm{Y}^*-\hat{\bm{Y}}=(\bm{Y}^*-\bm{\mu})+(\bm{\mu}-H\bm{\mu})-H(\bm{Y}-\bm{\mu}).$$

第1項は $\bm{Y}$ と独立で、平均は $\bm{0}$ である。第3項も平均 $\bm{0}$、第2項は定数である。よって交差項の期待値はすべて0になり、

$$E\lVert \bm{Y}^*-\hat{\bm{Y}}\rVert ^2=n\sigma^2+\lVert (I-H)\bm{\mu}\rVert ^2+E\lVert H\bm{\varepsilon}\rVert ^2,\qquad E\lVert H\bm{\varepsilon}\rVert ^2=\operatorname{Tr}(H\sigma^2I H^\top)=p\sigma^2 .$$

$E[S_e]$ は $E[\bm{Y}^\top(I-H)\bm{Y}]=\sigma^2\operatorname{Tr}(I-H)+\bm{\mu}^\top(I-H)\bm{\mu}$ から出る。

予測誤差を $\sigma^2$ で割って $n$ を引いた量は、次のとおりである。

$$\frac{E[S_e]}{\sigma^2}-n+2p .$$

この $\sigma^2$ に、全変数を入れたモデルの $\hat\sigma^2_{\mathrm{full}}$ を代入して推定したもの

$$C_p=\frac{S_e}{\hat\sigma^2_{\mathrm{full}}}-n+2p\tag{9.11}$$

を Mallows の $C_p$ と呼び、これが最小のモデルを選ぶ。偏りのないモデル（$(I-H)\bm{\mu}=\bm{0}$）では $E[S_e]=(n-p)\sigma^2$ である。よって $C_p\approx p$ となる。$C_p$ が $p$ より大きく離れたモデルは、必要な変数が抜けていることを示す。

AIC

第6章 §6.4 で導入した赤池情報量規準は、$k$ 個の母数をもつモデルの最大対数尤度 $\ell(\hat\theta)$ を使って次で定義される。

$$\mathrm{AIC}=-2\ell(\hat\theta)+2k .$$

将来のデータに対する期待対数尤度を $\ell(\hat\theta)$ だけで推定すると、$k$ だけ上に偏る。AIC はこれを補正したものである。

正規線形モデルの母数は、係数ベクトルの $p$ 個と $\sigma^2$ の1個である。よって $k=p+1$ である。最大対数尤度は、命題 9.1 の式に $\hat\sigma^2_{\mathrm{ML}}=S_e/n$ を入れて

$$\ell(\hat{\bm{\beta}},\hat\sigma^2_{\mathrm{ML}})=-\frac n2\log\Bigl(2\pi\frac{S_e}{n}\Bigr)-\frac{1}{2(S_e/n)}S_e=-\frac n2\log\Bigl(2\pi\frac{S_e}{n}\Bigr)-\frac n2 .$$

よって

$$\mathrm{AIC}=n\log\frac{S_e}{n}+n(\log2\pi+1)+2(p+1) .\tag{9.12}$$

同じデータでモデルを比べるときは $n(\log 2\pi+1)$ は共通なので、$n\log(S_e/n)+2(p+1)$ を比べれば足りる。$\log n$ を罰則に使う BIC は $-2\ell(\hat\theta)+k\log n$ である。$n\ge8$ では AIC より強く変数の追加を罰する。

交差検証

観測 $i$ を除いた $n-1$ 個で当てはめる。その当てはめで観測 $i$ を予測した誤差を、次のように書く。

$$e_{(i)}=Y_i-\bm{x}_i^\top\hat{\bm{\beta}}_{(i)} .$$

その2乗和 $\mathrm{PRESS}=\sum e_{(i)}^2$ を規準にするのが、1つ抜き交差検証である。$n$ 回当てはめ直す必要はなく、次の公式で1回の当てはめから計算できる。

証明

$A=X^\top X$ とおく。第 $i$ 行を除くと次のようになる。

$$X_{(i)}^\top X_{(i)}=A-\bm{x}_i\bm{x}_i^\top,\qquad X_{(i)}^\top\bm{Y}_{(i)}=X^\top\bm{Y}-\bm{x}_iY_i .$$

Sherman–Morrison の公式

$$(A-\bm{x}\bm{x}^\top)^{-1}=A^{-1}+\frac{A^{-1}\bm{x}\bm{x}^\top A^{-1}}{1-\bm{x}^\top A^{-1}\bm{x}}$$

を使うと

$$\bm{x}_i^\top(A-\bm{x}_i\bm{x}_i^\top)^{-1}=\bm{x}_i^\top A^{-1}+\frac{h_{ii}\,\bm{x}_i^\top A^{-1}}{1-h_{ii}}=\frac{\bm{x}_i^\top A^{-1}}{1-h_{ii}} .$$

したがって

$$\bm{x}_i^\top\hat{\bm{\beta}}_{(i)}=\frac{\bm{x}_i^\top A^{-1}(X^\top\bm{Y}-\bm{x}_iY_i)}{1-h_{ii}}=\frac{\hat Y_i-h_{ii}Y_i}{1-h_{ii}},$$

$$e_{(i)}=Y_i-\frac{\hat Y_i-h_{ii}Y_i}{1-h_{ii}}=\frac{Y_i-\hat Y_i}{1-h_{ii}}=\frac{e_i}{1-h_{ii}} .$$

$h_{ii}$ の平均は $p/n$ なので、次の近似が成り立つ。

$$\mathrm{PRESS}\approx\frac{S_e}{(1-p/n)^2}\approx S_e\Bigl(1+\frac{2p}{n}\Bigr).$$

$\sigma^2$ を $S_e/n$ で置き換えてみる。すると $C_p$ と同じ罰則 $2p$ が現れる。3つの規準は、予測誤差の推定という同じ考えを別の近似で実現したものである。

確認テスト

オリジナル

1. 9.3-1 自由度調整済み決定係数

$$\bar R^2=1-\frac{S_e/(n-p)}{S_{YY}/(n-1)}$$

を考える。説明変数を1つ足したとき、$\bar R^2$ が増えるための必要十分条件は、足した変数の係数の $F$ 統計量が1より大きいことである。これを示せ（$F=t^2$ である）。2. 9.3-2 $n=30$ で、入れ子の3つのモデル（$p=2,3,4$）の残差平方和が 120、96、92 であった。定数を除いた AIC で選ばれるモデルを答えよ。3. 9.3-3 説明変数を1つ足して $p$ から $p+1$ にしたとき、AIC（式 (9.12)）が下がる条件を求めたい。足した変数の $F$ 統計量 $F=(S_e-S_e')/\{S_e'/(n-p-1)\}$ で表せ。また $n$ が大きいとき、この条件がおよそ $F>2$ になることを示せ。4. 9.3-4 確認テスト 9.2-1 のてこ比を使い、単回帰で $h_{ii}$ が最大になるのはどの観測か答えよ。また $h_{ii}\to1$ のとき1つ抜き残差 $e_{(i)}$ がどうなるかを、命題 9.9 から説明せよ。