残差平方和に関連する標本分布(カイ2乗分布) #統計検定

はじめに

統計検定1級の統計数理・理工学において，回帰分析の関連テーマ(最小二乗推定・線形単回帰・線形重回帰)は2012年～2022年の間で複数回登場した頻出分野の1つである。
本記事では，線形回帰モデルにおける3つの主な課題のうち，誤差項 $\boldsymbol{\varepsilon}$ が正規分布にしたがうという仮定を加えた状態で，残差平方和の標本分布に関する性質を整理する。
特に，射影子 $P_X = \textbf{X} (\textbf{X}^T \textbf{X})^{-1} \textbf{X}^T$ と正規直交基底との関係や，正準形に関する説明を行なう。

改訂履歴

特異値分解を用いたカイ2乗分布の導出が，久保川達也著　「現代数理統計学の基礎」に説明されていることについて，参考文献に追記した(2024/10/31)。

回帰分析の関連テーマ(最小二乗推定・線形単回帰・線形重回帰)

統計検定1級での出題

統計検定1級の統計数理・理工学において，2012年～2020年の10年間で，最小二乗推定は5回出題されている。また関連テーマである線形単回帰は4回，線形重回帰は3回出題されており，割と出題頻度が多いテーマである。

線形回帰モデル

線形回帰モデルは，

$\begin{align} \boldsymbol{y} = \textbf{X} \boldsymbol{\beta} + \boldsymbol{\varepsilon} \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$ は誤差項である。

線形回帰モデルにおける主な課題は以下の通りである。

最小二乗推定量とこれにかかわる数値の性質の把握
標本分布の性質の把握
検定

本記事では，2つ目の課題である「標本分布の性質の把握」について紹介する。

誤差項における正規性の仮定

実は，最小二乗推定量の導出においては，誤差項 $\boldsymbol{\varepsilon}$ がしたがう分布については仮定を置いていなかった。
誤差項 $\boldsymbol{\varepsilon}$ について正規性の仮定

$\begin{align} \boldsymbol{\varepsilon} \sim N(\boldsymbol{0}_n, \sigma^2 I_n) \end{align}$

を置くと，最小二乗推定量などがしたがう分布を導出することができるようになる。

残差平方和の標本分布 (カイ2乗分布)

残差平方和とは，残差ベクトルの絶対値の二乗のことである。残差平方和を $Q(\beta)$ で表すと，

$\begin{align} Q(\beta) = \boldsymbol{e}^T \boldsymbol{e} &= ( \boldsymbol{y} - \boldsymbol{\hat{y}} )^T( \boldsymbol{y} - \boldsymbol{\hat{y}} ) \\ &= ( \boldsymbol{y} - \textbf{X} \boldsymbol{\hat{\beta}})^T( \boldsymbol{y} - \textbf{X} \boldsymbol{\hat{\beta}} ) \end{align}$

となる。

では，この残差平方和がしたがう確率分布を考えてみよう。

射影子を用いた残差平方和の変形

残差平方和を少し式変形する。射影子 $P_X = \textbf{X} (\textbf{X}^T \textbf{X})^{-1} \textbf{X}^T$ を用いると，

$\begin{align} Q(\beta) &= ( \boldsymbol{y} - \boldsymbol{\hat{y}} )^T( \boldsymbol{y} - \boldsymbol{\hat{y}} ) \\ &= ( (I_n - \textbf{X} (\textbf{X}^T \textbf{X})^{-1} \textbf{X}^T ) \boldsymbol{y} )^T (I_n - \textbf{X} (\textbf{X}^T \textbf{X})^{-1} \textbf{X}^T ) \boldsymbol{y} \\ &= ( (I_n - P_X ) \boldsymbol{y} )^T (I_n - P_X) \boldsymbol{y} \\ &= \boldsymbol{y} ^T (I_n - P_X)^2 \boldsymbol{y} \\ &= \boldsymbol{y} ^T (I_n - P_X) \boldsymbol{y} \quad (\because (I-P_X)^2=I-P_X) \\ & \equiv \boldsymbol{y} ^T M_X \boldsymbol{y} \quad (I-P_X \equiv M_X) \end{align}$

となる。

更に変形して，

$\begin{align} Q(\beta) &= \boldsymbol{y} ^T M_X \boldsymbol{y} \\ &= (\textbf{X} \boldsymbol{\beta} + \boldsymbol{\varepsilon} )^T M_X (\textbf{X} \boldsymbol{\beta} + \boldsymbol{\varepsilon} ) \\ &= \boldsymbol{\varepsilon} ^T M_X \boldsymbol{\varepsilon} \quad (M_X \textbf{X} = (I - P_X) \textbf{X} = 0)\\ \end{align}$

