統計学 学習ノート

はじめに

「$\hat{\beta} = (X^{\top}X)^{-1}X^{\top}y$ は覚えた。でも、なぜそれで $F$ 分布が出てくるの？」

準1級までの回帰分析は、公式を使って係数を求め、出力を読むことが中心でした。ところが1級では、「残差平方和を $\sigma^2$ で割ると $\chi^2$ 分布に従うことを示せ」のように、公式の裏側を問われます。

例えば、地図アプリの経路案内を思い浮かべてください。案内どおりに歩くことはできても、なぜその道が最短なのかを説明するには、地図全体を見る必要があります。回帰分析の場合、その地図にあたるのが、$n$ 次元の空間の図形です。

実は、回帰分析の定理はほとんどすべて、ある1つの事実から出てきます。それは、「最小二乗法はデータを平面へ直角に落とす操作である」という事実です。三平方の定理が平方和の分解になり、直角であることが独立性になります。

そこで登場するのが、今回の主役 最小二乗推定 と ハット行列 です。

前回の記事では、平方和を要因ごとに分ける分散分析を学びました。回帰分析も同じ線形モデルの仲間で、道具はそのまま使えます。今回の記事では、最小二乗推定の導出、残差の性質と回帰の分散分析、決定係数と重相関係数、変数変換、そして平均への回帰について解説します。

最小二乗推定を行列で導く

はじめに、最小二乗推定の導出について解説します。

説明変数が1つのモデルを 線形単回帰、$p$ 個のモデルを 線形重回帰 と呼びます。

$$y_i = \beta_0 + \beta_1 x_i + \varepsilon_i, \qquad y_i = \beta_0 + \beta_1 x_{i1} + \cdots + \beta_p x_{ip} + \varepsilon_i$$

どちらも行列で $y = X\beta + \varepsilon$ と書けます。$X$ は1列目がすべて1の $n \times (p+1)$ 行列で、列は一次独立とします。誤差には $E(\varepsilon) = 0$、$V(\varepsilon) = \sigma^2 I$ を仮定します。

最小二乗推定 は、$Q(b) = \|y - Xb\|^2$ を最小にする $b$ を $\beta$ の推定値とする方法です。$b$ で微分して0とおくと、次の 正規方程式 が得られます。

$$X^{\top}X\hat{\beta} = X^{\top}y, \qquad \hat{\beta} = (X^{\top}X)^{-1}X^{\top}y$$

ただし、微分が0というだけでは、最小であることの証明になりません。答案では、次の分解を示します。任意の $b$ について、$y - Xb = (y - X\hat{\beta}) + X(\hat{\beta} - b)$ です。

$$\|y - Xb\|^2 = \|y - X\hat{\beta}\|^2 + \|X(\hat{\beta} - b)\|^2 + 2(\hat{\beta} - b)^{\top}X^{\top}(y - X\hat{\beta})$$

正規方程式から $X^{\top}(y - X\hat{\beta}) = 0$ なので、最後の項は消えます。第2項は0以上で、$X$ の列が一次独立なら $b = \hat{\beta}$ のときだけ0です。したがって、$\hat{\beta}$ が唯一の最小解です。

ここで、図形の見方をします。$X$ の列が張る平面を、部屋の床だと思ってください。データ $y$ は空中の1点です。真上から光を当てると、床に影が落ちます。この影が当てはめ値 $\hat{y} = X\hat{\beta}$ で、影と本体を結ぶ線分が残差です。床の中で本体に最も近いのは、真下の点ですよね。

影を作る操作を行列で書いたものが ハット行列 です。$y$ に帽子をかぶせるので、この名前があります。

$$H = X(X^{\top}X)^{-1}X^{\top}, \qquad \hat{y} = Hy$$

$H^{\top} = H$ と $H^2 = H$ は、代入すれば確かめられます。一度影にしたものを、もう一度影にしても変わらない、というのが $H^2 = H$ の意味です。

具体的に見てみましょう。5人の1日の勉強時間 $x$ と模試の点の伸び $y$ が、次のようでした。

$x$	1	2	3	4	5
$y$	2	3	5	4	6

【正規方程式を解く手順】

1. $X^{\top}X$ を作る： 成分は $n = 5$、$\sum x_i = 15$、$\sum x_i^2 = 55$
2. $X^{\top}y$ を作る： 成分は $\sum y_i = 20$、$\sum x_i y_i = 69$
3. 逆行列を掛ける

