回帰診断法：1点で傾きが半分になる仕組みと、てこ比の正体【第17回】

はじめに

第17章は回帰診断法です。前回（第16回・重回帰分析）で「最小二乗法は射影である」「偏回帰係数の正体は残差回帰である」といった話を扱いました。あの回はモデルを作って係数を解釈する話でした。

今回は逆向きです。そのモデルの前提が本当に成り立っているのかを確かめる回です。

前回の最後にガウス・マルコフの定理が出てきました。「誤差が期待値ゼロ・等分散・無相関なら、最小二乗推定量は不偏推定量の中で最も分散が小さい」という定理です。よく考えるとこれは条件付きの保証で、条件が崩れたら何が起きるかは何も言っていません。今回はその「崩れたとき」を全部見ていきます。

そして今回、私の中でいちばん大きかった疑問はこれでした。

最小二乗法は全部の点の残差を2乗して足して最小化するのに、なぜ1点だけが結果を支配できるのか。 平等に扱っているように見えるのに、実際には1点で傾きが半分になる。この矛盾がずっと引っかかっていました。

答えは意外なところにありました。「二乗しているから、まさにそうなる」のです。最小二乗法は平等ではありません。式の上でそれが見えます。

この記事で使う言葉

先に5つ整理します。この記事はこの5語で回ります。

てこ比（leverage）：$h_{ii}$ と書きます。ある観測が説明変数の空間で「みんなの中心からどれだけ離れているか」の指標。てこ（レバー）の柄の長さにあたります。重要なのは$y$ の値を一切見ずに計算できることです。

外れ値（outlier）：$y$ の方向にズレている観測。てこ比とは別の概念で、こちらは残差を見て判定します。

影響力（influence）：その観測が実際に結果をどれだけ動かしたか。上の2つの掛け算で決まります。代表的な指標が次のクック距離です。

クック距離（Cook's distance）：$D_i$ と書きます。「その1点を抜いたら予測値全体がどれだけ動くか」を測った量。影響力の指標です。

不均一分散（heteroscedasticity）：誤差の分散が観測ごとに違う状態。反対語は等分散（homoscedasticity）。「hetero＝異なる」「skedasis＝散らばり」というギリシャ語由来で、字面のまま「散らばりが異なる」です。

略語も先に開いておきます。OLS（Ordinary Least Squares＝通常の最小二乗法）、WLS（Weighted Least Squares＝加重最小二乗法）、SE（Standard Error＝標準誤差）、DW（Durbin-Watson＝ダービン・ワトソン比）、BP検定（Breusch-Pagan test＝ブロイシュ・ペーガン検定）、RSS（Residual Sum of Squares＝残差平方和）、TSS（Total Sum of Squares＝全平方和）、RMSE（Root Mean Squared Error＝二乗平均平方根誤差）。

なお $p$ という記号が2つの意味で使われる問題は前回と同じです。この記事ではパラメータ数を $p$、有意確率は「p値」と書き分けます。

TL;DR

• 最小二乗法は平等ではない。 傾きは $\hat\beta_1=\sum w_i y_i$ という $y$ の加重和で、重みは $w_i=(x_i-\bar x)/S_{xx}$。実測で $\sum w_i y_i$ と OLS の傾きが 1.07709520 で一致
• $x_i=\bar x$ の点は傾きに1票も持っていない（重み $-6\times10^{-18}$）。だから $y$ を14動かしても傾きは変化 $-0.000$
• 1点が回帰を壊すには「$x$ が遠い」と「$y$ がズレている」の両方が必要。 $D_i=\frac{r_i^2}{p}\cdot\frac{h_{ii}}{1-h_{ii}}$ が掛け算だから、どちらかが0なら積は0
• $x$ が遠く $y$ もズレた1点で傾きが 1.024 → 0.350 に崩壊（$D_i=9.115$、$R^2$ は 0.945 → 0.181）
• てこ比は $y$ を見ずに計算できる（$X$ だけの関数）。$h_{ii}=\partial\hat y_i/\partial y_i$ ＝自分の予測への発言力で、5点すべてで実測の傾きが $h_{ii}$ と $10^{-16}$ 精度で一致
• 重回帰では「単独では平凡なのに組み合わせが異常」な点が生まれる。 相関0.639のとき、重心から同じ距離3.111でも $h$ が 0.465 対 0.213 と2.2倍違う
• てこ比が高い点は生の残差が小さく出る（$V[e_i]=\sigma^2(1-h_{ii})$）。$h=0.42$ の点の残差は 0.030 しかない
• 外れ値は自分でものさしを伸ばして自分を隠す（マスキング）。$\hat\sigma$ が 0.925 → 2.768 と3倍に膨らみ、標準化残差は $\sqrt{n-p}=4.243$ で飽和して止まる
• だから外れ値判定には外部（スチューデント化残差）を使う。 こちらは 11.79 まで伸び、$t(n-p-1)$ に厳密に従う（97.5%点が実測 2.2490 対 理論 2.2622）
• 不均一分散が壊すのは標準誤差だけ。 推定値は不偏のまま（平均 1.500377、真値 1.5）だが、SE は本当のSDの 0.944倍で被覆率 93.74%、偽陽性 6.55%
• ただし「必ず過小評価」は誤り。 分散が中央で大きい形だと比が 2.262 で過大評価になる。鍵はてこ比との相関
• $4/n$ は閾値ではなく目印。 完全に正常なデータで $n$=200 のとき誤検出率 100%（$E[D_i]\approx1/(n-p)$ なので $4/n$ は常に上位数%の線。原因は多重比較）
• 削除は区間を 0.721倍に縮めるが被覆率を下げる（95.73% → 91.07%、繰り返すと 80.55%）。一方で点推定のRMSEは改善する（0.102 → 0.084）
• 捨てるより重みを下げるほうが強い。 Huber ロバスト回帰は RMSE 0.080 で削除より良く、標本を1つも捨てない。ただし守るのは $y$ 方向だけ（汚染点が $x$=60 なら OLS と同様に壊れる）
• DW が小さい原因は自己相関でないことがある。 誤差が完全に独立な曲線データで DW=0.368、$x^2$ を入れたら 2.092 に戻った
• 自己相関は不均一分散よりずっと危険。 $\rho$=0.9 で被覆率 34.03%。ただしこれも向きは $x$ の並びで変わる
• アンスコムの4データは係数・$R^2$・$\hat\sigma$・残差SDが全部一致するのに中身が別物。 うち1つは $h_{ii}$ がちょうど 1.000
• おまけ：残差の平均は常に厳密にゼロなので、診断の材料にならない

