はじめに

千葉大学・株式会社Nospareの川久保です．前回の記事で，残差回帰と呼ばれる最小二乗推定量の部分ベクトルを求めるテクニックを紹介し，その分布について解説しました．今回は，前回紹介したテクニックを応用して，回帰係数のF検定統計量の導出や，パネルデータ分析における最小二乗ダミー変数（LSDV）推定量が，固定効果推定量やWG推定量と一致することを示します．

F検定

記法の準備

最小二乗推定量$\hat{\beta} = [\hat{\beta_1}^\top, \hat{\beta}_2^\top]^\top$に対して，その残差ベクトルを

\begin{align}
e &= y - X_1 \hat{\beta}_1 - X_2 \hat{\beta}_2 \\
&= y - X \hat{\beta} \\
&= (I_n - P)y
\end{align}

と書きます．ただし，$P = X(X^\top X)^{-1}X^\top$は$X$の列空間への射影行列です．繰り返しの宣伝になりますが，射影行列については，私の以前の記事を参考にしてください．

誤差分散$\sigma^2$の不偏推定量は，

\begin{align}
\hat{\sigma}^2 &= \frac{1}{n-k}(y - X\hat{\beta}) \\
&= \frac{1}{n-k}e^\top e,
\end{align}

すなわち残差二乗和を$n-k$で割った形です．ただし$k$は説明変数の数（$X$の列数）です．また，誤差項が正規分布にしたがう場合（この記事ではこれを仮定します），
$$
\tag{1}
\frac{e^\top e}{\sigma^2} = \frac{(n-k)\hat{\sigma}^2}{\sigma^2} \sim \chi_{n-k}^2
$$
が成り立ち，$\hat{\sigma}^2$と$\hat{\beta}$は独立です．

F検定統計量の導出

以下の仮説検定を考えます．
$$
H_0: \beta_2 = 0, \quad H_1: \beta_2 \not= 0
$$
このとき，帰無仮説$H_0$のもとで，以下が成り立ちます．

\tag{2}
F = \frac{\hat{\beta}_2^\top \tilde{X}_2^\top \tilde{X}_2 \hat{\beta}_2 \ / \ k_2}{\hat{\sigma}^2} \sim F_{k_2, n-k}

このことを以下で簡単に説明します．まず，前回の記事の最後で紹介しましたが，$(\hat{\beta}_2 - \beta_2)^\top \tilde{X}_2^\top \tilde{X}_2 (\hat{\beta}_2 - \beta_2) \ / \ \sigma^2 \sim \chi_{k_2}^2$が成り立ちます．よって帰無仮説$H_0$のもとで，(2)式の分子を$\sigma^2$で割った形
$$
\frac{\hat{\beta}_2^\top \tilde{X}_2^\top \tilde{X}_2 \hat{\beta}_2}{\sigma^2 k_2}
$$
は自由度$k_2$のカイ二乗分布にしたがう確率変数を$k_2$で割った形であることが分かります．次に，分母を$\sigma^2$で割った形
$$
\frac{\hat{\sigma}^2}{\sigma^2} = \frac{(n-k)\hat{\sigma}^2}{\sigma^2} \frac{1}{n-k}
$$
は，(1)から自由度$n-k$のカイ二乗分布にしたがう確率変数を$n-k$で割った形です．さらに$\hat{\sigma}^2$と$\hat{\beta}_2$は独立であることから，(2)が成り立ちます．

F検定統計量の残差二乗和を使った表現

F検定統計量の別表現を導出します．帰無仮説$H_0$のもと（$\beta_2 = 0$という制約を課したモデル）での残差は，

\begin{align}
y - X_1(X_1^\top X_1)^{-1} X_1^\top y &= (I_n - P_1)y
\end{align}

です．ここで，$y = X_1\hat{\beta}_1 + X_2\hat{\beta}_2 + e$（$\hat{\beta}_1, \hat{\beta}_2, e$は，制約なしモデルの回帰係数の推定量と残差ベクトル）なので，

\begin{align}
(I_n - P_1)y &= (I_n - P_1)(X_1\hat{\beta}_1 + X_2\hat{\beta}_2 + e) \\
&= (I_n - P_1)X_2\hat{\beta}_2 + (I_n - P_1)e \quad (\because (I_n - P_1)X_1 = 0) \\
&= \tilde{X}_2\hat{\beta}_2 + (I_n - P_1)(I_n - P)y. \quad (\because (I_n - P_1)X_2 = \tilde{X}_2)
\end{align}

最後の等式の第二項については，

\begin{align}
(I_n - P_1)(I_n - P) &= I_n - P_1 - P + P_1P \\
&= I_n - P_1 - P + P_1 \\
&= I_n - P \tag{3}
\end{align}