$$\hat{\beta} = \begin{pmatrix} 5 & 15 \\ 15 & 55 \end{pmatrix}^{-1} \begin{pmatrix} 20 \\ 69 \end{pmatrix} = \begin{pmatrix} 1.1 & -0.3 \\ -0.3 & 0.1 \end{pmatrix} \begin{pmatrix} 20 \\ 69 \end{pmatrix} = \begin{pmatrix} 1.3 \\ 0.9 \end{pmatrix}$$

回帰直線は $\hat{y} = 1.3 + 0.9x$ です。単回帰の公式 $\hat{\beta}_1 = S_{xy}/S_{xx} = 9/10$ とも一致します。

最後に、推定量の性質です。$y = X\beta + \varepsilon$ を代入すると、$\hat{\beta} = \beta + (X^{\top}X)^{-1}X^{\top}\varepsilon$ と書けます。ここから、次の2つがすぐに出ます。

$$E(\hat{\beta}) = \beta, \qquad V(\hat{\beta}) = (X^{\top}X)^{-1}X^{\top}(\sigma^2 I)X(X^{\top}X)^{-1} = \sigma^2 (X^{\top}X)^{-1}$$

単回帰では、$(X^{\top}X)^{-1}$ の右下の成分が $1/S_{xx}$ なので、$V(\hat{\beta}_1) = \sigma^2/S_{xx}$ です。上の例でも、右下の成分は $0.1 = 1/10$ になっています。

試験のポイント：「最小である」ことの証明では、微分だけで終わらせず、平方和の分解で交差項が消えることを示します。

以上が、最小二乗推定の導出についての解説です。

残差と回帰の分散分析

次に、残差の性質と回帰の分散分析について説明します。

残差 は $e = y - \hat{y} = (I - H)y$ です。実測と当てはめ値の差で、誤差 $\varepsilon$ そのものではありません。$(I - H)X = 0$ なので $e = (I - H)\varepsilon$ と書けて、誤差のうち床に垂直な成分だけが残差として見えています。

【残差の性質】

1. $X^{\top}e = 0$： 残差は $X$ のすべての列と直交する
2. $\sum_i e_i = 0$： $X$ の1列目がすべて1だから（切片があるときだけ）
3. $\hat{y}^{\top}e = 0$： 当てはめ値とも直交する
4. $V(e) = \sigma^2 (I - H)$： 残差の分散は点ごとに違い、等分散ではない

次に、$\sigma^2$ の不偏推定量を導きます。$e^{\top}e = \varepsilon^{\top}(I - H)\varepsilon$ です。一般に $E(\varepsilon^{\top}A\varepsilon) = \sigma^2 \operatorname{tr}A$ で、トレースは積の順番を入れ替えられるので次が成り立ちます。

$$\operatorname{tr}H = \operatorname{tr}\{(X^{\top}X)^{-1}X^{\top}X\} = p + 1, \qquad E(e^{\top}e) = (n - p - 1)\sigma^2$$

したがって、$s^2 = e^{\top}e/(n - p - 1)$ が不偏です。推定に使った係数の個数だけ、分母が減るのです。$n$ 次元のうち $p + 1$ 次元は床に吸い取られ、垂直な方向に残る自由度が $n - p - 1$ だから、と言っても同じです。

続いて、回帰の分散分析 です。$\bar{y}$ を並べたベクトルを $\bar{y}\mathbf{1}$ とします。切片があるので $\mathbf{1}$ は床の上にあり、$\hat{y} - \bar{y}\mathbf{1}$ も床の上です。残差は床と直交するので、三平方の定理が使えます。

$$\underbrace{\sum_i (y_i - \bar{y})^2}_{S_T} = \underbrace{\sum_i (\hat{y}_i - \bar{y})^2}_{S_R} + \underbrace{\sum_i e_i^2}_{S_E}$$

全体のばらつき $S_T$ が、回帰で説明できる部分 $S_R$ と、説明できない部分 $S_E$ に分かれました。これを表にまとめます。

要因	平方和	自由度	平均平方	$F$
回帰	$S_R$	$p$	$S_R/p$	$\dfrac{S_R/p}{S_E/(n-p-1)}$
残差	$S_E$	$n-p-1$	$S_E/(n-p-1)$
全体	$S_T$	$n-1$