なぜ二乗の最小化なのに1点が支配できるのか

順序を変えます。教科書だと診断量の定義から入りますが、最初にこの疑問を片付けたほうが後が全部楽になると気づいたので、ここから始めます。

傾きは $y$ の加重和だった

単回帰の傾きの式を変形します。

ここで $S_{xx}=\sum_i(x_i-\bar x)^2$ です。$\sum_i(x_i-\bar x)\bar y = \bar y\sum_i(x_i-\bar x) = 0$ なので分子の $\bar y$ が消えて、この形になります。

つまり傾きは $y$ のただの加重和です。重み $w_i$ は $x$ だけで決まっていて、$y$ を見る前に確定しています。

実際に確かめました。

下2行は重みが満たす制約です。$\sum w_i=0$ かつ $\sum w_i x_i=1$ で、この2本が「不偏推定量になるための条件」にあたります。

そして重みの中身を見ると、疑問が一気に解けます。

(a) が核心です。 棒グラフは「どの点が傾きに何票入れているか」を表しています。$x_i=\bar x$ のところで棒がゼロになっているのが見えるでしょうか。実測値は $-6.22\times10^{-18}$、つまり真ん中の点は傾きについて1票も持っていません。端の点の票は重く、$x$ が $\bar x$ から離れるほど比例して重くなります（端から2番目 0.03000 対 中央から2番目 0.01000 でちょうど3倍）。

(b) は「なぜ二乗なのか」の答えです。 支点から距離 $d$ の点を考えると、回転への効きは $d$ に比例する一方、二乗損失の減り方は $d^2$ に比例します。だから最小化アルゴリズムは遠い点を優先して合わせに行く。

もし損失が絶対値（$|e_i|$）なら、残差の大きさに対する効きは一定になります（$y$ が10ズレても100ズレても効きが同じ）。これが後で出てくるロバスト回帰の発想の源です。

ただし注意が必要で、絶対値損失にしても $x$ 方向のてこ比は消えません。勾配を書くと明らかです。

どちらにも $(x_i-\bar x)$ が残っています。 絶対値損失が一定にするのは「$y$ 方向のズレに対する効き」だけで、「$x$ の遠さに対する効き」はそのままです。この事実は後でロバスト回帰の限界として戻ってきます。

(c) は別の角度からの説明です。 最小二乗法は次の2本を必ず満たします。

2本目が問題です。これは「残差を $x$ で重み付けした和がゼロ」という制約なので、$x_i$ が大きい点の残差には係数 $x_i$ が掛かって効きます。少し残るだけで和を崩すため、つじつまを合わせるために直線がその点へ引っぱられる。実データで確認すると、どちらのデータでも $\sum x_i e_i$ は $10^{-14}$ 台、つまり例外なく厳密に成立していました。

ここが今回の土台です。

最小二乗法が平等なのは残差の扱いだけで、係数への寄与は最初から不平等です。 $\hat\beta_1=\sum w_i y_i$、$w_i=(x_i-\bar x)/S_{xx}$ という式がそれを示しています。

• $x_i$ が $\bar x$ に近い点 … $w_i\approx0$。$y$ がどれだけ暴れても傾きは動かない
• $x_i$ が遠い点 … $w_i$ が大きい。$y$ が少し動くと傾きが動く

てこ比が高い点とは「$w_i$ が大きい点」＝1票の重みが大きい点のことです。 そして $y$ 方向のズレが「その1票の中身」。影響力＝票の重み × 票の中身。

「二乗誤差の最小化」という手続きが、勝手にこの重み付けを作り出している。ここが腑に落ちると、以降の診断量が全部「この不平等を測る道具」として読めるようになります。

重回帰でも同じことが起きる

「単回帰の話では？」と思ったので、$p=3$（切片＋説明変数2つ）で実演しました。上の図の (d)(e)(f) です。$x_1$ と $x_2$ の相関は 0.639 に設定しています。

1点を $(x_1,x_2)=(-2.2,\ 2.2)$ に置いて、その点の $y$ だけを下げていきます。

$\hat\beta_1$ が 2.05 から 3.34、$\hat\beta_2$ が −1.39 から −2.96（真値はそれぞれ 2.0 と −1.5）。1点で両方の係数が2倍近くまで動きました。 $h_{ii}$ が 0.465 のまま一定なのは、てこ比が $X$ だけで決まるので $y$ をいくら動かしても変わらないためです。図(f)で全点のクック距離を並べると、この点の 9.77 が他の40点の最大 0.114 の86倍でした。

