線形回帰モデルで制限付きモデルとの残差平方和の差がカイ二乗分布に従うことを証明！！！

この記事では線形回帰モデルで制限付きモデルとの残差平方和の差がカイ二乗分布に従うことを証明します。

線形回帰モデルの設定

\begin{align*} y = X \beta + \varepsilon, \quad \varepsilon \sim N(0, \sigma^2 I) \end{align*}
という線形回帰モデルを考えます。\(X\)は1列目が全て1の\(n \times {p+1}\)行列です。
また、\(\beta = (\beta_0, \beta_1, \ldots, \beta_p)^t \in \mathbb R^{p + 1}\)という添字の振り方をしているとします。

例えば帰無仮説
\begin{align*} h_0: \quad \beta_k = \beta_{k+1} = \cdots = 0 \end{align*}
を考えてみます。(後述のハット行列\(H\)と表記がややこしくなるので、小文字の\(h\)で帰無仮説を表記しています）。

\begin{align*} X_{f} &= X \end{align*}
と表記することにします。また、\(X\)の\(k\)列以降を全て0ベクトルとしたベクトルを\(X_r\)とします。つまり、
\begin{align*} \tilde X_r = X \begin{bmatrix}I_k & 0 \\0 & 0 \end{bmatrix} \end{align*}
です。ついでに、\(X\)を最初の\(k\)列とそれ以降を分けて表記します。つまり、
\begin{align*} X = \begin{bmatrix} X_r & X_e \end{bmatrix}\end{align*}
により\(X_r, X_e\)を定めます。

また、
\begin{align*} \beta_f = \beta \end{align*}
とし、
\begin{align*} \tilde \beta_r = \begin{bmatrix}I_k & 0 \\0 & 0 \end{bmatrix} \beta \end{align*}
とします。ついでに最初の\(k\)行とそれ以降を分けて表記します。つまり、
\begin{align*} \beta = \begin{bmatrix} \beta_r \\ \beta_e \end{bmatrix} \end{align*}

ここで、
\begin{align*} \tilde X_r = \begin{bmatrix} X_r & 0 \end{bmatrix}\end{align*}
と、
\begin{align*} \tilde \beta_r = \begin{bmatrix} \beta_r \\ 0 \end{bmatrix} \end{align*}
という関係になっています。

帰無仮説の下では、
\begin{align*} y = X_f \tilde \beta_r + \varepsilon, \quad \varepsilon \sim N(0, \sigma^2 I) \end{align*}
と仮定していることになります。
\begin{align*} H_f = X_f \left( X_f^t X_f\right)^{-1} X_f^t \end{align*}
とします。
\begin{align*} H_r = X_r \left( X_r^t X_r\right)^{-1} X_r^t \end{align*}
と定めます。一応\(H_f, H_r\)はともに\(n\)次正方行列です。

残差平方和の差分

一般に、線形回帰モデルの残差平方和
\begin{align*} \sum (y_i – \hat{y}_i)^2 \end{align*}
は行列で表記すると、
\begin{align*} y^t \left(I – H \right) y \end{align*}
でした。
\begin{align*} RSS_f =y^t \left(I – H_f \right) y \end{align*}
とし、
\begin{align*} RRS_r = y^t \left(I – H_r \right) y\end{align*}
とします。これらの差を考えます。
\begin{align*} RRS_f – RRS_r = y^t \left(I – H_f \right)y – y^t \left(I – H_r \right) y = y^t \left(H_f – H_r \right) y\end{align*}
です。

\begin{align*} X_f \beta_f &= X_f \begin{bmatrix}I_k & 0 \\0 & 0 \end{bmatrix} \beta \\&= X_f \begin{bmatrix}I_k & 0 \\0 & 0 \end{bmatrix} \begin{bmatrix}I_k & 0 \\0 & 0 \end{bmatrix} \beta \\&= \tilde X_r \tilde \beta_r \end{align*}
であり、
\begin{align*} \tilde X_r \tilde \beta_r = X_r \beta _r \end{align*}
ということを念頭に置いておきます（すみません極めて当たり前かもしれないで）。つまり、
\begin{align*} X_f \beta_f = X_r \beta _r \end{align*}

\begin{align*} y^t H_f y &= \beta_f X_f ^t H_f X_f \beta _f + \varepsilon^t H_f X_f \beta_f + \beta_f^t X_f^t H_f \varepsilon + \varepsilon^t H_f \varepsilon
\\&= \beta_f X_f ^t X_f \beta _f + \varepsilon^t X_f \beta_f + \beta_f^t X_f^t \varepsilon + \varepsilon^t H_f \varepsilon
\\&= \beta_r X_r ^t X_r \beta _r + \varepsilon^t X_r \beta_r + \beta_r^t X_r^t \varepsilon + \varepsilon^t H_f \varepsilon \end{align*}
です。ただし\(H_f X_f = X_f\)や\(X_f \beta_f = X_r \beta_r\)をところどころ用いています。
また同様に、
\begin{align*} y^t H_r y &= \beta_r X_r ^t H_r X_r \beta _r + \varepsilon^t H_r X_r \beta_r + \beta_r^t X_r^t H_r \varepsilon + \varepsilon^t H_r \varepsilon
\\&= \beta_r X_r ^t H_r X_r \beta _r + \varepsilon^t H_r X_r \beta_r + \beta_r^t X_r^t H_r \varepsilon + \varepsilon^t H_r \varepsilon
\\&= \beta_r X_r ^t X_r \beta _r + \varepsilon^t X_r \beta_r + \beta_r^t X_r^t \varepsilon + \varepsilon^t H_r \varepsilon \end{align*}
となります。というわけで、
\begin{align*}y^t \left(H_f – H_r \right) y
&= \varepsilon^t \left(H_f – H_r \right) \varepsilon \end{align*}
となります。

\begin{align*} \frac{1}{\sigma} \varepsilon \end{align*}
を考えることにします。つまり、
\begin{align*} \frac{1}{\sigma} \varepsilon \sim N(0,I)\end{align*}
ですので、多変量標準正規分布です。
\(\left(H_f – H_r \right)\)は階数\(p+1 – k\)の実対称冪等行列であることがわかります。
適当に直交行列\(S \in O_n \)で対角化することで、
\begin{align*} \left(H_f – H_r \right) = S^t \Lambda_{p+1 – k }S\end{align*}
とすることができます。ただし、
\begin{align*} \Lambda_p = \begin{bmatrix}I_{p+1-k} & 0 \\0 & 0 \end{bmatrix} \end{align*}
という表記を用いています。
\begin{align*} \frac{1}{\sigma} \varepsilon^t \left(H_f – H_r \right) \frac{1}{\sigma} \varepsilon
&= \frac{1}{\sigma} \varepsilon^t S^t \Lambda_{p + 1 – k} S \frac{1}{\sigma} \varepsilon \end{align*}
とできます。
\begin{align*} S \frac{1}{\sigma} \varepsilon \end{align*}
は多変量標準正規分布の直交行列による変換なので多変量標準正規分布に従います。
というわけで、
\begin{align*}\frac{1}{\sigma} \varepsilon^t S^t \Lambda_{p + 1 – k} S \frac{1}{\sigma} \varepsilon \end{align*}
は標準正規分布の2乗の\(p+1-k\)個の和なので自由度\(p + 1 -k\)のカイ二乗分布に従います。
というわけで、帰無仮説の下で
\begin{align*} \frac{1}{\sigma^2 } \left(y^t \left(I – H_f \right) y \,\, – \,\, y^t \left(I – H_r \right) y \right) \sim \chi_{p + 1 – k} \end{align*}
ということが分かります。

回帰平方和との関係

回帰平方和は
\begin{align*} ESS = \sum \left(\hat y_i – \bar y_i \right) ^2 \end{align*}
により定義されますが、これを行列の表記に一旦直してみます。
\begin{align*} e = (1, 1, \ldots, 1)^t \in \mathbb R^{n}\end{align*}
と表記することにすると、
\begin{align*} ESS = \sum \left(\hat y_i – \bar y_i \right) ^2 = y^t \left( H – \frac{1}{n} e e^t \right) y\end{align*}
であることが分かります。\(ESS_f = y^t \left( H_f – \frac{1}{n} e e^t \right) y\)と表記することにします。
というのも、
\begin{align*} \hat y = \frac{1}{n} e e^t y \end{align*}
だからです。ここで、帰無仮説として切片\(\beta_0\)以外の回帰係数が全て\(0\)であるという帰無仮説を考えてみます。
つまり、
\begin{align*} h_0 : \beta_1 = \beta_2 = \cdots = \beta_p = 0\end{align*}
ということです。
\begin{align*} X_r = \begin{bmatrix} 1 \\ 1 \\ \vdots \\ 1 \end{bmatrix}\end{align*}
と定めます。
\begin{align*} H_r &= X_r\left( X_r^t X_r \right)^{-1}X_r^t \\& = X_r \frac{1}{n} X_r^t \\&= \frac{1}{n}ee^t \end{align*}
となります。
つまり、帰無仮説において
\begin{align*} RSS_f – RSS_r = ESS_f \end{align*}
であることがわかります。