帰無仮説 $H_0: \beta_1 = \cdots = \beta_p = 0$ のもとで、$\varepsilon$ が正規分布なら、$F$ は自由度 $(p, n-p-1)$ の $F$ 分布に従います。この証明の筋道を押さえましょう。

【F分布を導く筋道】

1. $H_1 = \mathbf{1}\mathbf{1}^{\top}/n$ とおくと、$H_0$ のもとで $S_R = \varepsilon^{\top}(H - H_1)\varepsilon$、$S_E = \varepsilon^{\top}(I - H)\varepsilon$ である
2. $\mathbf{1}$ は床の上にあるので $HH_1 = H_1$ である。ここから $H - H_1$ と $I - H$ はどちらも対称で冪等で、$(H - H_1)(I - H) = 0$ となる
3. 対称で冪等な行列 $A$ は、直交行列で対角化すると、対角成分が1と0だけになる。1の個数は $\operatorname{tr}A$ である
4. 標準正規ベクトルは直交変換しても標準正規ベクトルなので、$\varepsilon^{\top}A\varepsilon/\sigma^2$ は $\operatorname{tr}A$ 個の標準正規変数の2乗の和になる
5. 積が0の2つの行列は、別々の座標を拾う。よって $S_R/\sigma^2 \sim \chi^2_p$ と $S_E/\sigma^2 \sim \chi^2_{n-p-1}$ は独立である

あとは $F$ 分布の定義どおりです。前回の分散分析の証明と、まったく同じ骨組みですね。

同じ理由で、$\hat{\beta}$ と $e$ も独立です。$\operatorname{Cov}(Hy, (I-H)y) = \sigma^2 H(I-H) = 0$ で、正規分布では無相関が独立を意味するからです。$c_{jj}$ を $(X^{\top}X)^{-1}$ の対角成分とすると、個々の係数の検定統計量が $t$ 分布に従います。

$$t_j = \frac{\hat{\beta}_j - \beta_j}{s\sqrt{c_{jj}}} \sim t_{n-p-1}$$

例題：単回帰の分散分析

先ほどの5人のデータについて、次に答えよ。〔1〕残差を求め、$\sum e_i = 0$ と $\sum x_i e_i = 0$ を確かめよ。〔2〕分散分析表を作り、$H_0: \beta_1 = 0$ を有意水準5%で検定せよ。ただし $F_{0.05}(1, 3) = 10.13$ とする。〔3〕単回帰では $t_1^2 = F$ であることを示せ。

【解答】

〔1〕当てはめ値は 2.2, 3.1, 4.0, 4.9, 5.8 だから、残差は $-0.2, -0.1, 1.0, -0.9, 0.2$ である。和は0で、$\sum x_i e_i = -0.2 - 0.2 + 3.0 - 3.6 + 1.0 = 0$ である。

〔2〕$S_T = 10$、$S_E = 0.04 + 0.01 + 1 + 0.81 + 0.04 = 1.9$、$S_R = 10 - 1.9 = 8.1$ である。

要因	平方和	自由度	平均平方	$F$
回帰	8.1	1	8.1	12.79
残差	1.9	3	0.633
全体	10.0	4

$F = 12.79 > 10.13$ なので、$H_0$ を棄却する。勉強時間の効果はあると判断する。

〔3〕$\hat{y}_i - \bar{y} = \hat{\beta}_1(x_i - \bar{x})$ より $S_R = \hat{\beta}_1^2 S_{xx}$ である。$\beta_1 = 0$ のとき $t_1 = \hat{\beta}_1/(s/\sqrt{S_{xx}})$ だから、次を得る。

$$t_1^2 = \frac{\hat{\beta}_1^2 S_{xx}}{s^2} = \frac{S_R}{S_E/(n-2)} = F$$

試験のポイント：$S_R$ と $S_E$ の独立性は「対称で冪等な行列の積が0」から示します。自由度は行列のトレースとして求めます。

以上が、残差と回帰の分散分析についての説明です。

決定係数と重相関係数

続いて、当てはまりの指標について解説します。

決定係数 は、全体のばらつきのうち回帰で説明できた割合です。

$$R^2 = \frac{S_R}{S_T} = 1 - \frac{S_E}{S_T}$$