そして重回帰のほうが厄介です。 (e) を見てください。$x_1$ と $x_2$ に相関があると、てこ比の等高線は傾いた楕円になります。

同じ距離なのに $h$ が2.2倍違います。 前者は「$x_1$ が小さいのに $x_2$ が大きい」、つまり相関 0.639 という全体の傾向に逆らっている点だからです。

ここが実務で怖いところでした。$x_1$ も $x_2$ も単独では $\pm2.2$ という値で、標準偏差1の変数ならよくある値です。単独では平凡なのに、組み合わせとして異常だから高てこ比になる。つまり1変数ずつヒストグラムを描いて外れ値チェックをしても絶対に見つかりません。$h_{ii}$ を計算して初めて見えます。これが $h_{ii}$ という量をわざわざ使う実用上の理由でした。

なお、てこ比が「重心からの楕円距離」であることは式でも確認できます（後で出てくる図の (d) が等高線の全体像です）。マハラノビス距離 $d_M$ を使うと

で、実測との最大誤差は $2.8\times10^{-16}$ でした（$d_M$ は共分散行列を $n-1$ で割って測った場合。$n$ で割る流儀なら分母も $n$ になります）。楕円が傾く理由は、マハラノビス距離が共分散行列の逆行列で測る距離だからです。

1点を動かすと直線はどう動くか

ここからが今回いちばん見たかったものです。20点を固定して、21番目の点だけを動かします。動かし方を3通り変えると、まったく違う挙動になります。

A：$x$ は真ん中のまま、$y$ を大きく動かす

$x=5$（ちょうど $\bar x$）に置いたまま、$y$ を14上へ動かします。右側の数値に注目してください。

$y$ を14も動かしたのに、傾きは 1.046 から1ミリも動きません。 動くのは切片だけ（1.749 → 2.415）。スチューデント化残差は +15.48 という異常な値まで行き、$R^2$ は 0.902 から 0.391 へ落ちるのに、傾きだけが無傷です。

理由はもう分かります。$x_i=\bar x$ の点は $w_i=0$、つまり傾きについて1票も持っていない。てこの支点の真上を押しても回転モーメントが生じないのと同じです。

B：$y$ は直線上のまま、$x$ を遠ざける

今度は $y$ を直線の延長上に置いたまま、$x$ を 9.5 から 16 へ遠ざけます。

てこ比が 0.182 から 0.518 までぐんぐん上がります。危険な位置に移動しているわけです。ところが傾きの変化はごくわずか（1.040 → 1.024）、クック距離も 0.081 止まり。

しかも $R^2$ は 0.913 から 0.945 に「良くなって」しまいます。$x$ の範囲が広がって $S_{xx}$ が増えたためで、あてはまりの指標としては改善に見える。てこ比が高いだけでは何も起きていないというのがこの図の主張です。

C：$x$ が遠く、$y$ もズレる

本命です。$x=16$ に置いたまま（$h_{ii}=0.518$ に固定）、$y$ を15下へ動かします。

傾きが 1.024 から 0.350 へ崩壊しました（変化 −0.696）。クック距離は 9.115 で、よく使われる閾値1をはるかに超えています。$R^2$ は 0.945 から 0.181。

Bと違うのは $y$ 方向にもズレたことだけ、Aと違うのは $x$ が遠いことだけです。

3つのアニメが言っていること。

1点が回帰を壊すには、「$x$ が遠い」と「$y$ がズレている」の両方が必要です。

• A … $y$ がズレているだけ → 傾きは無傷
• B … $x$ が遠いだけ → 何も起きない（むしろ $R^2$ は改善）
• C … 両方そろって初めて崩壊

式で見ると理由が明快です。

$$D_i = \frac{r_i^2}{p}\cdot\frac{h_{ii}}{1-h_{ii}}$$

$y$ 方向のズレ（$r_i$）と $x$ 方向の遠さ（$h_{ii}$）の掛け算になっています。 掛け算だから、どちらかが0なら積は0。これがCだけが壊れる理由です。

外れ値・てこ比・影響力を2×2で整理する

この3つが別物だという話を、同じデータで並べます。21番目の点の置き場所を4通り変えただけです（$h_{ii}$ の下限が $1/21=0.048$ になっているのがその証拠です）。

大事な工夫が1つあります。$y$ の位置を「その点を除いた直線」からの上下で指定しました。こうしないとてこ比だけを孤立させられません（最初は「真の直線」から指定して失敗しました。それだと③に残差が残ってしまい、影響力がゼロになりません）。

灰色の破線が「問題の点を入れない直線」、赤が「入れた直線」です。2本が重なっていれば影響ゼロと読めます。

②と③が対比のポイントです。

② は残差が +9.40 という明確な外れ値なのに、傾きの変化は −0.000。動いたのは切片だけ（+0.429）です。クック距離 0.395 は「中くらい」で、閾値1には届きません。

③ はてこ比 0.560 という危険な位置にいるのに、クック距離が 0.000、傾きの変化も +0.000。図でも2本の線が完全に重なっています。

④ は③と同じ位置から $y$ を12下げただけで、クック距離が 9.681、傾きが 1.123 から 0.585 へ半減しました。

3つの言葉の定義（これで区別が付きます）。

• 外れ値 … $y$ 方向にズレている。残差を見て判定する。データを取ったあとで初めて分かる
• てこ比 … $x$ 方向に遠い。$y$ を一切見ずに計算できる（$h_{ii}$ は $X$ だけの関数）。つまり「まだ何も起きていないが、起きたら大きい」という潜在的な危険度
• 影響力 … 実際に結果を動かした量。上の2つの掛け算

