残差回帰とその応用テクニック（1）

# はじめに
千葉大学・株式会社Nospareの[川久保](https://sites.google.com/site/ykawakubostat/)です．今回と次回で，線形回帰モデルにおける残差回帰（residual regression，またはpartitioned regression）と呼ばれるテクニックを紹介します．残差回帰とは，最小二乗推定量の部分ベクトルを求めるためのテクニックです．残差回帰のテクニックを用いると，回帰係数のF検定統計量の導出や，私が以前に紹介した[パネルデータ分析](https://qiita.com/ykawakubo/items/0dcfacaa00188b48737d)における最小二乗ダミー変数（LSDV）推定量が，固定効果推定量やWG推定量と一致することを示せます．


今回は残差回帰の説明と，そこから導かれる最小二乗推定量の部分ベクトルの分布に関する性質を紹介します．次回は，今回紹介したテクニックを用いた応用テクニックを紹介します．


# 残差回帰とFWL定理
$n\times 1$目的変数ベクトル$y$と，$n\times k$説明変数行列について，以下の線形回帰モデルを考えます．

$$
\tag{1}
y = X \beta + \varepsilon, \quad \varepsilon \sim \mathrm{N}(0,\sigma^2 I_n)
$$

ここで，$X = [X_1: X_2]$と分割され，$X_1$が$n\times k_1$，$X_2$が$n\times k_2$，$k = k_1 + k_2$とします．すなわち，モデル(1)が

$$
y = X_1\beta_1 + X_2\beta_2 + \varepsilon,
$$
と書けるとします．ここで$\beta = [\beta_1^\top, \beta_2^\top]^\top$です．モデル(1)の最小二乗推定量を$\hat{\beta}$と書き，対応する部分ベクトルの分割を$\hat{\beta} = [\hat{\beta}_1^\top, \hat{\beta}_2^\top]^\top$と書くとします．

残差回帰とは，$\hat{\beta}_2$を以下の3つのステップから求める方法です．

1. $y$を$X_1$に回帰して残差ベクトル$\tilde{y}$を求める．
1. $X_2$を$X_1$に回帰して，残差行列$\tilde{X}_2$を求める．
1. $\tilde{y}$を$\tilde{X}_2$に回帰する．

上記のステップ3で求まった回帰係数の推定量は，モデル(1)の最小二乗推定量$\hat{\beta}$の部分ベクトル$\hat{\beta}_2$と等しくなり，FWL定理（Frisch-Waugh-Lovell定理）として知られています．


# 残差回帰から求まる推定量のかたち
残差回帰の3つのステップを，1つ1つ数式を追ってみたいと思います．

### ステップ1
まずステップ1ですが，$y$を$X_1$に回帰させたときの最小二乗推定量は$(X_1^\top X_1)^{-1} X_1^\top y$であることに気をつけると，
$$
\tilde{y} = y - X_1 (X_1^\top X_1)^{-1} X_1^\top y = (I_n - P_1) y,
$$
と表されます．ただし，$I_n$は$n\times n$単位行列，$P_1 = X_1 (X_1^\top X_1)^{-1} X_1^\top$は$X_1$の列空間への射影行列です．射影行列については，私の[以前の記事](https://qiita.com/ykawakubo/items/dd352bc0921f4c1047b2)を参考にしてください．

### ステップ2
$X_2$の各列を$X_1$に回帰させたときの回帰係数ベクトル（$k_1 \times 1$）を，$k_2$列並べた行列（$k_1 \times k_2$）は，
$$
(X_1^\top X_1)^{-1} X_1^\top X_2
$$
であることから，残差ベクトルは，
$$
\tilde{X}_2 = X_2 - X_1(X_1^\top X_1)^{-1} X_1^\top X_2 = (I_n - P_1)X_2
$$
と表されます．

### ステップ3
$\tilde{y}$を$\tilde{X}_2$に回帰したときの最小二乗推定量を$\tilde{\beta}_2$とすると，

```math
\tag{2}
\begin{align}
\tilde{\beta}_2 &= (\tilde{X}_2^\top \tilde{X}_2)^{-1} \tilde{X}_2^\top \tilde{y} \\
&= \{ X_2^\top (I_n - P_1)X_2 \}^{-1} X_2^\top (I_n - P_1) y
\end{align}
```

と表されます．上式の2つ目の等式には，$I_n - P_1$が射影行列であり，対称なベキ等行列であるという性質を用いています．この$\tilde{\beta}_2$が，モデル(1)の最小二乗推定量の部分ベクトル$\hat{\beta}_2$と等しくなるというのが，FWL定理です．


# FWL定理の証明
モデル(1)の最小二乗推定量$\hat{\beta}$は，$X^\top X \hat{\beta} = X^\top y$を満たします．この式を，$X = [X_1:X_2]$，$\hat{\beta} = [\hat{\beta}_1^\top, \hat{\beta}_2^\top]^\top$の分割に注意して書き直すと，
```math
\tag{3}
\begin{bmatrix}
X_1^\top X_1 & X_1^\top X_2 \\
X_2^\top X_1 & X_2^\top X_2
\end{bmatrix}
\begin{bmatrix}
\hat{\beta}_1 \\
\hat{\beta}_2
\end{bmatrix}
=
\begin{bmatrix}
X_1^\top y \\
X_2^\top y
\end{bmatrix}
```
となります．(3)式の1行目を$\hat{\beta}_1$について解くと，

```math
\begin{align}
\hat{\beta}_1 &= (X_1^\top X_1)^{-1}X_1^\top y - (X_1^\top X_1)^{-1} X_1^\top X_2 \hat{\beta}_2 \\
&= (X_1^\top X_1)^{-1} X_1^\top (y - X_2\hat{\beta}_2 )
\end{align}
```

となります．これを(3)式の2行目
$$
X_2^\top X_1 \hat{\beta}_1 + X_2^\top X_2 \hat{\beta}_2 = X_2^\top y
$$
に代入すると，
$$
X_2^\top P_1 y - X_2^\top P_1 X_2 \hat{\beta}_2 + X_2^\top X_2 \hat{\beta}_2 = X_2^\top y
$$
となり，これを$\hat{\beta}_2$について解くと，

```math
\begin{align}
\hat{\beta}_2 &= \{ X_2^\top (I_n - P_1) X_2 \}^{-1} X_2^\top (I_n - P_1)y \\
&= \tilde{\beta}_2
\end{align}
```

となり，たしかに(2)式の残差回帰推定量$\tilde{\beta}_2$と一致しました．


# 最小二乗推定量の部分ベクトルの分布
最小二乗推定量の部分ベクトル$\hat{\beta}_2$が，(2)式の残差回帰推定量$\tilde{\beta}_2$と一致することから，以下のような$\hat{\beta}_2$の分布の性質を導くことができます．


まず$\hat{\beta}_2$の平均は，

```math
\begin{align}
\hat{\beta}_2 &= (\tilde{X}_2^\top \tilde{X}_2)^{-1} \tilde{X}_2^\top \tilde{y} \\
&= (\tilde{X}_2^\top \tilde{X}_2)^{-1} X_2^\top (I_n - P_1) y \\
&= (\tilde{X}_2^\top \tilde{X}_2)^{-1} X_2^\top (I_n - P_1)(X_1\beta_1 + X_2\beta_2 + \varepsilon) \\
&= (\tilde{X}_2^\top \tilde{X}_2)^{-1} X_2^\top (I_n - P_1) (X_2\beta_2 + \varepsilon) \quad (\because (I_n - P_1)X_1 = 0) \\
&= (\tilde{X}_2^\top \tilde{X}_2)^{-1} \tilde{X}_2^\top \tilde{X}_2 \beta_2 + (\tilde{X}_2^\top \tilde{X}_2)^{-1} \tilde{X}_2^\top \varepsilon \\
&= \beta_2 + (\tilde{X}_2^\top \tilde{X}_2)^{-1} \tilde{X}_2^\top \varepsilon \tag{4}
\end{align}
```

という変形から，$E[\varepsilon] = 0$より，$E[\hat{\beta}_2] = \beta_2$です．もちろん，$\hat{\beta}_2$は最小二乗推定量の部分ベクトルで，最小二乗推定量は不偏推定量であることからも当然の結果です．

次に分散共分散行列は，

```math
\begin{align}
V(\hat{\beta}_2) &= E[(\hat{\beta}_2 - \beta_2)(\hat{\beta}_2 - \beta_2)^\top] \\
&= (\tilde{X}_2^\top \tilde{X}_2)^{-1} \tilde{X}_2^\top E[\varepsilon \varepsilon^\top] \tilde{X}_2 (\tilde{X}_2^\top \tilde{X}_2)^{-1} \quad (\because (4)) \\
&= \sigma^2 (\tilde{X}_2^\top \tilde{X}_2)^{-1} \tilde{X}_2^\top \tilde{X}_2 (\tilde{X}_2^\top \tilde{X}_2)^{-1} \quad (\because E[\varepsilon \varepsilon^\top] = V(\varepsilon) = \sigma^2 I_n) \\
&= \sigma^2 (\tilde{X}_2^\top \tilde{X}_2)^{-1}
\end{align}
```

と求めることができます．よって，$\varepsilon$に正規分布を仮定しているとき，
$$
\hat{\beta}_2 \sim \mathrm{N}(\beta_2, \sigma^2(\tilde{X}_2^\top \tilde{X}_2)^{-1})
$$
が成り立ちます．さらに，このことから，

```math
\begin{align}
& (\hat{\beta}_2 - \beta_2)^\top V(\hat{\beta}_2)^{-1} (\hat{\beta}_2 - \beta_2) \\
&= (\hat{\beta}_2 - \beta_2)^\top \frac{\tilde{X}_2^\top \tilde{X}_2}{\sigma^2} (\hat{\beta}_2 - \beta_2) \sim \chi_{k_2}^2
\end{align}
```

すなわち自由度$k_2$（$\beta_2$の次元）のカイ二乗分布にしたがうことがわかります．この性質は，次回紹介するF検定の導出で用います．


# おわりに

株式会社Nospareには，統計学の様々な分野を専門とする研究者が所属しております．統計アドバイザリーやビジネスデータの分析につきましては[株式会社Nospare](https://bit.ly/38WtGZH) までお問い合わせください．
<a href="https://bit.ly/38WtGZH">
<img src="https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1023506/d5c625c0-d6f4-680b-d702-3c40f11d6c39.jpeg" alt="株式会社Nospare" width=400>
</a>