先ほどの例では $8.1 \div 10 = 0.81$ です。

重相関係数 は、実測値 $y$ と当てはめ値 $\hat{y}$ の相関係数です。これは決定係数の正の平方根に一致します。証明しましょう。切片があれば $\sum e_i = 0$ なので、$\hat{y}$ の平均も $\bar{y}$ です。$y_i - \bar{y} = (\hat{y}_i - \bar{y}) + e_i$ と分けると、$\hat{y}$ と $e$ の直交から次が成り立ちます。

$$\sum_i (y_i - \bar{y})(\hat{y}_i - \bar{y}) = \sum_i (\hat{y}_i - \bar{y})^2 = S_R$$

したがって、相関係数は $S_R/\sqrt{S_T S_R} = \sqrt{S_R/S_T}$ となります。ここがポイントです。「決定係数は相関の2乗」という言い方は、単回帰だけの話ではありません。重回帰でも、$y$ と $\hat{y}$ の相関の2乗として成り立ちます。単回帰では $\hat{y}$ が $x$ の1次式なので、$R^2 = r_{xy}^2$ となるわけです。

綱引きの例えで言えば、$p$ 人の選手が力を合わせて最も強く引いたときの力が、重相関係数です。実際、$\hat{y}$ は、説明変数のあらゆる1次結合のうちで、$y$ との相関が最も高いものになっています。

補足（ここは読み飛ばしてOK）： 中心化した床の上の任意のベクトル $v$ について、中心化した $y$ と $v$ の内積は、中心化した $\hat{y}$ と $v$ の内積に等しくなります。コーシー・シュワルツの不等式から、相関が最大になるのは $v$ が $\hat{y}$ の方向のときです。

分散分析の $F$ は、分子分母を $S_T$ で割ると、$R^2$ だけで書けます。

$$F = \frac{R^2/p}{(1 - R^2)/(n - p - 1)}$$

一方で、説明変数を増やすと、$R^2$ は決して下がりません。床が広がると、$y$ に近い点を探せる範囲が広がるので、$S_E$ は増えないからです。意味のない変数を入れても $R^2$ は上がるので、変数の数を比べるときは自由度調整済み決定係数やAICを使います。

例題：相関係数から決定係数を求める

$n = 30$ の標本で、相関係数が $r_{y1} = 0.8$、$r_{y2} = 0.6$、$r_{12} = 0.5$ だった。$y$ を $x_1, x_2$ に回帰する。〔1〕標準化した変数の回帰係数 $b_1, b_2$ を求めよ。〔2〕決定係数と重相関係数を求めよ。〔3〕$H_0: \beta_1 = \beta_2 = 0$ を検定せよ。ただし $F_{0.05}(2, 27) = 3.35$ とする。

【解答】

〔1〕すべての変数を平均0、分散1に標準化すると、正規方程式は相関行列で書ける。

$$\begin{pmatrix} 1 & 0.5 \\ 0.5 & 1 \end{pmatrix} \begin{pmatrix} b_1 \\ b_2 \end{pmatrix} = \begin{pmatrix} 0.8 \\ 0.6 \end{pmatrix}$$

解くと、$b_1 = (0.8 - 0.3)/0.75 = 2/3$、$b_2 = (0.6 - 0.4)/0.75 = 4/15$ である。

〔2〕標準化した変数では $S_T = y^{\top}y$、$S_R = \hat{y}^{\top}y = b_1 x_1^{\top}y + b_2 x_2^{\top}y$ で、$x_j^{\top}y / y^{\top}y = r_{yj}$ である。よって $R^2 = b_1 r_{y1} + b_2 r_{y2} = 0.533 + 0.160 = 0.693$ である。重相関係数は $R = \sqrt{0.693} = 0.833$ である。

〔3〕$F = (0.693/2)/(0.307/27) = 30.5 > 3.35$ なので、$H_0$ を棄却する。

$r_{y1}^2 + r_{y2}^2 = 1.00$ にはなりません。$x_1$ と $x_2$ に相関があるので、説明の重なりの分だけ差し引かれるのです。

以上が、決定係数と重相関係数についての解説です。

変数変換で線形に持ち込む

さらに、曲線の関係を扱う変数変換について説明します。