たとえで言うと、てこ比はハンマーの柄の長さ、外れ値は振り下ろす力、影響力は実際に開いた穴です。 柄が長くても振らなければ穴は開かない（③）。思い切り振っても柄が短ければ穴は小さい（②）。

てこ比が「$y$ を見ずに計算できる」というのは、実験計画の話とも繋がります。どこで観測するかを自分で決められる場合、てこ比はデータを取る前に設計できるのです。前回、説明変数がよく散らばっていると標準誤差が下がるという話をしましたが（単回帰なら $\mathrm{SE}(\hat\beta_1)=\hat\sigma/\sqrt{S_{xx}}$ で、分母が「てこの長さ」にあたります）、あれは「てこ比の高い点をあえて作ると精度が上がる」という意味でもありました。精度を上げる工夫と、1点に支配されるリスクは同じ源から来ているわけです。

てこ比とハット行列の幾何

てこ比の意味をもう少し掘ります。名前の由来から入ります。

(a) ハット行列という名前の由来。 $H = X(X^\top X)^{-1}X^\top$ を使うと $\hat y = Hy$ と書けます。$y$ に帽子（hat）をかぶせる行列なのでハット行列です。やっていることは射影、つまり影を落とす操作だけです。$X$ の列が張る平面が「作れる予測値の全体」で、$H$ は $y$ をそこに垂直に落とす。落とした先が $\hat y$、落ちた差が残差 $e=(I-H)y$ です。

$n=3$ の実例で確認しました。$y=(2,1,5)$ に対して $\hat y=(1.1429,\ 2.2857,\ 4.5714)$、残差 $(0.8571,\ -1.2857,\ 0.4286)$。

ピタゴラスの定理が成立しています。残差と列空間の直交も $e^\top\mathbf{1}=1.6\times10^{-15}$、$e^\top x=2.9\times10^{-15}$ で確認できました。前回の「最小二乗法は射影」がそのまま再登場した形です。

(b) てこ比は $\bar x$ からの距離で決まる。 単回帰なら閉じた式が書けます。

図の×印がこの式の値で、実測との最大誤差は $1.9\times10^{-16}$ でした。ここから3つの性質が読めます。

合計が $p$ になるのは前回の $\mathrm{tr}(I-H)=n-p$ と同じ話です。だから平均は $p/n$ で、目安の $2p/n$ や $3p/n$ は「平均の2〜3倍」という発想から来ています。

(c) ここがてこ比の意味そのものです。 $y_i$ を1だけ動かすと $\hat y_i$ はどれだけ動くか、を実験しました。

完全に一致しました。つまり

で、自分の予測値を自分でどれだけ決めているかの割合です。$h_{ii}=1$ なら回帰線はその点を必ず通る（他の点の意見が一切入らない）。$h_{ii}=1/n$ なら、他の点と平等に $1/n$ しか発言していない。これが「てこ比」という訳語の意味でした。

(e)(f) 実務上の落とし穴。 残差の分散は一定ではありません。

6000回のシミュレーションで、生の残差の実測SDが理論の $\sigma\sqrt{1-h_{ii}}$ と一致することを確認しました（最大誤差 0.0166）。標準化すると 0.986〜1.018 で平らになります。

これが意味するのは、てこ比が高い点は生の残差が小さく出るということです。実測では $h=0.425$ の点の残差が 0.030 しかありませんでした。$1/\sqrt{1-h}=1.318$ 倍に補正して初めて他の点と比較できます。

生の残差プロットだけを見ていると、いちばん危険な点が「よくあてはまっている」ように見える。 これがてこ比を別途チェックすべき理由です。

3つの残差と、クック距離の定義

診断量の定義をまとめます。前提として $e_i=y_i-\hat y_i$（生の残差）、$p$ はパラメータ数、$\hat\sigma^2=\mathrm{RSS}/(n-p)$ です。

① てこ比

$y$ を使わないのが特徴です。目安は $2p/n$（甘め）または $3p/n$（厳しめ）。

② 標準化残差（内部スチューデント化残差）

$V[e_i]=\sigma^2(1-h_{ii})$ で割って点どうしを比較できるようにしたものです。$t$ 分布には従いません。分母の $\hat\sigma$ の中に分子の $e_i$ 自身が入っているので、分子と分母が独立でないためです。

では何に従うかというと、$r_i^2/(n-p)$ が $\mathrm{Beta}\left(\frac12, \frac{n-p-1}{2}\right)$ に従います。ベータ分布は $[0,1]$ に台を持つので、$r_i^2 \le n-p$、つまり $|r_i|\le\sqrt{n-p}$ という上限がここから出ます。$n$=12、$p$=2 で6万回まわしたときの $\max|r_i|$ は 3.0723 で、上限 3.1623 に届きませんでした。

③ スチューデント化残差（外部・削除残差）

$\hat\sigma_{(i)}$ はその点 $i$ を抜いて計算した誤差の標準偏差です。こうすると分子と分母が独立になり、ちょうど $t(n-p-1)$ に従います。だから「この点は外れ値か」を検定できます。

右側の等式が実用上ありがたいところで、$n$ 回の再フィットをしなくても $r_i$ から変換できることを意味します。

④ クック距離

左が定義（$\hat y_{j(i)}$ は点 $i$ を抜いて再推定した予測値）、右が計算式です。「その1点を抜いたら、全部の予測値がどれだけ動くか」を $\hat\sigma^2$ で割って無次元化したものです。