となりますが，$P_1P = P_1$は，$P$が$X$の列空間への射影行列，$P_1$が$X_1$の列空間への射影行列，そして$X_1$の列空間は$X$の列空間の部分空間であることから成り立ちます（より小さい空間への射影になるということ）．よって，制約つきモデルの残差は，

\begin{align}
(I_n - P_1)y &= \tilde{X}_2\hat{\beta}_2 + (I_n - P)y \\
&= \tilde{X}_2 \hat{\beta}_2 + e
\end{align}

と表されます．

制約つきモデルの残差二乗和（residual sum of squares, RSS）は，

\begin{align}
\mathrm{RSS}_\mathrm{R} &= (\tilde{X}_2\hat{\beta}_2 + e)^\top (\tilde{X}_2\hat{\beta}_2 + e) \\
&= \hat{\beta}_2^\top \tilde{X}_2^\top \tilde{X}_2 \hat{\beta}_2 + e^\top e +2e^\top \tilde{X}_2\hat{\beta}_2
\end{align}

となりますが，第三項（交差項）については，

\begin{align}
e^\top \tilde{X}_2 &= y^\top (I_n - P) (I_n - P_1) X_2 \\
&= y^\top (I_n - P)X_2 \quad (\because (3)) \\
&= 0
\end{align}

より0となります．$e^\top e$は制約なしモデルの残差二乗和$\mathrm{RSS}_\mathrm{UR}$なので，

\hat{\beta}_2^\top \tilde{X}_2^\top \tilde{X}_2 \hat{\beta}_2 = \mathrm{RSS}_\mathrm{R} - \mathrm{RSS}_\mathrm{UR}

が成り立ちます．さらに，$\hat{\sigma}^2 = \mathrm{RSS}_\mathrm{UR} \ / \ (n-k)$より，F検定統計量(2)は，以下のように書き直すことができます．

F = \frac{(\mathrm{RSS}_\mathrm{R} - \mathrm{RSS}_\mathrm{UR}) \ / \ k_2}{\mathrm{RSS}_\mathrm{UR} \ / \ (n-k)}

この形の方が教科書などでおなじみかと思います．

定数項のみモデルのF検定

上記の検定問題において，$X_1 = \mathbb{1}_n$（1が$n$個縦に並んだベクトル）である場合，すなわち「定数項を除いたその他の説明変数にかかる回帰係数が全て0である」という帰無仮説を検定する問題を考えます．このとき，制約つきモデル（定数項のみモデル）の回帰係数の最小二乗推定量は，

\begin{align}
(X_1^\top X_1)^{-1} X_1^\top y &= (\mathbb{1}_n^\top \mathbb{1}_n)^{-1} \mathbb{1}_n^\top y \\
&= \frac{1}{n}(y_1+\dots+y_n) \\
&=: \bar{y},
\end{align}

すなわち$y$の標本平均になります．よって，制約つきモデルの残差二乗和は，

\begin{align}
\mathrm{RSS}_\mathrm{R} &= (y - \mathbb{1}_n\bar{y})^\top (y - \mathbb{1}_n\bar{y}) \\
&= \sum_{i=1}^n (y_i - \bar{y})^2
\end{align}

となり，F検定統計量は，

F = \frac{ \{ \sum_{i=1}^n (y_i - \bar{y})^2 - \sum_{i=1}^n (y_i - x_i^\top \hat{\beta})^2 \} \ / \ k_2 }{ \sum_{i=1}^n (y_i - x_i^\top \hat{\beta})^2 \ / \ (n-k_2 - 1) }

と表されます．

LSDV推定量の導出

パネルデータ分析における，固定効果モデルの3つの推定量は全て一致することを以前の記事で紹介しました．3つの推定量のうちの1つである最小二乗ダミー変数（LSDV）推定量は，残差回帰のテクニックを用いて導出することができます．

固定効果モデル

$\{ y_{it}, x_{it} \}$を，$i$番目の個体（$i=1,\dots,N$）の$t$期目（$t=1,\dots,T$）の観測とし，
$$
y_{it} = x_{it}^\top \beta + b_i + \varepsilon_{it}, \quad \varepsilon_{it} \overset{\mathrm{iid}}{\sim} (0,\sigma^2) \quad (i = 1,\dots,N; \ t = 1,\dots,T),
$$
なるモデルを考えます．$X_i = (x_{i1},\dots,x_{iT})^\top, y_i = (y_{i1},\dots,y_{iT})^\top, \varepsilon_i = (\varepsilon_{i1},\dots,\varepsilon_{iT})^\top$とすると，このモデルは個体$i$ごとに

y_i = X_i \beta + \mathbb{1}_T b_i + \varepsilon_i, \quad (i=1,\dots,N),

と書け，さらに$y = (y_1^\top,\dots,y_N^\top)^\top, X = (X_1^\top,\dots,X_N^\top)^\top, b = (b_1,\dots,b_N)^\top, \varepsilon = (\varepsilon_1^\top,\dots,\varepsilon_N^\top)^\top$とし，また$Z = I_N \otimes \mathbb{1}_T$とおくと，