「線形」モデルの線形とは、$x$ についてではなく、母数 $\beta$ について1次 という意味です。ですから、$x$ や $y$ を前もって変換しても、最小二乗推定の理論はそのまま使えます。曲がった針金をアイロンでまっすぐに伸ばしてから、定規を当てるようなものです。

【変数変換で線形にできる例】

• べき乗型 $y = ax^{b}u$： 両辺の対数をとると $\log y = \log a + b\log x + \log u$
• 指数型 $y = ae^{bx}u$： 対数をとると $\log y = \log a + bx + \log u$
• 多項式 $y = \beta_0 + \beta_1 x + \beta_2 x^2 + \varepsilon$： $x$ と $x^2$ を別の説明変数として重回帰

具体的に見てみましょう。店の広さ $x$ と売上 $y$ が、(1, 3)、(2, 5.2)、(4, 8.8)、(8, 15.5) だったとします。両対数をとって回帰すると、傾きは約0.787、切片は約1.098で、$\hat{y} = 3.00\, x^{0.787}$ となります。傾きは弾力性と呼ばれ、「$x$ が1%増えると $y$ が約0.79%増える」と読めます。

ただし、変換には注意点が3つあります。

【変数変換の注意点】

1. 誤差の入り方が変わる： $y = ax^{b}u$ は掛け算の誤差で、$y = ax^{b} + \varepsilon$ は対数で線形にならない
2. 逆変換すると偏る： $\log y \sim N(m, s^2)$ なら $e^{m}$ は $y$ の中央値で、平均は $e^{m + s^2/2}$
3. 決定係数を比べられない： $y$ と $\log y$ では $S_T$ が別物である

2番目を数値で見ると、$s = 0.2$ なら $e^{0.02} \approx 1.02$ で、$e^{m}$ は平均より約2%小さくなります。

試験のポイント：回帰の仮定（等分散、正規性）が求められるのは、変換後の誤差です。元の尺度での予測には逆変換の偏りが入ることも、答案に一言添えます。

以上が、変数変換についての説明です。

平均への回帰：回帰の語源

最後に、平均への回帰を見ていきましょう。

平均への回帰（回帰効果、regression to the mean）とは、1回目に平均から離れた値をとった個体が、2回目には平均の近くの値をとりやすい現象です。19世紀のゴルトンが、背の高い親の子は親ほどには高くならないことを見つけ、これが「回帰」という言葉の語源になりました。

理論を導きましょう。$(X, Y)$ が2変量正規分布に従い、平均 $\mu_X, \mu_Y$、標準偏差 $\sigma_X, \sigma_Y$、相関係数 $\rho$ をもつとします。$\beta = \rho\sigma_Y/\sigma_X$ とおくと、$\operatorname{Cov}(X, Y - \beta X) = 0$ です。正規分布では無相関が独立なので、$Y - \beta X$ は $X$ と独立です。$Y = \beta X + (Y - \beta X)$ の条件付き期待値をとると、次が得られます。

$$E(Y \mid X = x) = \mu_Y + \rho\frac{\sigma_Y}{\sigma_X}(x - \mu_X), \qquad V(Y \mid X = x) = \sigma_Y^2 (1 - \rho^2)$$

標準化すると、$E(Z_Y \mid Z_X = z) = \rho z$ です。$|\rho| < 1$ なので、標準偏差の単位で見て、2回目の予想はいつも1回目より平均に近づきます。

なぜでしょうか？ 成績を「実力 + 運」で考えると分かります。$X = T + e_1$、$Y = T + e_2$ として、実力 $T$ と運 $e_1, e_2$ は独立とします。$\operatorname{Cov}(X, Y) = V(T)$ なので、相関係数は次のようになります。

$$\rho = \frac{V(T)}{V(T) + V(e)}$$

つまり、$\rho$ は「ばらつきのうち実力が占める割合」です。1回目の高得点者には、実力に加えて運のよかった人が多く含まれます。運は2回目に持ち越されないので、平均に近づくのです。ガチャの大当たりの翌日に、同じ運を期待できないのと同じです。

例題：2回のテストの回帰効果

全員が2回のテストを受ける。どちらの点も平均60点、標準偏差10点で、2変量正規分布に従い、相関係数は0.7とする。〔1〕1回目が80点だった人の2回目の点の条件付き期待値と条件付き標準偏差を求めよ。〔2〕その人の2回目が80点を超える確率を求めよ。〔3〕2回目の点の分散が1回目より小さくならない理由を説明せよ。