自作モジュールで両者が一致することを確認しました（最大誤差 $1.2\times10^{-15}$）。定義どおり $n$ 回再フィットした結果と、右の式が完全に同じ値になります。

4つの関係を1行で。

$h_{ii}$ は $x$ だけの量、$t_i$ は $y$ のズレを $t$ 分布に載せた量、$D_i$ はその2つの掛け算。 $r_i$ は $t_i$ を計算する途中に出てくる中間量、という位置づけです。

なぜ「自分を抜く」というひと手間が必要なのか

定義だけ見ると③の $\hat\sigma_{(i)}$ は面倒に思えます。実験して理由を確かめました。

(a) 1点を上へずらしていくと、内部と外部で挙動が正反対になります。

青（内部）は 4.005 で飽和して止まります。この実験は $n$=20・$p$=2 なので、理論上限 $\sqrt{n-p}=4.243$ を超えられないためです（上限は $n$ で変わります。次の (c) は $n$=12 なので 3.162 が壁になります）。赤（外部）は 11.79 までまっすぐ伸びます。

(b) 原因はこれです。 $\hat\sigma$ が 0.925 から 2.768 へ、つまり3倍に膨らんでいます。一方 $\hat\sigma_{(i)}$ は 0.940 のまま不動。

つまり外れ値は自分自身でものさしを伸ばして、自分を平凡に見せてしまう。これをマスキング（自己隠蔽）と呼びます。分子（残差）が伸びても分母（ものさし）が同じ倍率で伸びるので、比が伸びないわけです。

(c) 分布も確認しました。 $n=12$、$p=2$ で2万回まわして分位点を照合します。

外部は $t(9)$ にきれいに従っています。内部の実測範囲は $[-2.895,\ 3.020]$ で、$\sqrt{n-p}=3.162$ の壁で切られていて裾が足りません。

マスキングが実務で意味すること。

$r_i$ には $\sqrt{n-p}$ という上限があるので、$n$ が小さいほど天井が低い。 $n$=12、$p$=2 なら $|r_i|$ は 3.162 を絶対に超えないので、 「$|r|>3$ を外れ値とする」という基準はこのサイズだとほぼ機能しません。

→ 外れ値の判定には必ず外部（スチューデント化残差）を使う。 これが「自分を抜く」ひと手間を払う理由です。

(d)(e) クック距離の2つの因子。 等高線で見ると、①〜④が前半の2×2にそのまま対応します。そして $h/(1-h)$ という形が非線形であることが効いてきます。

$h=0.5$ で1倍、0.9 で9倍、0.99 で99倍。$h\to1$ で発散します。だからてこ比が 0.8 を超えたあたりから急に危険になるわけです。

(f) そして3つの指標は別の点を指しました。 これが「3つとも見る必要がある」ことの実演です。

$h$ 最大は $i$=25、$|t_i|$ 最大は $i$=6、$D_i$ 最大は $i$=23。3つとも違う点です。 どれか1つだけ見ていると取りこぼします。

回帰の4つの仮定と、確認する方法

ここから後半です。診断の対象になる仮定を先に並べます。

この表の3列目が今回の要点です。 崩れ方が2種類に分かれます。線形性の破れは $\hat\beta$ を偏らせるのでモデルを直すしかない。一方、独立性と等分散性の破れは $\hat\beta$ を偏らせずSEだけを壊すので、SEの計算法を替えれば済みます。

正規性だけ性質が違って、$n$ が大きければ中心極限定理で救われます。第13回で検定の頑健性を扱ったときと同じ構図です。「正規性がいちばん重要そう」と思いがちですが、実は4つの中でいちばん優先度が低いというのが今回の学びでした。

残差プロットの判断表

過去問では図の読み取りが問われやすいので、6パターンを同じ形式で並べました。すべて横軸に予測値 $\hat y$、縦軸に標準化残差です。オレンジの線は移動平均で、傾向を見やすくするために引きました。

同じデータの診断量も出しました。

この表を作っていて2つ気づきました。

1つめ。③の $R^2$ が 0.641 でそれなりに高く、⑥の $R^2$ は 0.896 で①の 0.893 とほぼ同じです。$R^2$ という数字だけ見ていたら、③も⑥も見逃します。

2つめが自分で驚いたところで、③の DW が 0.31 になっています。③は誤差を完全に独立に生成しているのに、自己相関の指標が「強い正の自己相関」を示しました。

理由を考えると納得できます。$x$ でソートされたデータに曲線があると、直線をあてはめた残差は「左は正、中央は負、右は正」となめらかに変化します。すると隣同士が似るので、自己相関として検出されてしまう。実際に $x^2$ を入れて曲線を吸収したら DW が 2.09 に戻りました（後述）。

DW が小さいときは、まず線形性を疑う。 自己相関と決めつけて時系列モデルに走ると誤診します。手順は 「①残差 vs $\hat y$ で曲線がないか確認 → ②なければ自己相関を疑う」。

横軸に何を取るかで見えるものが変わる

なぜ横軸に $\hat y$ を取るのかというと、残差と予測値は直交している（前回の話）ので、仮定が正しければ無相関に見えるはずだからです。構造が見えたらそれは仮定の破れです。$y$ を横軸にすると残差と相関してしまうので間違いです。

不均一分散：壊れるのは推定値か標準誤差か

ここは区別が大事だと思ったので、実測で切り分けました。$\sigma(x)=0.30x$（0.3 から 3.3 まで11倍の差）というかなり強い不均一分散を作り、8000回まわします。

まず結論。推定値はまったく壊れていません。

