はじめに
千葉大学・株式会社Nospareの川久保です．今回と次回で，計量経済学の分野で頻出するパネルデータ分析法の超入門編を書きたいと思います．今回は，パネルデータに対する基本的なモデルと推定法を紹介します．次回は，各種の検定法と，データ解析例を紹介する予定です．この記事で紹介する種々の推定量の名称や例は，千木良弘朗・早川和彦・山本拓「動学的パネルデータ分析」を参考にしました．
パネルデータとは
パネルデータとは，いくつかの個体（家計，企業など）について複数期間にわたって観測されたデータのことをいいます．記法としては，$\{ y_{it}, x_{it} \}$を，$i$番目の個体（$i=1,\dots,N$）の$t$期目（$t=1,\dots,T$）の観測として表します．このようなデータ構造は，経時データとも呼ばれますが，慣習として計量経済学ではパネルデータ，生物・医学統計では経時データと呼ばれている印象です．
そして，経時データに対する統計モデルとしては，$x_{it}$では捕らえられない個体$i$の個別効果を変量効果（random effect）として回帰モデルに入れた，線形混合モデルがよく用いられます．私の過去の記事「線形混合モデルとその応用例」もご参照ください．
パネルデータ分析の3つの基本的なモデル
パネルデータに対する基本的な分析法として，以下の3つが考えられます．
- Pooled OLS
- 変量効果モデル
- 固定効果モデル
それぞれ以下で説明していきます．
Pooled OLS
これは，パネルデータの構造を使わず，サンプルサイズが$NT$のデータだと思って，通常の線形回帰モデルを当てはめる手法です．つまり，
$$
\tag{1}
y_{it} = x_{it}^\top \beta + \varepsilon_{it}, \quad \varepsilon_{it} \overset{\mathrm{iid}}{\sim} (0,\sigma^2), \quad (i=1,\dots,N; \ t = 1,\dots,T)
$$
というモデルを想定します．ここで$\varepsilon_{it}$は独立同一に平均$0$，分散$\sigma^2$の分布にしたがうとし，以下最小二乗法などのモーメント法で推定するため，パラメトリック分布の仮定はおきません．(1)式の回帰係数を，最小二乗法（ordinary least square, OLS）で以下のように推定します．
\hat{\beta}_{\mathrm{POLS}} = \left( \sum_{i=1}^N \sum_{t=1}^T x_{it} x_{it}^\top \right)^{-1} \sum_{i=1}^N \sum_{t=1}^T x_{it} y_{it}
OLSはあくまで推定法ですが，慣習としてPooled OLSをモデル名のように言ってしまうことが多いです．
変量効果モデル
(1)のモデルに，$x_{it}$では捕らえられない個体$i$の個別効果を変量効果として入れた
$$
\tag{2}
y_{it} = x_{it}^\top \beta + b_i + \varepsilon_{it}, \quad b_i \overset{\mathrm{iid}}{\sim} (0,\tau^2), \quad \varepsilon_{it} \overset{\mathrm{iid}}{\sim} (0,\sigma^2)
$$
というモデルを，変量効果モデルと呼びます．お気づきの方もおられると思いますが，これはまさに線形混合モデルのrandom intercept modelのことです．計量経済学では変量効果モデルという名称が一般的なのですが，それは(2)式の$b_i$を変量効果としてでなく，固定効果（fixed effect）とするモデルも存在するからです．
変量効果モデルに対しては，$u_{it} = b_i + \varepsilon_{it}$という誤差項だと思い，
$$
y_{it} = x_{it}^\top \beta + u_{it}
$$
を一般化最小二乗法（generalized least squares, GLS）で推定します．$(u_{i1},\dots,u_{iT})^\top$の共分散行列$\Omega$は，対角成分が$\sigma^2 + \tau^2$，非対角成分が$\tau^2$です． よって，この$\Omega$の推定量$\hat{\Omega}$を用いると，$\beta$のGLS推定量は，
\hat{\beta}_{\mathrm{RE}} = \left( \sum_{i=1}^N X_i^\top \hat{\Omega}^{-1} X_i \right)^{-1} \sum_{i=1}^N X_i^\top \hat{\Omega}^{-1} y_i
となります．ただし，$X_i = (x_{i1},\dots,x_{iT})^\top, y_i = (y_{i1},\dots,y_{iT})^\top$です．$\hat{\beta}_{\mathrm{RE}}$は，変量効果推定量とも呼ばれます．この変量効果推定量は，線形混合モデルの記事で紹介したBLUPにおける$\hat{\beta}$と同じものです（$\hat{\Sigma} = \mathrm{diag}(\hat{\Omega},\dots,\hat{\Omega})$です）．
固定効果モデル
モデルの表現としては(2)式と同じですが，$E[b_i \mid X_i] \not= 0$のときは，$\hat{\beta}_\mathrm{RE}$ も $\hat{\beta}_\mathrm{POLS}$も，$\beta$の一致推定量になりません．具体的には，個別効果$b_i$と説明変数$X_i$が相関してしまうと，$E[b_i \mid X_i] \not= 0$になります．計量経済学の実証分析上，個別効果と説明変数は往々にして相関すると考えられます．
例えば，生産関数を推定するために，企業の生産$y_{it}$を説明する投入$x_{it}$をモデリングするとします．このとき個別効果$b_i$の解釈として経営能力が考えられますが，投入$x_{it}$は経営能力と相関すると考える方がもっともらしいです．
そこで，$E[b_i \mid X_i] \not= 0$のもとで，$\beta$の一致推定量を得る方法を考えます．
固定効果モデルの推定
今から3つの推定法を見ますが，得られる推定量は実はすべて同じです！
固定効果推定量（階差をとってGLS）
(2)のモデルの階差をとってモデルを変換することで，$b_i$を消し去ってしまおうというアイデアです．つまり，
$$
\tag{3}
\Delta y_{it} = y_{it} - y_{i,t-1} = \Delta x_{it}^\top \beta + \Delta \varepsilon_{it} \quad (i=1,\dots,N; \ t=2,\dots,T)
$$
というモデルに変換すると，$b_i$は消えています．これを行列表記します．$(T-1) \times T$行列として，
D =
\begin{bmatrix}
-1 & 1 & 0 & \cdots & 0 & 0 \\
0 & -1 & 1 & \cdots & 0 & 0 \\
\vdots & & \ddots & & \vdots \\
0 & 0 & 0 & \cdots & -1 & 1
\end{bmatrix}
を定義すると，
$$
D y_i = D X_i \beta + D\varepsilon_i \quad (i=1,\dots,N)
$$
が(3)のモデルの行列表記になります．このモデルの誤差項$D \varepsilon_i$には系列相関が生じているため，共分散行列を計算してGLSで推定します．$V(D \varepsilon_i \mid X_i) = \sigma^2 DD^\top$より，GLS推定量は，
\begin{align}
\hat{\beta}_\mathrm{FE} &= \left\{ \sum_{i=1}^N (DX_i)^\top (DD^\top)^{-1} DX_i \right\}^{-1} \sum_{i=1}^N (DX_i)^\top (DD^\top)^{-1} Dy_i \\
&= \left( \sum_{i=1}^N X_i^\top Q X_i \right)^{-1} \sum_{i=1}^N X_i^\top Q y_i
\end{align}
となります．ただし$Q = D^\top(DD^\top)^{-1}D$です．この$\hat{\beta}_\mathrm{FE}$は，固定効果推定量と呼ばれています．
Within-Group推定量
同じく，モデルを変換して$b_i$を消し去るというアイデアですが，階差をとるのではなく個体$i$の平均からの偏差を考えます．つまり，
y_{it} - \bar{y}_i = (x_{it} - \bar{x}_i)^\top \beta + (\varepsilon_{it} - \bar{\varepsilon}_i)
というモデルに変換します．この変換をWithin-Group（WG）変換といい，WG変換されたモデルのOLS推定量，つまり，
\hat{\beta}_\mathrm{WG} = \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)
をWG推定量といいます．実は，固定効果推定量の式中の$Q$は，$J_T$を全要素が1の$T \times T$行列とすると，$Q = I_T - T^{-1} J_T$であることが示せ，固定効果推定量とWG推定量は同じであることがわかります．
最小二乗ダミー変数推定量
$b_i$を個体$i$に対するダミー変数の回帰係数だと思って，これをOLSで推定して得られる$\beta$の推定量を，最小二乗ダミー変数（Least squares dummy variables, LSDV）推定量といいます．補助回帰のテクニックを用いると，LSDV推定量も，固定効果推定量とWG推定量と同じになることが示せます．
どのアプローチの場合も，結果として得られる推定量は同じなので気にしなくていいのですが，名称が入り乱れている上に，GLSは固定効果モデルでも変量効果モデルでも出てくるので，非常にややこしいです．OLSやGLSは，あくまで推定法であってモデルの名称ではないと既に述べましたが，どの（変換）モデルに対するOLSなのかGLSなのかをきちんと整理して理解できていることが大事だと思います．
補足
モデル？推定法？
固定効果と変量効果は，その名称から$b_i$を非確率変数として扱うか，確率変数として扱うかの違いのような印象を受けますが，実際は推定法の違いです．
固定効果のメリットとデメリット
メリット
固定効果の推定は，$b_i$が$X_i$と相関していても，$\beta$を一致推定できるというメリットがあります．そして多くの計量経済学の実証分析上は，$b_i$と$X_i$は相関していると考えられるため，固定効果モデルが採用されることが多いです．
デメリット
一方で固定効果の推定は，階差をとったりWG変換をすることで，サンプルサイズが$NT$から$N(T-1)$に減少しています．LSDV推定の観点から見ると，ダミー変数に対する回帰係数$N$個を余計に多く推定しなければなりません．つまり，サンプルサイズが$N$個減るか，推定しなければならないパラメータ数が$N$個増えるかで，実質同じことですが，推定の効率が落ちてしまいます．
また，$x_{it}$に時間不変な変数を入れることができない（階差やWG変換で時間不変な変数も消えてしまうから，LSDVだと多重共線性が生じるから）というデメリットもあります．
変量効果のメリットとデメリット
メリット
固定効果の推定のデメリットと裏返しの関係で，推定しなければならないパラメータは分散パラメータ$\tau^2$の1つが増えるのみです．また$x_{it}$に時間不変な変数も含めることができます．
デメリット
これも固定効果の推定のメリットと裏返しの関係で，$b_i$と$X_i$が相関している場合，$\beta$を一致推定できないというデメリットがあります．
定数項のあつかい
説明の便宜上，どのモデルにおいても$x_{it}$は定数項を含まないものとしてきましたが，通常はPooled OLSと変量効果モデルにおいては定数項を含めて推定します．一方，固定効果の推定の場合は，階差をとったりWG変換したりすると定数項は消えてしまうので（LSDV推定では定数項とダミー変数たちが多重共線性をおこすため），$x_{it}$に定数項は含めません（ただし一部のソフトウェアでは，$b_i$の推定値の平均が0になるような変換を施し，定数項の推定値を提示しているものもあります）．
つまり，定数項ありのPooled OLSをベースとした推定すべきパラメータ数の増加分は，変量効果モデルは1つ（変量効果の分散$\tau^2$），固定効果モデルは$(N-1)$個（ダミー変数の回帰係数$N$個から定数項$1$つを引いたもの）となります．
おわりに
株式会社Nospareには，統計学の様々な分野を専門とする研究者が所属しております．統計アドバイザリーやビジネスデータの分析につきましては株式会社Nospare までお問い合わせください．