【解答】

〔1〕$E(Y \mid X = 80) = 60 + 0.7 \times 20 = 74$ である。条件付き標準偏差は $10\sqrt{1 - 0.49} = 7.14$ である。

〔2〕$P(Y > 80 \mid X = 80) = P(Z > 6/7.14) = P(Z > 0.84) = 0.20$ である。

〔3〕全分散の公式より、次のようになる。

$$V(Y) = V\{E(Y \mid X)\} + E\{V(Y \mid X)\} = \rho^2\sigma^2 + (1 - \rho^2)\sigma^2 = \sigma^2$$

条件付き期待値のばらつきは $\rho^2\sigma^2$ に縮むが、そのまわりに $(1 - \rho^2)\sigma^2$ のばらつきが加わるので、全体の分散は変わらない。

答えで大切なのは〔3〕です。回帰効果は、集団全体が平均へ縮んでいく現象ではありません。2回目に飛び抜けた点をとる人は、別の顔ぶれで現れます。

せっかくなので、向きを逆にした例も見てみましょう。2回目が80点だった人の1回目の予想も、同じ式で74点です。時間を巻き戻しても平均に近づくことから、これは時間の流れの効果ではなく、相関が1未満であることの効果だと分かります。

ここから、実務で最も大事な注意が出てきます。血圧が特に高い人だけを集めて健康法を試すと、何もしなくても次回の血圧は下がりやすくなります。成績の悪かった生徒を叱ると次は伸びるので、叱るのが効いたと勘違いするのも同じです。効果を確かめるには、同じ条件の対照群と比べる必要があります。

試験のポイント：回帰効果は $|\rho| < 1$ の統計的な性質で、因果関係ではありません。前後の差だけで効果を判断する分析を批判する問題で問われます。

以上が、平均への回帰についての解説です。

細かい部分を忘れてしまっても問題ありません。ここでは、「$H$ が影を作り、直交が三平方と独立性を生む」という点だけ押さえておきましょう。

今回のまとめ

それでは、今回の内容を振り返ります。

今回は最小二乗推定の理論を、ハット行列を軸に紹介しました。公式が多く見えますが、最小二乗法は床に影を落とす操作で、残差はその垂線です。平方和の分解も、$F$ 分布も、$R^2$ の性質も、「直角である」というただ1つの事実の言い換えでした。

今回のポイントを振り返りましょう。

振り返り

• 最小二乗推定は正規方程式から求まり、平方和の分解で最小性を証明する。
• ハット行列は影を作る対称で冪等な行列で、そのトレースが自由度を決める。
• 三平方の定理が $S_T = S_R + S_E$ になり、行列の積が0であることが $F$ 分布の独立性になる。
• 重相関係数は $y$ と $\hat{y}$ の相関係数で、その2乗が決定係数である。
• 平均への回帰は $|\rho| < 1$ から生じ、因果関係と取り違えてはいけない。

キーワードのおさらい

• 線形単回帰: 説明変数が1つの線形モデル
• 線形重回帰: 説明変数が複数の線形モデル。$y = X\beta + \varepsilon$ と書く
• 最小二乗推定: 残差平方和を最小にする推定。$\hat{\beta} = (X^{\top}X)^{-1}X^{\top}y$
• 回帰の分散分析: $S_T$ を $S_R$ と $S_E$ に分け、$F$ 検定を行うこと
• 重相関係数: $y$ と $\hat{y}$ の相関係数
• 決定係数: 全体の平方和のうち回帰で説明できた割合。$R^2 = S_R/S_T$
• 残差: 実測値と当てはめ値の差。$e = (I - H)y$
• 変数変換: 対数などで関係を母数について線形にすること
• 平均への回帰（回帰効果）: 平均から離れた値の次の値が、平均に近づきやすい現象

さて、ここまでは量的な変数を扱ってきました。では、「合格・不合格」と「教材A・教材B」のように、人数を数えるだけのデータでは、どのような導出が待っているのでしょうか？ しかも度数が小さいと、カイ二乗分布の近似も怪しくなります。

次回は、分割表の解析について紹介します。

それでは、次の記事へ進んでいきましょう。