壊れているのは標準誤差です。

その結果どうなるか。

答え：壊れるのは標準誤差です。推定値は不偏のまま。

式で追えます。$\hat\beta = \beta + (X^\top X)^{-1}X^\top\varepsilon$ なので、 $E[\varepsilon]=0$ でありさえすれば分散がどうなっていようと $E[\hat\beta]=\beta$。 不偏性は「誤差の平均がゼロ」しか使っていません。

一方 $V[\hat\beta]=\sigma^2(X^\top X)^{-1}$ は「全部の点で $\sigma^2$ が同じ」を使って導いた式なので、 そこが崩れると使えません。正しくは

$$V[\hat\beta]=(X^\top X)^{-1}X^\top \Omega X (X^\top X)^{-1}$$

（$\Omega$ は誤差の分散共分散行列）。この形をサンドウィッチ型と呼び、 $\Omega$ の対角だけを残差の2乗で置き換えたのがロバスト標準誤差です。 逆に言えば、HC系は $\Omega$ が対角（＝不均一分散だけ）を前提にしているので 自己相関には効きません（後述の自己相関には HAC＝Newey-West 標準誤差を使います）。

→ 実務的には「係数の値は信じていい、p値と信頼区間は信じるな」という判断になります。

「必ず過小評価」は間違いだった

教科書でよく「不均一分散はSEを過小評価する」と書かれます。私もそう覚えていたのですが、実測すると向きが変わりました。

鍵はてこ比との相関でした。 てこ比が高い点（$x$ が端）の分散が大きいとSEは過小評価になり、逆にてこ比が低い点（中央）の分散が大きいと過大評価になります。

考えてみれば当然で、$\hat\sigma^2$ は全部の残差を平均した「代表値」です。代表値が、効きの強い点（高てこ比）の実際の分散より小さければSEは小さすぎ、大きければ大きすぎになる。「不均一分散があると危険」は正しいですが、「必ず甘くなる」は誤りでした。

3つの対処法の使い分け

WLSが効率を上げる理由は、ばらつきの小さい点（信頼できる観測）を重く、大きい点を軽く扱うからです。分散が 55.5% 減るというのは、同じ精度を得るのに必要な標本が半分以下になるということです。

ただし WLS が最良になるのは重み $1/\sigma_i^2$ が（比だけでも）分かっているときです。分散の形をデータから推定して重みにすると、その推定の不確実性が入るので保証は消えます。だから第一選択がロバスト標準誤差なのでした。

ここで前回のガウス・マルコフの定理と繋がります。あの定理は「等分散のもとでOLSが最良」と言っていました。つまり裏を返せば、等分散でなければOLSより良いものがあるという宣言でもあったわけです。WLSがまさにそれでした。

対数変換すると解釈はどう変わるか

変数変換は不均一分散の対処として出てきますが、係数の意味が変わるのが厄介なところです。ここを数値で確かめました。

真のモデルを $y = 3.0 \cdot x^{1.4} \cdot e^{\varepsilon}$、つまり誤差が掛け算で入る形にします。売上や PV のように「何割増える」という性質の量はこうなりがちです。

(a)(b) 変換が効く仕組み。 生のスケールだと右へ行くほど散らばりが大きく、そのまま回帰すると BP検定のp値が 0.0029 で不均一分散が検出されます。両辺の対数を取ると

で、直線かつ等分散になります。推定値は傾き 1.3976（真値 1.4）、切片 1.1191（真値 $\log 3 = 1.0986$）、BP検定のp値は 0.546 で等分散を棄却しません。

掛け算の誤差は $\log$ を取ると足し算になる。 これが変換の原理でした。逆に言えば、誤差が最初から足し算で入っているデータに対数変換をしても意味がない（むしろ壊す）ことになります。

(c) 両対数の傾きは弾性値です。

1%の変化なら近似がほぼ完璧ですが、100%の変化だと厳密 163.5% 対 近似 139.8% で大きく外れます。「1%増えると $\beta_1$% 増える」は小さな変化に限った近似でした。

(d) 片対数（$y$ だけ対数）の場合。

$\beta_1=0.2001$ なので「$x$ が1増えると $y$ が約20%増える」と読みますが、$\Delta x=5$ では厳密 171.9% 対 近似 100.0% とまったく合いません。複利で効くので当然です。

4つのモデルを整理します。

対数変換の副作用（3つとも実務で踏みます）。

① $R^2$ を変換前後で比べてはいけない。 上の例では 0.772 から 0.903 に上がりましたが、 $y$ と $\log y$ では被説明変数が別物なので TSS の意味が違います。比較は無意味です。

② 予測値を戻すときに偏りが出る。 $e^{\widehat{\log y}}$ は $E[y]$ ではなく中央値の推定になります。 理由は2段階です。(1) $\exp$ は単調増加なので分位点をそのまま移すため、 $\log y$ の中央値（正規なら平均と同じ）を戻すと $y$ の中央値になる。 (2) 一方 $E[y]$ とは一致せず必ず下回る——ここで $\exp$ が凸であることから イェンセンの不等式が効きます。平均が欲しいなら $\exp(\hat\mu + \hat\sigma^2/2)$ のような補正が必要です。

③ ゼロや負の値には使えない。 $\log(y+1)$ で逃げる手はありますが、 係数の解釈が「%」から崩れるので安易に使うと意味を失います。

なお Box-Cox 変換は

という族の中から $\lambda$ を最尤法で選ぶ方法です。$\lambda=1$ なら変換不要、0.5 なら平方根、0 なら対数、$-1$ なら逆数。「どの変換がよいか」をデータに決めさせる枠組みで、$\lambda$ の信頼区間が1を含むなら変換しないという判断ができます。

