線形回帰モデルのF検定 #統計検定

はじめに

統計検定1級の統計数理・理工学において，回帰分析の関連テーマ(最小二乗推定・線形単回帰・線形重回帰)は2012年～2022年の間で複数回登場した頻出分野の1つである。
本記事では，線形回帰モデルにおける3つの主な課題の1つである検定問題のうち，線形回帰モデルのF検定について整理する。

線形回帰モデル

統計検定1級の統計数理・理工学において，2012年～2020年の10年間で，線形回帰モデルのF検定は2014年に出題されている。

線形回帰モデルは，

$\begin{align} \boldsymbol{y} = \textbf{X} \boldsymbol{\beta} + \boldsymbol{\varepsilon} \tag{1} \end{align}$

であらわされる。
ここで $\boldsymbol{y} \in \mathbb{R}^n$ は目的変数， $\textbf{X} \in \mathbb{R}^{n \times d}$ は計画行列(または説明変数行列)， $\boldsymbol{\beta} \in \mathbb{R}^d$ は回帰係数， $\boldsymbol{\varepsilon} \in \mathbb{R}^n$ は誤差項である。

また誤差項は正規分布にしたがい， $\boldsymbol{\varepsilon} \sim N(\boldsymbol{0}_n, \sigma^2 I_n)$ とする。

F検定の問題設定

F検定と変数選定

線形回帰モデルに対してF検定を行なうモチベーションは，回帰係数のうち，特定の $\beta_i$ が0になるかどうかを検定するというものである。これは，機械学習でもおなじみの変数選択に通じるものがある。

(1)式に示した線形回帰モデルについて，計画行列を列で区切って2つの計画行列に分割する。これに対応して，回帰係数も2つに分割する。

$\begin{align} \textbf{X} &= [\textbf{X}_{(1)}, \textbf{X}_{(2)} ] \\ \boldsymbol{\beta} ^T &= [ \boldsymbol{\beta}_{(1)}^T, \boldsymbol{\beta}_{(2)}^T ] \end{align}$

ただし，

$\begin{align} \textbf{X} &\in \mathbb{R}^{n \times d}, \textbf{X}_{(1)} \in \mathbb{R}^{n \times s}, \textbf{X}_{(2)} \in \mathbb{R}^{n \times (d-s)} \\ \boldsymbol{\beta} &\in \mathbb{R}^d, \boldsymbol{\beta}_{(1)} \in \mathbb{R}^s, \boldsymbol{\beta}_{(2)} \in \mathbb{R}^{(d-s)} \end{align}$

である(なお， $d \gt s$ である)。

このように分割すると，線形回帰モデルは以下のようになる。

$\begin{align} \boldsymbol{y} = \textbf{X} \boldsymbol{\beta} + \boldsymbol{\varepsilon} = \textbf{X}_{(1)} \boldsymbol{\beta}_{(1)} + \textbf{X}_{(2)} \boldsymbol{\beta}_{(2)} + \boldsymbol{\varepsilon} \tag{2} \end{align}$

検定問題

F検定の検定問題は以下のようになる。

\begin{align}
\begin{cases}
H_0 &: \boldsymbol{\beta}_{(2)} = \boldsymbol{0} \\
H_1 &: \boldsymbol{\beta}_{(2)} \neq \boldsymbol{0} \text{(ただし$\boldsymbol{\beta} \in \mathbb{R}^d$)}
\end{cases}
\end{align}

部分空間による表現

残差平方和に関連する標本分布(カイ2乗分布) #統計検定 - jiku logでは，線形回帰モデルの部分空間による解釈のもと，残差の性質を考察した。今回も同様に，部分空間による解釈を行なう。

線形回帰モデルでは，回帰係数は計画行列 $\textbf{X}$ が張る $d$ 次元部分空間にあると解釈できる。帰無仮説 $H_0$ のもと，回帰係数はさらに次元数の少ない部分空間( $s$ 次元部分空間)にあると解釈できる。
具体例として，サンプル数が $n=3$ ，回帰係数の次元数が $d=2$ ，帰無仮説 $H_0$ のもとにおける回帰係数の次元数が $s=1$ ，という状況が挙げられる。

図示すると下図の様になる。

また，この $d$ 次元部分空間を真上から眺めると，下図のようになる。

残差平方和の算出

残差平方和に関連する標本分布(カイ2乗分布) #統計検定 - jiku logで説明したように，残差平方和はカイ2乗分布にしたがう。そのため，計画行列 $\textbf{X}$ を使ったときの残差と，帰無仮説のもとの計画行列 $\textbf{X}_{(1)}$ を使った残差の比を取るとF分布にしたがうので，帰無仮説の検定を行なうことができる。
ただし具体的なF分布を求めるためには，残差平方和がしたがうカイ2乗分布の自由度を求める必要がある。

F統計量の分母

はじめに，残差平方和に関連する標本分布(カイ2乗分布) #統計検定 - jiku log を参照しつつ，計画行列 $\textbf{X} \in \mathbb{R}^{(n \times d)}$ を用いたときの残差平方和 $\boldsymbol{e}^T \boldsymbol{e}$ を求めよう。
この残差平方和は，最小二乗推定量 $\boldsymbol{\hat{\beta}} = (\textbf{X}^T \textbf{X} )^{-1} \textbf{X}^T \boldsymbol{y}$ を用いて，