となり，正規分布にしたがう誤差項 $\boldsymbol{\varepsilon}$ の二次形式になっている。

特異値分解を用いた射影子の性質の確認

ここで，射影子 $P_X$ および $M_X$ の性質を，特異値分解を用いて確認していく。

$\textbf{X}$ を以下のように特異値分解する。

$\begin{align} \textbf{X} &= \textbf{U} \textbf{S} \textbf{V}^T \end{align}$

ただし， $\textbf{U} \in \mathbb{R}^{n \times d}$ は $d$ 個の $n$ 次元正規直交基底からなる行列， $\textbf{V} \in \mathbb{R}^{n \times d}$ は $d$ 個の $d$ 次元正規直交基底からなる行列， $\textbf{S}$ は $d$ 個の特異値を対角成分に持つ対角行列である。

$\begin{align} \textbf{U} &= \begin{bmatrix} \boldsymbol{u}_1 & \boldsymbol{u}_2 & \cdots & \boldsymbol{u}_d \\ \end{bmatrix} \in \mathbb{R}^{n \times d} \\ \textbf{V} &= \begin{bmatrix} \boldsymbol{v}_1 & \boldsymbol{v}_2 & \cdots & \boldsymbol{v}_d \\ \end{bmatrix} \in \mathbb{R}^{d \times d} \\ \textbf{S} &= \mathrm{diag} [ \sigma_1 \cdots \sigma_d ] \in \mathbb{R}^{d \times d} \end{align}$

$\textbf{U}, \textbf{V}$ は上記のように正規直交基底からなる行列なので， $\textbf{U}^T \textbf{U} = I_d$ および $\textbf{V}^T \textbf{V} = I_d$ を満たす。

$\textbf{X}$ の特異値分解を用いると，

$\begin{align} \textbf{X}^T \textbf{X} = \textbf{V} \textbf{S} \textbf{U}^T \cdot \textbf{U} \textbf{S} \textbf{V}^T = \textbf{V} \textbf{S}^2 \textbf{V}^T \in \mathbb{R}^{d \times d} \end{align}$

および

$\begin{align} (\textbf{X}^T \textbf{X})^{-1} = \textbf{V} \textbf{S}^{-2} \textbf{V}^T = \textbf{V} \begin{bmatrix} 1/\sigma_1^2 & & \\ & \ddots & \\ & & 1/\sigma_d^2 \\ \end{bmatrix} \textbf{V}^T \end{align}$

なので，

$\begin{align} \textbf{X} (\textbf{X}^T \textbf{X})^{-1} \textbf{X}^T &= \textbf{U} \textbf{S} \textbf{V}^T \cdot \textbf{V} \textbf{S}^{-2} \textbf{V}^T \cdot \textbf{V} \textbf{S} \textbf{U}^T \\ &= \textbf{U} \textbf{S} \textbf{S}^{-2} \textbf{S} \textbf{U}^T \\ &= \textbf{U} \textbf{U}^T \in \mathbb{R}^{n \times n} \end{align}$

が得られる。すなわち射影子 $P_X$ は，正規直交基底ベクトルの積になるということが分かる。

また $\textbf{U} \textbf{U}^T$ は，

$\begin{align} \textbf{U} \textbf{U}^T &= \begin{bmatrix} \boldsymbol{u}_1 & \boldsymbol{u}_2 & \cdots & \boldsymbol{u}_d \\ \end{bmatrix} \begin{bmatrix} \boldsymbol{u}_1^T \\ \boldsymbol{u}_2^T \\ \vdots \\ \boldsymbol{u}_d^T \\ \end{bmatrix} \\ &= \sum_{i=1}^d \boldsymbol{u}_i \boldsymbol{u}_i ^T \end{align}$

となる。

$\textbf{U}$ は $d$ 個の正規直交基底からなるが，これらと直交するような正規直交基底ベクトル $\boldsymbol{u}_{d+1},..., \boldsymbol{u}_n$ を追加した行列 $\textbf{G}$ を考える。
すなわち，

$\begin{align} \textbf{G} &= \begin{bmatrix} \boldsymbol{u}_1 & \boldsymbol{u}_2 & \cdots & \boldsymbol{u}_d & | & \boldsymbol{u}_{d+1} \cdots & \boldsymbol{u}_n\\ \end{bmatrix} \end{align}$

である。これを用いると，

$\begin{align} \textbf{G} \textbf{G}^T &= \sum_{i=1}^d \boldsymbol{u}_i \boldsymbol{u}_i ^T + \sum_{i=d+1}^n \boldsymbol{u}_i \boldsymbol{u}_i ^T \\ &= P_X + \sum_{i=d+1}^n \boldsymbol{u}_i \boldsymbol{u}_i ^T \\ &= I_n \end{align}$

