以下の線形モデルを考えます.
$$
y=X\beta+\varepsilon,\quad \varepsilon\sim N(0,\sigma^{2}I)
$$
記号の説明をします.$n$ をサンプルサイズ,$p$ を説明変数の数とします.$y\in\mathbb{R}^{n}$ は目的変数ベクトルです.$X\in\mathbb{R}^{n\times p}$ は説明変数からなる計画行列です.$\varepsilon$ は誤差項で,平均 $0$ で分散 $\sigma^{2}I$ の多変量正規分布に従います.$\beta\in\mathbb{R}^{p}$ は回帰係数で,我々は $\beta$ の推定に興味があります.この記事を通して,$X^{\top}X$ は正則であることを仮定します.このとき,$\beta$ の最尤推定量 $\hat\beta$ が以下の形で与えられることはよく知られています.
$$
\hat\beta=(X^{\top}X)^{-1}X^{\top}y
$$
いま,回帰係数 $\mathbf\beta$ に関する帰無仮説
$$
H_0 : H^\intercal\mathbf\beta = \mathbf\xi
$$
を考えます.ただし,$H\in\mathbb{R}^{p\times q}$ は $\operatorname{rank} H = q$ を満たすとします.また,$\xi\in\mathbb{R}^{q}$ です.例えば,$H:=(1,0,\dots,0)^{\top},\xi:=0$ とすると,$H_{0}:\beta_{1}=0$ という帰無仮説になるため,説明変数が目的変数に寄与しているかを検証できます.本記事の目的は,帰無仮説の下で
F_{0}:=\frac{n-p}{q}\frac{(H^{\top}\hat\beta - \xi )^\intercal \{H^{\top}(X^{\top}X)^{-1}H\}^{-1}( H^{\top}\hat\beta - \xi)}{\|y-X\hat\beta\|^2} \sim F(q,n-p)
を示すことです.ただし,$F(q,n-p)$ は自由度 $(q,n-p)$ の F 分布です.$F_{0}$ を検定統計量として,帰無仮説 $H_{0}$ に関する仮説検定が行えます.$\sigma^{2}$ に依らないこと,検定統計量が従う分布がよく知られた分布であることがポイントです.
残差平方和が従う分布
$e:=y-X\hat\beta$ とします.$\|e\|^{2}$ が従う分布を導出します.はじめに,$P=I-X(X^{\top}X)^{-1}X^{\top}$ とします.$P$ は $P^{2}=P$ を満たす対称行列であることが確認できます.$PX=O$ に注意すると,$e$ は次のように整理できます.
\begin{align*}
e&=y-X\hat\beta\\
&=y-X(X^{\top}X)^{-1}X^{\top}y\\
&=Py\\
&=PX\beta+P\varepsilon\\
&=P\varepsilon
\end{align*}
したがって,$\|e\|^{2}=\varepsilon^{\top}P^{\top}P\varepsilon=\varepsilon^{\top}P\varepsilon$ が成り立ちます.ここから,$\sigma^{-2}\|e\|^{2}$ が自由度 $n-p$ の $\chi^{2}$ 分布 $\chi^{2}(n-p)$ に従うことを頑張って示します.
冪等行列の性質
$A^{2}=A$ を満たす正方行列を冪等行列といいます.冪等行列 $A$ は以下の性質をもちます.
- $A$ の任意の固有値 $\lambda$ は $\lambda\in{0,1}$ を満たす
- $\operatorname{rank}A = \operatorname{tr}A$
- $I-A$ は冪等行列
多変量正規分布の二次形式の積率母関数
$A\in\mathbb{R}^{n\times n}$ は対称行列とします.$\varepsilon^{\top}A\varepsilon$ の積率母関数 $M_{\varepsilon^{\top}A\varepsilon}(t)$を求めます.
\begin{align*}
M_{\varepsilon^{\top}A\varepsilon}(t)
&=E\left[\exp\left(t\varepsilon^{\top}A\varepsilon\right)\right]\\
&=\frac{1}{(2\pi\sigma^{2})^{n/2}}\int_{\mathbb{R}^{n}}
\exp\left[t\varepsilon^{\top}A\varepsilon-\frac{1}{2\sigma^{2}}\varepsilon^{\top}\varepsilon\right]\,\mathrm{d}\varepsilon\\
&=\frac{1}{(2\pi\sigma^{2})^{n/2}}\int_{\mathbb{R}^{n}}
\exp\left[-\frac{1}{2}\varepsilon^{\top}\left(\frac{1}{\sigma^{2}}I-2tA\right)\varepsilon\right]\,\mathrm{d}\varepsilon\\
\end{align*}
ここで,絶対値が十分小さい $t$ については,$B_{t}:=\sigma^{-2}I-2tA$ は正定値行列になるので,
\begin{align*}
M_{\varepsilon^{\top}A\varepsilon}(t)
&=\frac{1}{(2\pi\sigma^{2})^{n/2}}\int_{\mathbb{R}^{n}}
\exp\left[-\frac{1}{2}\varepsilon^{\top}B_{t}\varepsilon\right]\,\mathrm{d}\varepsilon\\
&=\frac{1}{(2\pi\sigma^{2})^{n/2}}(2\pi)^{n/2}(\det{B_{t}})^{-1/2}
\end{align*}
となります.最後に,$\det{B_{t}}=\sigma^{-2n}\det(I-2t\sigma^{2}A)$ を用いると,
$$
M_{\varepsilon^{\top}A\varepsilon}(t)
=\det(I-2t\sigma^{2}A)^{-1/2}
$$
が得られました.
残差平方和が従う分布
$\chi^{2}(m)$ の積率母関数は $(1-2t)^{-m/2}$ で与えられます.行列式は固有値の積で表されることと,冪等行列の性質から,$A=\sigma^{-2}P$ に関する積率母関数は
$$
M_{\varepsilon^{\top}A\varepsilon}(t)
=\det(I-2tP)^{-1/2}=(1-2t)^{-\operatorname{rank}{P}/2}
$$
となります.
$$
\operatorname{rank}{P}=\operatorname{tr}{P}=\operatorname{tr}{I}-\operatorname{tr}(X(X^{\top}X)^{-1}X^{\top})=n-p
$$
が成り立つので,積率母関数の一意性から,$\sigma^{-2}\|e\|^{2}=\sigma^{-2}\varepsilon^{\top}P\varepsilon\sim\chi^{2}(n-p)$ がわかりました.
検定統計量の導出
Cochran の定理
以下に示す Cochran の定理を用いて,$\sigma^{-2}\|e\|^{2}$ と $\sigma^{-1}(\hat\beta-\beta)$ が独立であることを示します.
Cochran の定理
$z\sim N(0,I)$ とする.冪等行列 $A$ と行列 $B$ について,$z^{\top}Az$ と $Bz$ が独立であるための必要十分条件は $AB^{\top}=O$ である.
$A:=P,B:=(X^{\top}X)^{-1}X^{\top}$ とします.このとき,
$$
AB^{\top}
=\left(I-X(X^{\top}X)^{-1}X^{\top}\right)X(X^{\top}X)^{-1}=O
$$
が成り立ちます.また,$\hat\beta-\beta=(X^{\top}X)^{-1}X^{\top}\varepsilon$ なので,Cochran の定理より,$\sigma^{-2}\|e\|^{2}$ と $\sigma^{-1}(\hat\beta-\beta)$ は独立であることがわかりました.
検定統計量の導出
$H^{\top}\beta=\xi$ を仮定して,検定統計量 $F_{0}$ を導出します.はじめに,$F_{0}$ の分子 $(H^{\top}\hat\beta - \xi )^\intercal {H^{\top}(X^{\top}X)^{-1}H}^{-1}( H^{\top}\hat\beta - \xi)$ が従う分布を導出します.
$$
\hat\beta=(X^{\top}X)^{-1}X^{\top}y=\beta+(X^{\top}X)^{-1}X^{\top}\varepsilon
$$
が成り立つので,$\hat\beta\sim N(\beta,\sigma^{2}(X^{\top}X)^{-1})$ です.
また,$\hat\xi:=H^{\top}\hat\beta\sim N(H^{\top}\beta,\sigma^{2}H^{\top}(X^{\top}X)^{-1}H)$ となります.$D:=H^{\top}(X^{\top}X)^{-1}H$ とすると,$\chi^{2}$ 分布と正規分布の関係から, $v_{1}:=(\hat\xi-\xi)^{\top}\sigma^{-2}D^{-1}(\hat\xi-\xi)\sim\chi^{2}(q)$ が従います.$v_{1}$ は $\sigma^{-1}(\hat\beta-\beta)$ の関数の形で表されるので,$v_{1}$ と $v_{2}:=\sigma^{-2}\|e\|^{2}$ は独立です.
最後に,F 分布の性質から,
F_{0}=\frac{v_{1}/q}{v_{2}/(n-p)}=\frac{n-p}{q}\frac{(\hat\xi-\xi)^{\top}D^{-1}(\hat\xi-\xi)}{\|e\|^{2}}\sim F(q,n-p)
が従います.
以上で,線形回帰の回帰係数の線形制約に関する検定統計量の導出ができました.
参考文献
佐和隆光 (1979). "回帰分析(新装版)". 朝倉書店.
付録
相関係数に関する検定統計量として,
$$
T:=\frac{|r|\sqrt{n-2}}{\sqrt{1-r^2}} \sim t(n-2)
$$
はよく知られています.ここで,$t(n-2)$ は自由度 $n-2$ の t 分布で,$r$ は標本相関係数です.この検定統計量が回帰係数に関する検定統計量として導出できることを示します.
$(x_{1},y_{1})^{\top},\dots,(x_{n},y_{n})^{\top}\overset{i.i.d}{\sim} N(\mu,\Sigma)$ とします.ただし,$\mu:=(\mu_{X},\mu_{Y})^{\top}\in\mathbb{R}^{2},\Sigma:=\begin{pmatrix}\sigma_{X}^{2}&\rho\sigma_{X}\sigma_{Y}\\\rho\sigma_{X}\sigma_{Y}&\sigma_{Y}^{2}\end{pmatrix},\sigma_{X}>0,\sigma_{Y}>0,\rho\in(0,1)$ です.
$$
y_{i}\mid x_{i}\sim N\left(\mu_{Y}+\rho\frac{\sigma_{Y}}{\sigma_{X}}(x_{i}-\mu_{X}),\left(\sigma_{Y}\sqrt{1-\rho^{2}}\right)^{2}\right)
$$
より,$\beta_{2}:=\rho\sigma_{Y}/\sigma_{X},\beta_{1}:=\mu_{Y}-\beta_{2}\mu_{X},\sigma:=\sigma_{Y}\sqrt{1-\rho^{2}}$ とすると,
$$
y\mid x\sim\beta_{1}\mathbb{1}+\beta_{2}x+\varepsilon
$$
と表せます.ただし,$\mathbb{1}$ は要素がすべて $1$ のベクトル,$\varepsilon\sim N(0,\sigma^{2}I)$ です.
$H_{0}:\beta_{2}=0$ に関する検定統計量から,$T$ を導出できることを示します.
まず,$H=(0,1)^{\top},\xi=0$ の場合の $F_{0}$ を求めます.$\overline{x}:=n^{-1}\sum_{i=1}^{n}x_{i},\overline{y}:=n^{-1}\sum_{i=1}^{n}y_{i}$ とします.
\begin{align*}
S_{XX}&:=n^{-1}\sum_{i=1}^{n}(x_{i}-\overline{x})^{2},\\
S_{YY}&:=n^{-1}\sum_{i=1}^{n}(y_{i}-\overline{y})^{2},\\
S_{XY}&:=n^{-1}\sum_{i=1}^{n}(x_{i}-\overline{x})(y_{i}-\overline{y})
\end{align*}
とすると,
\begin{align*}
(\hat\xi-\xi)^{\top}(\hat\xi-\xi)&=\hat\beta_{2}^{2}=\frac{S_{XY}^{2}}{S_{XX}^{2}},\\
D&=S_{XX}^{-1},\\
\|e\|^{2}&=\frac{S_{XY}^{2}}{S_{XX}}\frac{1-r^{2}}{r^{2}}
\end{align*}
と表せます(行間です.そこそこの量の計算を省略しました.).したがって,
\begin{align*}
F_{0}
&=\frac{n-p}{q}\frac{(\hat\xi-\xi)^{\top}D^{-1}(\hat\xi-\xi)}{\|e\|^{2}}\\
&=(n-2)\frac{S_{XY}^{2}}{S_{XX}^{2}}S_{XX}\frac{S_{XX}}{S_{XY}^{2}}\frac{r^{2}}{1-r^{2}}\\
&=(n-2)\frac{r^{2}}{1-r^{2}}\sim F(1,n-2)
\end{align*}
が得られます.F 分布と t 分布の関係から,
$$
T=\sqrt{F_0}\sim t(n-2)
$$
が成り立ちます.以上より,相関係数の検定統計量が導出できました.