\begin{align}
\boldsymbol{e}^T \boldsymbol{e} &= (\boldsymbol{y} - \textbf{X} \boldsymbol{\hat{\beta}}) \\
&= \boldsymbol{y} ^T (I_n - P_X) \boldsymbol{y} \\
&= \boldsymbol{\varepsilon}^T \left( \sum_{i=d+1}^{n} \boldsymbol{u}_i \boldsymbol{u}_i ^T \right) \boldsymbol{\varepsilon}
= \sum_{i=d+1}^{n} z_i^2
\tag{3}
\end{align}

となる。

ただし， $P_X = \textbf{X} (\textbf{X}^T \textbf{X} )^{-1} \textbf{X}^T = \sum_{i=d+1}^{n} \boldsymbol{u}_i \boldsymbol{u}_i ^T$ は $d$ 次元部分空間への射影子である。また $z_i (i=d+1,\cdots n)$ は，正準形にしたときに表れる，正規分布 $N(0, \sigma^2)$ にしたがう確率変数である。

よって，

$\begin{align} \frac{\boldsymbol{e}^T \boldsymbol{e} }{ \sigma^2} \sim \chi^2 (n-d) \end{align}$

となり，自由度 $n-d$ のカイ2乗分布が得られた。

F統計量の分子

次に，帰無仮説 $H_0$ のもと計画行列 $\textbf{X}_{(1)} \in \mathbb{R}^{(n \times d)}$ を用いたときの残差平方和 $\boldsymbol{\tilde{e}}^T \boldsymbol{\tilde{e}}$ を求めよう。これは，上記の残差平方和において， $d$ 次元の部分を $s$ 次元(ただし，先に説明したとおり $d \gt s$ )に書き換えればよいので，

\begin{align}
\boldsymbol{\tilde{e}}^T \boldsymbol{\tilde{e}}
&= \boldsymbol{y} ^T (I_n - P_{X_{(1)}}) \boldsymbol{y} \\
&= \boldsymbol{\varepsilon}^T \left( \sum_{i=s+1}^{n} \boldsymbol{u}_i \boldsymbol{u}_i ^T \right) \boldsymbol{\varepsilon}
= \sum_{i=s+1}^{n} z_i^2
\end{align}

となる。

F検定に用いるF統計量のうち，分母に来るのは式(3)の左辺を自由度で割ったものである。一方，分子に来るものは残差平方和の差を自由度で割ったものであり，以下の式で表される。

$\begin{align} \boldsymbol{\tilde{e}}^T \boldsymbol{\tilde{e}} - \boldsymbol{e}^T \boldsymbol{e} &= \sum_{i=s+1}^{n} z_i^2 - \sum_{i=d+1}^{n} z_i^2 = \sum_{i=s+1}^{d} z_i^2 \\ &= \boldsymbol{y} ^T P_{X} \boldsymbol{y} - \boldsymbol{y} ^T P_{X_{(1)}} \boldsymbol{y} \end{align}$

よって，

$\begin{align} \frac{ \boldsymbol{\tilde{e}}^T \boldsymbol{\tilde{e}} - \boldsymbol{e}^T \boldsymbol{e} }{ \sigma^2} \sim \chi^2 (d-s) \end{align}$

となり，自由度 $d-s$ のカイ2乗分布が得られた。

残差平方和の振り分け

残差平方和は，正準形で表しておくと見通しがよい。上に示した，F統計量の分子・分母における残差平方和を図示すると以下のようになる。

※ただし $z_i$ は，過去記事で示した $\boldsymbol{z} = \textbf{G}^T \boldsymbol{\varepsilon}$ という変換ではなく，竹村彰通著「現代数理統計学」に倣って $\boldsymbol{z} = \textbf{G}^T \boldsymbol{y}$ という変換であることに注意する必要がある。

線形回帰モデルのF検定で一番つまずくところは，「F統計量の分子と分母ってなんだっけ…？」ということだが，これについて竹村彰通著「現代数理統計学」には以下のような説明があった。

回帰分析における係数ベクトルのF検定はさまざまの同値な形に表すことができる

つまり，検定に用いるF分布の自由度さえ押さえておけば，分母・分子の選び方には自由度があるということである。自由度を把握しておくためには，正準形を念頭に置いておくとわかりやすい。

(補足)正準形の注意点

前節で，「上図の $z_i$ は $\boldsymbol{z} = \textbf{G}^T \boldsymbol{y}$ という変形によるものである」，という説明を行なった。 $\boldsymbol{y}$ と $\boldsymbol{\varepsilon}$ の違いは以下の通りである。

$\begin{align} \boldsymbol{y} &\sim N(\textbf{X} \boldsymbol{\beta}, \sigma^2 I_n) \\ \boldsymbol{\varepsilon} &\sim N(\boldsymbol{0}, \sigma^2 I_n) \end{align}$

「残差平方和がカイ2乗分布にしたがう」という説明をしているときは，「計画行列 $\textbf{X}$ が張る空間の直交補空間」を考えているため，誤差の平均は0になっている。
ただし $\boldsymbol{z} = \textbf{G}^T \boldsymbol{y}$ という変換をすると， $d$ 次元部分空間中にある $z_i(i=1,\cdots,d)$ の平均は0とは限らなくなる。