\begin{bmatrix}
y_1 \\
\vdots \\
y_N
\end{bmatrix}
=
\begin{bmatrix}
X_1 \\
\vdots \\
X_N
\end{bmatrix}
\beta +
\begin{bmatrix}
\mathbb{1}_T & & \\
 & \ddots & \\
 & & \mathbb{1}_T
\end{bmatrix}
\begin{bmatrix}
b_1 \\
\vdots \\
b_N
\end{bmatrix}
+
\begin{bmatrix}
\varepsilon_1 \\
\vdots \\
\varepsilon_N
\end{bmatrix},

つまり

\tag{4}
y = X\beta + Zb + \varepsilon,

と行列表記できます．

残差回帰を用いたLSDV推定量の導出

LSDV推定では，(4)式のモデルを，

y = [X:Z]
\begin{bmatrix}
\beta \\
b
\end{bmatrix}
+ \varepsilon,

という大きな線形回帰モデルだと見て，$b$を$N$個のダミー変数$Z$にかかる回帰係数だとみなします．回帰係数ベクトル$(\beta^\top, b^\top)^\top$を最小二乗推定し，最小二乗推定量の$\beta$に関する部分ベクトルをLSDV推定量と呼んでいます．そこで，最小二乗推定量の部分ベクトルを求める残差回帰のテクニックが用いられるわけです．

(4)式における$Zb$の項を，残差回帰の解説における$X_1\beta_1$，(4)式における$X\beta$の項を，残差回帰の解説における$X_2\beta_2$だとみなして，残差回帰を実行します．$Z$の列空間への射影行列を$P_Z = Z(Z^\top Z)^{-1}Z^\top$とすると，

\hat{\beta} = \{ X^\top (I_{NT} - P_Z) X \}^{-1} X^\top (I_{NT} - P_Z)y

が，$(\beta^\top, b^\top)^\top$の最小二乗推定量の部分ベクトル（つまり$\beta$のLSDV推定量）です．これが，固定効果推定量およびWG推定量に一致することを以下で示していきます．$P_Z$は以下のように変形できます．

\begin{align}
P_Z &= (I_N \otimes \mathbb{1}_T) \{ (I_N \otimes \mathbb{1}_T)^\top (I_N \otimes \mathbb{1}_T) \}^{-1} (I_N \otimes \mathbb{1}_T)^\top \\
&= (I_N \otimes \mathbb{1}_T) (I_N \otimes \mathbb{1}_T^\top \mathbb{1}_T)^{-1} (I_N \otimes \mathbb{1}_T^\top) \\
&= I_N \otimes \mathbb{1}_T(\mathbb{1}_T^\top \mathbb{1}_T)^{-1}\mathbb{1}_T^\top \\
&= I_N \otimes (T^{-1}J_T)
\end{align}

ここで，$J_T$は，要素が全て1の$T\times T$正方行列です．さらに，

\begin{align}
I_{NT} - P_Z &= I_N \otimes (I_T - T^{-1}J_T) \\
&= I_N \otimes Q
\end{align}

となります．ただし，$Q = I_T - T^{-1}J_T$はWithin-Group変換（各観測から個体内標本平均を引く操作）を施す行列です．よって，

\begin{align}
\hat{\beta} &= \{ X^\top (I_N \otimes Q) X \}^{-1} X^\top (I_N \otimes Q) y \\
&= \left\{
\begin{bmatrix}
X_1^\top, \cdots, X_N^\top
\end{bmatrix}
\begin{bmatrix}
Q & & \\
 & \ddots & \\
 & & Q
\end{bmatrix}
\begin{bmatrix}
X_1 \\
\vdots \\
X_N
\end{bmatrix}
\right\}^{-1}
\begin{bmatrix}
X_1^\top, \cdots, X_N^\top
\end{bmatrix}
\begin{bmatrix}
Q & & \\
 & \ddots & \\
 & & Q
\end{bmatrix}
\begin{bmatrix}
y_1 \\
\vdots \\
y_N
\end{bmatrix}
\\
&= \left( \sum_{i=1}^N X_i^\top Q X_i \right)^{-1} \left( \sum_{i=1}^N X_i^\top Q y \right) \\
&= \left\{ \sum_{i=1}^N \sum_{t=1}^T (x_{it} - \bar{x}_i)(x_{it} - \bar{x}_i)^\top \right\}^{-1} \sum_{i=1}^N \sum_{t=1}^T (x_{it} - \bar{x}_i)(y_{it} - \bar{y}_i)
\end{align}

と式変形でき，LSDV推定量が固定効果推定量およびWG推定量と一致することが示されました．

おわりに

株式会社Nospareには，統計学の様々な分野を専門とする研究者が所属しております．統計アドバイザリーやビジネスデータの分析につきましては株式会社Nospare までお問い合わせください．