クック距離の閾値は 4/n か 1 か

諸説あって混乱していたので、決着をつけました。完全に仮定を満たす正常なデータを作り、誤検出率を測ります（$p$=3、各条件3000回）。

$n$=200 の正常なデータで誤検出率100%。 つまり必ず誰かが引っかかります。

理由を最初は「$n$ が増えると $4/n$ は0に近づくのに $D_i$ の分布は縮まないから」だと考えたのですが、これは間違いでした。確かめると $D_i$ の分布も同じ速さで縮んでいます。

$E[D_i]\approx1/(n-p)$ なので、$4/n$ は $n$ によらず常に「平均の約4倍」＝上位5〜9%あたりの位置に張り付いています。つまり $4/n$ は絶対的な基準ではなく、最初から「上位数%を拾う線」なのです。

なお厳密には $E[D_i \mid X] = \frac{1}{p}\overline{\left(\frac{h_{ii}}{1-h_{ii}}\right)}$ で、$1/(n-p)$ に一致するのは全部のてこ比が $p/n$ に揃っているときだけです。$h/(1-h)$ が凸関数なのでイェンセンの不等式により常に $1/(n-p)$ 以上になり、$n$ が小さいほど上振れします（$n$=15 で 1.14倍、$n$=200 で 1.007倍）。

そして本当の原因は多重比較でした。1点あたり5〜9%の判定を $n$ 個くり返すので、少なくとも1個が引っかかる確率はほぼ $1-(1-0.05)^n$ で、$n$=200 なら実質1になります。実測でも平均10.71個が引っかかっていました。

第16回で「無意味な変数100本のうち5.55本が偽陽性になる」という話をしましたが、まったく同じ構造です。検定を $n$ 回くり返せば、5%水準なら $0.05n$ 件は必ず出ます。

逆側（本物を見逃さないか）も測りました。$n$=30 で1点を汚染します。この実験では汚染点を $x$ 方向にも遠い位置に置いています（前半のケースCと同じ状況）。

置く位置を重心に変えると結果が一変します。

重心に置いた場合、$y$ を6ずらしても $D_i>1$ の検出率は0%です。これは欠陥ではなく仕様で、クック距離は影響力の指標なので「影響していない外れ値」は当然拾いません（前半のケースAそのものです）。$4/n$ なら 99.6% 拾えるので、目印としての $4/n$ が役に立つ場面でもあります。いずれにせよ $t_i$ との併用が必要になります。

結論：$4/n$ は「閾値」ではなく「並べ替えて上位を見るための目印」です。

• $4/n$ … 感度は高い（小さい汚染も89%拾う）が偽陽性だらけ。 「$D_i$ を降順に並べて上から数点を目で見る」出発点として使う。 「$4/n$ を超えたら異常」と読んではいけない
• $D_i>1$ … 超えたらほぼ確実に本物（$n\ge30$ で誤検出1%未満）。ただし小さい影響は見逃す
• 実務的な答え … 閾値で自動判定しようとしないこと。全点プロットして「他から突き抜けている点があるか」を相対的に見るのが本来の使い方

前半の重回帰の例で、1点だけが他の86倍だったように、 本物の高影響点は相対的に明らかに浮きます。

表の3列目にある $F_{0.5}(p, n-p)$ について補足します。$D_i$ は形の上では「$\hat\beta$ の信頼楕円体の上をどれだけ移動したか」を測っていて、F統計量と同じスケールになります。そこから Cook 自身は「$D_i > F_{0.5}$ なら、1点抜くと $\hat\beta$ が50%信頼領域の外に出る」という読み方を提案しました。実測すると $F_{0.5}(3, n-3)$ は 0.835 から 0.791 と1に近い値で、誤検出率も $D_i>1$ とほぼ同じでした。よく使われる「$D_i>1$」という目安は、この $F_{0.5}\approx0.8$ を丸めたものだと理解できます。

$n$ が大きいときは $D_i$ より DFBETAS（係数ごとの変化を標準化した量、目安 $2/\sqrt{n}$）や DFFITS（目安 $2\sqrt{p/n}$）のほうが「どの係数が動いたか」まで分かって実用的です。

誤差の自己相関とダービン・ワトソン比

DW比の定義はこれです。

(a) 実測で近似式を確認しました。

範囲は 0〜4 で、2付近なら無相関、0に近ければ正の自己相関、4に近ければ負の自己相関。

$\rho$ が大きいところで理論とズレる（$\rho$=0.95 で実測 0.354 対 理論 0.100）理由は2つあります。1つは回帰が誤差の滑らかな成分を直線に吸収してしまうので、残差の自己相関が真の誤差より弱まること（$\rho$=0.95 のとき真の誤差では 0.855 なのに残差では 0.791 でした）。もう1つは近似式が落としている端点項で、厳密には

です（$\hat\rho = \sum_{i=2}^n e_ie_{i-1}/\sum_i e_i^2$ と定義したとき厳密に成立します。実測で差は $10^{-16}$ 台）。

この2つは逆向きに働きます。吸収効果は DW を上へ押し、端点項は必ず正なので下へ引く。$\rho$=0.95 の系列で内訳を出すと、残差の自己相関 $\hat\rho = 0.791$ から $2(1-\hat\rho) = 0.418$、そこから端点項 0.053 を引いて 0.365。上の表の 0.354 は多数回の平均なので完全一致はしませんが、理論値 0.100 との差の大半が吸収効果によるものだと分かります。正味では吸収効果が勝つので、正の自己相関のときは実測が理論値より大きく出ます。