よって，

$\begin{align} M_X &= I_n - P_X = \sum_{i=d+1}^n \boldsymbol{u}_i \boldsymbol{u}_i ^T \end{align}$

となり， $M_X$ もまた，正規直交基底ベクトルの積になるということが分かる。

正準形への変形

前節で計算した正規直交基底 $\textbf{G}$ を用いて，誤差項 $\boldsymbol{\varepsilon}$ の見通しを良くするために変形する。

$\begin{align} \boldsymbol{z} = \textbf{G}^T \boldsymbol{\varepsilon} \end{align}$

で定義される $\boldsymbol{z}$ を導入する。 $\boldsymbol{z}$ は，正規分布にしたがう確率変数 $\boldsymbol{\varepsilon}$ に正規直交基底 $\textbf{G}$ をかけたものなので， $\boldsymbol{z}$ も正規分布にしたがう。

$\begin{align} E [ \boldsymbol{z} ] &= \textbf{G}^T E [ \boldsymbol{\varepsilon} ] = \boldsymbol{0} _n \\ V [ \boldsymbol{z} ] &= \textbf{G}^T V [ \boldsymbol{\varepsilon} ] \textbf{G} = \sigma^2 I_n \\ \end{align}$

であるので，

$\begin{align} \boldsymbol{z} &\sim N(\boldsymbol{0}_n, \sigma^2 I_n) \\ z_i &\sim N(0, \sigma^2) \quad (i=1, ..., n) \end{align}$

となる。このように $\textbf{X}$ が張る空間上の正規直交基底と，これに直交する空間上の正規直交基底で変換した形を正準形と呼ぶ。

さらに， $\boldsymbol{z} = \textbf{G}^T \boldsymbol{\varepsilon}$ より，

$\begin{align} \begin{bmatrix} z_1 \\ \vdots \\ z_n \\ \end{bmatrix} &= \begin{bmatrix} \boldsymbol{u}_1^T \\ \vdots \\ \boldsymbol{u}_n^T \\ \end{bmatrix} \boldsymbol{\varepsilon} = \begin{bmatrix} \boldsymbol{u}_1^T \boldsymbol{\varepsilon} \\ \vdots \\ \boldsymbol{u}_n^T \boldsymbol{\varepsilon} \\ \end{bmatrix} \\ \end{align}$

よって，

$\begin{align} \therefore z_i &= \boldsymbol{u}_i^T \boldsymbol{\varepsilon} \end{align}$

となる。

さらに，

$\begin{align} z_i^2 = (\boldsymbol{u}_i^T \boldsymbol{\varepsilon})^T \boldsymbol{u}_i^T \boldsymbol{\varepsilon} = \boldsymbol{\varepsilon}^T \boldsymbol{u}_i \boldsymbol{u}_i^T \boldsymbol{\varepsilon} \end{align}$

となるので，

$\begin{align} \sum_{i=d+1}^{n} z_i^2 &= \sum_{i=d+1}^{n} \boldsymbol{\varepsilon}^T \boldsymbol{u}_i \boldsymbol{u}_i^T \boldsymbol{\varepsilon} \\ &= \boldsymbol{\varepsilon}^T M_X \boldsymbol{\varepsilon} \\ &= \boldsymbol{e}^T \boldsymbol{e} \end{align}$

となり，残差平方和 $Q(\beta) = \boldsymbol{e}^T \boldsymbol{e}$ が，正準形で表現されることが確認できた。

カイ2乗分布の導出

$z_i \sim N(0, \sigma^2)$ より， $(z_i/\sigma) \sim N(0, 1)$ となるので，

$\begin{align} \sum_{i=d+1}^{n} \left( \frac{z_i}{\sigma} \right) ^2 \sim \chi^2 (n - d) \end{align}$

となる。よって，

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

となり，残差平方和は自由度 $(n-d)$ のカイ2乗分布にしたがうということが確認できた。

誤差項の分散の不偏推定量

残差平方和は自由度 $(n-d)$ のカイ2乗分布にしたがうので，誤差項 $\boldsymbol{\varepsilon} \sim N(\boldsymbol{0}_n, \sigma^2 I_n)$ における $\sigma^2$ の不偏推定量を求めてみよう。

自由度 $(n-d)$ のカイ2乗分布の平均は $(n-d)$ なので，

$\begin{align} E \left [ \frac{\boldsymbol{e}^T \boldsymbol{e} }{\sigma^2} \right ] &= n - d \\ \therefore E \left [ \frac{\boldsymbol{e}^T \boldsymbol{e} }{n -d} \right ] &= \sigma^2 \end{align}$

となり， $\sigma^2$ の不偏推定量が得られた。