(b)(c) 自己相関の害は不均一分散よりはるかに深刻でした。

$\rho$=0.9 では被覆率が 34% まで崩壊します。不均一分散が 93.7% だったのと比べると桁違いです。

理由は直感的に言えて、自己相関があると「実質的な標本サイズ」が減るからです。隣の観測がほぼ同じ情報しか持っていないなら、$n$ 個あっても独立な情報は $n$ 個ない。それなのに $\sqrt{n}$ で割ってしまうので過小評価になります。

ただし、ここでも「必ず過小評価」ではありませんでした。 不均一分散のときと同じ罠です。向きは$x$ の並び方で決まります（$\rho$=0.7 で実測）。

「実質的な標本サイズが減る」が効くのは、$x$ 自身も滑らかに動いていて、誤差の滑らかさと噛み合うときです。$x$ が1つおきに反転する設計だと、正の自己相関は隣同士の差を安定させるので傾きの分散が下がり、SEは過大評価になります。

不均一分散では「てこ比との相関」で向きが決まりましたが、自己相関では「$x$ の変動と誤差の変動の噛み合わせ」で決まる。同じ構図でした。

(d) 予告した罠の実演です。 誤差を完全に独立に生成した曲線データで DW = 0.368。$x^2$ を入れて曲線を吸収したら DW = 2.092 に戻りました。DW が小さい原因は自己相関でなく線形性の破れだったわけです。

DW を使うときの注意点。

• 観測に意味のある順序がなければ計算しても無意味。 DW は「隣」を使うので、 並び順を変えると値が変わります。時系列・空間・ロット順など順序が実在するときだけ意味がある
• $x$ でソートしたデータでは曲線が自己相関に化ける（上の (d)）
• 検定表は $d_L$ と $d_U$ の2つの臨界値を持ち、間に入ると「判定不能」になります。 これは分布が $X$ に依存して正確な臨界値が計算できないという事情から来ています
• 1次の自己相関しか見ません。季節性などはブロイシュ・ゴッドフレイ検定や Ljung-Box 検定を使います
• 見つかったときの対処は HAC（Newey-West）標準誤差が手軽です。不均一分散のロバストSEと同じ発想で、$\Omega$ の非対角成分まで拾います。構造をモデル化したいなら時系列モデルや一般化最小二乗法へ

外れ値は削除していいのか

いちばん実務的な問いです。数値で答えを出しました。$t(3)$ という重い裾の誤差（外れ値が自然に発生する）で $n$=30、4000回の実験です。真のモデルは正常で、外れ値は分布の性質として出ているだけという設定です。

削除すると区間幅は 0.721倍に縮むので「精度が上がった」ように見えます。 ところが実際に当たる率は 95.73% から 91.07% へ下がり、繰り返すと 80.55% まで落ちます。

理由は同じデータで「外れ値を選ぶ」と「推定する」の両方をやっていることです。残差が大きい点を捨てる操作は $\hat\sigma$ を人為的に小さくするので、自分で自分のものさしを縮めていることになります。そしてこの選択の不確実性は、最後に報告するSEにまったく反映されません。

第16回で「有意でない変数を消す」という手順が偽陽性を生む話をしましたが、構造が同じです。同じデータで選抜と推定を繰り返すと、報告される不確実性が嘘になります。

ただし点推定の精度は改善する

ここは意外でした。同じ実験で $\hat\beta_1$ の精度を見ると、削除は有効なのです。

つまり「削除すると係数の当たりは良くなるが、その係数に付ける誤差の見積もりが嘘になる」。この2つを分けて理解する必要がありました。

そしてHuber ロバスト回帰は削除よりRMSEが小さく、しかも標本を1つも捨てていません。やることは重み付けで、まず残差を頑健なスケール $s$（$\mathrm{MAD}/0.6745$、MADは中央絶対偏差）で割って $u_i = e_i/s$ を作り、$|u_i|\le k$ なら重み1、超えたら $k/|u_i|$ に下げます（$k=1.345$ は誤差が正規のとき効率95%になる定数）。

これは冒頭の (b) の図に戻る話です。Huber は残差が小さいうちは二乗、大きくなったら絶対値に切り替える損失を使うので、$y$ 方向に大きくズレた点の効きが一定で打ち止めになります。

捨てるより重みを下げるほうが強い。 削除は「重み1か0」という極端なロバスト化で、境界のすぐ内外で扱いが不連続に変わります。 Huber はなめらかに下げるので情報を捨てすぎません。

ただし Huber は $y$ 方向しか守らない

さきほどの勾配の話がここで効いてきます。絶対値損失でも $(x_i-\bar x)$ は残るので、てこ比が高い点に対しては無力です。前半のケースCと同じデータ（20点＋汚染1点、真の傾き 1.0、$y$ は15下）で、汚染点の $x$ だけを動かしました。

$x$=16 なら Huber は 0.906 でかなり守り切りますが（OLS は 0.350 まで崩壊）、$x$=60 まで遠ざけると OLS の 0.740 と Huber の 0.735 がほぼ同じになります。$x$ 方向の破断点は0、つまり1点でも十分遠ければロバスト回帰も壊れるわけです。

記事の中心命題「壊れるには $x$ 遠と $y$ ズレの両方が必要」に照らすと、Huber は後者だけを守る道具でした。$x$ 方向の外れにも対処するには、てこ比も見る手法（MM推定、最小刈り込み二乗法＝LTS）が必要です。