操作変数法

**操作変数法**（そうさへんすうほう、）とは、統計学、計量経済学、疫学、また関連分野において、統制された実験が出来ない時、もしくは処置がランダムに割り当てられない時に、因果関係を推定するための方法である。直感的に言えば、操作変数は説明変数と被説明変数の間の相関が二変数間の因果関係をもっともらしく反映していない時に用いられる。妥当な操作変数は説明変数に影響を与えるが被説明変数に独立的な影響を持たず、研究者が被説明変数に対する説明変数の因果効果を明らかにすることを可能とする。

操作変数法は説明変数（共変数）が回帰モデルにおける誤差項と相関している時にすることを可能とする。このような相関は、被説明変数の変化が共変数の少なくとも一つの値を変化させる時（"逆"の因果）、説明変数と被説明変数の双方に影響を与える除外変数が存在する時、共変数に測定誤差がある時（error-in-variables models）に起こるだろう。回帰の文脈において一つないしは複数の問題を持つ説明変数は時折、内生性として言及される。この状況下では、最小二乗法はバイアスを持ち一致性を持たない推定量を生み出す。しかし、もし*操作変数*が利用可能ならば、一致推定量を得ることができる。操作変数とはそれ自身は説明すべき方程式には依存していないが、内生的な説明変数とほかの共変数の値による条件の下で相関している変数のことである。線形モデルにおいては操作変数法を用いるために二つの必要な仮定がある。
* 操作変数は他の共変数で条件付けた時に、内生的な説明変数と相関しなくてはならない。もしこの相関が統計的に有意なほど高ければ、その操作変数は**強い第一段階**（）を持つと言う。相関が弱いとパラメータの推定値と標準誤差について間違った推論を導きかねない。
* 操作変数は説明方程式の誤差項と他の共変数で条件付けた時に相関してはならない。言い換えると、操作変数は元の予測変数と同じ問題に直面することがない。もしこの条件が満たされているならば、その操作変数は**除外制約**（）を満たすと言う。

## 導入
操作変数法の概念はフィリップ・ライト（）と共著者で息子のシューアル・ライトにより、1928年に出版された著書*The Tariff on Animal and Vegetable Oils*においての文脈で導出された。1945年、Olav Reiersøl は彼の学位論文において、errors-in-variables modelsの文脈で同じ手法を用い、その手法に名前を与えた。

操作変数法の背後にあるアイデアは広いモデルのクラスに拡張できるが、操作変数法についての非常に一般的な文脈は線形回帰にある。伝統的に、操作変数は Z で定義され、操作変数と相関を持つ独立変数は X 、操作変数と無相関な誤差項は U として定義され、以下のような方程式を満たす。

: Y = \beta X + U \, 

ここで X は通常、1のみからなる列と他の共変数からなる追加的な列を持つ行列である。この場合において操作変数が解くことのできる問題について考えよう。すると、操作変数法がいかにして問題を解くかを示すことができる。最小二乗法(OLS)が、cov(X,U) = 0 の下で  \beta  について問題を解くことを思い出そう（これは簡単である。誤差の二乗和を最小化する時、min (Y-\beta  X)^2 、一階条件はまさしく  X' (Y-\beta X) = X' U = 0  である）。もし、上でリストアップした理由の一つのために、本当のモデルでは cov(X,U) \neq 0 であるとしたら、例えばもし X と Y の両方に別々に影響を与える除外変数が存在するならば、OLSの手続きは、Y に対する X の因果的な効果を生み出さないだろう。OLSはただ単純に X と相関しないように結果的になる誤差を生み出すパラメータを取り出すであろう。

一変数の場合を用いるとより明確になる。一変数と定数についての回帰を考えているとしよう（ひょっとしたら他の共変数は必要ないかもしれない、またひょっとしたらすでに他の関連する共変数をしているかもしれない）。

:y=\alpha + \beta x + u

この場合、興味のある説明変数に対する係数は  \beta = \frac{cov(x,y)}{var(x)}  として与えられる。y について代入すると以下のようになる。

:  \beta = \frac{cov(x,y)}{var(x)} = \frac{cov(x,\alpha + \beta x + u)}{var(x)} = \beta + \frac{cov(x,u)}{var(x)} \rightarrow \beta 

もし想定しているモデル上において cov(x,u) \neq 0 ならば、OLSは興味のある因果効果を反映*していない*係数を推定する。操作変数法は、x が u と無相関か否かというより、他の変数 z が u と無相関か否かに基づいてパラメータ \vec{\beta} を識別するので問題を解決することが可能になる。もし理論上 z が x と関係し（第一段階）、u と無相関（除外制約）ならば、操作変数法は最小二乗法が失敗した興味のある因果パラメータを識別するだろう。線形の場合には操作変数推定量を使い導出する複数の特定の方法が存在するので、さらなる議論は推定の節で行う。

もちろん、操作変数法はより広い非線形モデルのクラスにも適用されてきた。操作変数の一般的な定義は、反事実的かつグラフィカルな形式論を用いることで、ジューディア・パールによって与えられた。グラフィカルな操作変数の定義は *Z* が次の条件を満たすことで与えられる。

: (Z \perp\!\!\!\perp Y)_{G_{\overline{X}}} \qquad(Z \not\!\!{\perp\!\!\!\perp} X)_G 

ここで \perp\!\!\!\perp はベイジアン・ネットワークにおける であり、G_{\overline{X}} は X に入る矢印がすべてカットオフされるようなベイジアン・ネットワークにおけるグラフである。

操作変数の反事実的（）な定義は操作変数 *Z* が以下を満たすことである。
: (Z \perp\!\!\!\perp Y_x)\qquad (Z \not\!\!{\perp\!\!\!\perp} X)
ここで Y_x は X が x であった時の Y の取りうる値であり、 \perp\!\!\!\perp は独立を表している。

もし追加的な共変数 W があるのならば、W で条件付けた下での操作変数 Z として定義が変更される。

パールの定義のエッセンスは
# 興味のある方程式は"構造的"なものであり、単なる"回帰"ではない。
# 誤差項 U は X が定数である時に Y に影響を与えるすべての外生的要因を表している。
# 操作変数 Z は U と独立でなければならない。
# 操作変数 Z は X が定数ならば、Y に影響を与えてはならない。（除外制約）
# 操作変数 Z は X と独立ではない。
ということである。

これらの条件は方程式の特定の関数形に依存しておらず、ゆえに非線形方程式、つまり誤差項 U が非加法的である場合にも適用できる。これらの条件はまた複数の方程式からなるシステムにも適用でき、そこでは X（と他の要因）がいくつかの中間的な変数を通して Y に影響を与える。操作変数は X の原因である必要はない。そのような原因の代理変数は、条件1から5を満たすならば、操作変数としてまた使えるだろう。また、除外制約（条件4）は条件2と3から導けるので省略できる。

## 例
カジュアルな言い方では、ある変数 X が他の変数 Y に与える因果効果を推定しようとする時、操作変数は X における効果のみを通して Y に影響を与える変数 Z のことである。例えば、研究者は喫煙の一般的な健康における因果効果を推定したいとしよう。健康と喫煙の相関は、他の変数が健康と喫煙の両方に影響を与えた、もしくは健康状態が喫煙に影響を与えたと考えることもできるので、喫煙が健康を悪化させる原因であるということは意味しない。一般の母集団において喫煙状態を制御した実験を行うのはとても難しく費用がかかる。研究者は、因果分析における喫煙の操作変数としてタバコ製品の税率の時系列を用いて観測データから健康における喫煙の因果効果を推定しようとするだろう。研究者はタバコ製品についての税率は喫煙に与える効果のみを通して健康に影響を与えると仮定するので、タバコ製品についての税率は操作変数としては合理的な選択である。もし研究者がタバコ税と健康状態が相関しているのを発見できたならば、それは喫煙が健康状態の変化の原因であるという証拠と見なしうる。

ヨシュア・アングリストとアラン・クルーガーは操作変数法の使用と歴史についてのサーベイを行っている。

## 適切な操作変数の選択
U は観測できないので、Z が U から独立であるという仮定が満たされるかどうかはデータからは導かれず、かわりにモデルの構造、つまりデータの生成過程を決定しなくてはならない。因果のグラフはこの構造を表現していて、上で与えられているグラフの定義は変数 Z が共変数 W が与えられた下で操作変数として機能するかどうかを手早く決めるために用いることが出来る。この点を見るために、以下の例を考えよう。

ランダムに学生が寮に割り当てられるような大学において、大学のチュータープログラム（Tutoring Program）がGPAに与える影響を推定したいとしよう。チュータープログラムへの出席とGPAの間の関係は数多くの要因によって混同される。チュータープログラムに出席した学生は自分の成績に注意を払うかもしれないし、努力をするかもしれない（この混同は図.1-3における右側のチュータープログラムとGPAの間の双方向に伸びる弧を通して描写される）。学生が寮にランダムに割り当てられたとして、チュータープログラムが行われる場所からの学生の寮の近さは操作変数の自然な候補になる。

しかしながら、チュータープログラムが大学の図書館で行われたらどうなるだろうか。チュータープログラムが行われる場所と寮の近さは学生がより図書館で時間を費やそうとする原因になりうるだろうし、今度はそれがGPAを改善するだろう（図1を参照）。図2において描写されている因果グラフを用いると、チュータープログラムの行われる場所の近さは、グラフ G_{\overline{X}} において「場所の近さ（Proximity）  \rightarrow  図書館の使用時間（Library Hours）  \rightarrow  GPA」という経路を通してGPAと繋がっているため、操作変数としてはふさわしくない。しかしながら、共変数として図書館の使用時間を加えてコントロールすれば、G_{\overline{X}} において図書館の使用時間が与えられた下で場所の近さはGPAから分離されるため、場所の近さは操作変数になる。

今、図3のように学生"そのものが持つ能力"（Ability）が学生のGPAと同じくらい学生の図書館にいる時間に影響を与えるとしよう。因果グラフを用いれば、図書館の使用時間はコライダーとなり、図書館で条件付けることは「場所の近さ（Proximity） \rightarrow 図書館の使用時間（Library Hours） \leftrightarrow GPA」という経路を作るだろう。結果として、場所の近さは操作変数として用いることは出来ない。

最後に、図4で示されているように、図書館で勉強しない学生は単に他の場所で勉強するので、図書館の使用時間はGPAに何の影響も与えないとしよう。この場合、図書館の使用時間を制御することは場所の近さからGPAまでの信用できない経路を開く。しかしながら、図書館の使用時間をコントロールせず、共変数から除外すると、場所の近さは操作変数として用いることが出来る。

## 推定
今またここで操作変数法の詳細を考えよう。以下のような形でデータが生成されるとする。
:  y_i = X_i' \beta  + e_i, 
ここで
* i は観測値の添え字、
* y_i は被説明変数、
* X_i は説明変数と定数のベクトル、
* e_i は X_i とは異なる y_i に影響を与えるすべての要因を表す観測できない誤差項、
* \beta は観測できないスカラー値のパラメータ、
* 上付きの添え字 ' は行列ないしはベクトルの転置、
とする。

パラメータ \beta は X_i の各要素が一単位動き、他の y_i に変動を与えるすべての要因が一定である時に y_i が受ける因果効果を表している。計量経済学的な目的は \beta を推定することである。単純化のために、e は互いに無相関で、同じ分散である分布から生成されるものとする。つまり誤差項は自己相関がなく分散均一である。

また同じ形の回帰モデルも導出できるとする。観測値のランダムなサンプルのサイズを T とすると、最小二乗法による推定量は以下のようになる。
:  \widehat{\beta}_\mathrm{OLS} = (X'X)^{-1} X'y = (X'X)^{-1} X'(X \beta + e) = \beta + (X'X)^{-1} X'e
ここで X はそれぞれの X_i' を並べた T 行の行列、y と e は長さ T の列ベクトルを表している。この方程式は導入節における  cov(X,y)  についての方程式と類似している（ここでの方程式は行列バージョンである）。X と e が無相関であるとき、ある正則条件の下で第2項を X で条件付けた期待値は0となり、さらに極限において0に収束する。よってこの推定量は不偏かつ一致推定量である。X と e に含まれる測定されない因果変数が相関すると、しかしながら、最小二乗法による推定量は一般的に \beta についてバイアスを持ち、一致性もない。この場合、X の値が与えられた場合の y の値を予測するための推定量としては妥当であるが、X の y に対する因果効果はこの推定量では分からない。

パラメータ  \beta  を正しく推定するために、それぞれの内生的な X と強く相関するが、y とは相関しない（言い換えれば、e とは相関しない）変数 Z を導入する。簡単化のために、X は 定数と内生変数の列からなる T 行2列の行列であるとし、Z は 定数と操作変数の列からなる T 行2列の行列であるとする。しかしながら、この方法は X が定数と、例えば、5つの内生変数からなる行列であり、Z が定数と5つの操作変数からなる場合といった時にも拡張できる。以下の議論においては X は T 行 K 列の行列であり、K は未定のままであると仮定する。X と Z が共に T 行 K 列の行列である時の推定量はと言われる。

それぞれの内生的要素 x_i と操作変数の間の関係が以下のように与えられると仮定する。

:  x_i = Z_i \gamma + v_i, 

最も一般的な操作変数による特定化は以下の推定量を用いる。

: \widehat{\beta}_\mathrm{IV} = (Z'X)^{-1} Z'y 

この特定化は、真のモデルにおいて  Z' e = 0  が満たされる限り、サンプルサイズが大きくなれば真のパラメータへと近づいていく。

: \widehat{\beta}_\mathrm{IV} = (Z'X)^{-1} Z'y = (Z'X)^{-1} Z'X \beta + (Z'X)^{-1} Z'e  \rightarrow \beta    

データを生成する過程において  Z' e = 0  が満たされる限り、操作変数推定量の適切な使用により、パラメータ \beta が識別される。操作変数法は  Z'e = 0  を満たす一意なパラメータについて解くので、これは機能し、そしてゆえにサンプルサイズが大きくなるにつれ真のパラメータに近づいていく。

今、拡張を行う。興味のある方程式における共変数の数より操作変数の数の方が大きいとする。つまり Z は T 行 M 列行列で M>K であるとする。これはしばしば**過剰識別**のケースと呼ばれる。この場合、一般化モーメント法(GMM)を用いることができる。GMM推定量は以下のようになる。

: \widehat{\beta}_\mathrm{GMM} = (X' P_Z X)^{-1}X' P_Z y,

ここで P_Z は射影行列であり、P_Z=Z(Z' Z)^{-1}Z' を満たす。

この表現は、操作変数の数と興味のある方程式における共変数の数が一致する時に、最初の表現にまとめることができる。過剰識別の操作変数法はそれゆえに、適切に識別された場合の操作変数法の一般化の一つである。

\widehat{\beta}_\mathrm{GMM} を以下のように展開する。
: \widehat{\beta}_\mathrm{GMM} = (X' Z(Z' Z)^{-1}Z' X)^{-1}X' Z(Z' Z)^{-1}Z' y

適切に識別されている時、操作変数の数と共変数の数は同じである。つまり X の次元と Z の次元は同じである。ゆえに、X' Z, Z' Z と Z' X は同じ次元の正方行列である。任意の n 行 n 列の行列 A と B について (AB)^{-1} = B^{-1}A^{-1} であるという事実（逆行列を参照のこと。）を用いて、（\widehat{\beta}_\mathrm{GMM} の中の）逆行列を展開すると、

:
\begin{align}
\widehat{\beta}_\mathrm{GMM} &= (Z' X)^{-1}(Z' Z)(X' Z)^{-1}X' Z(Z' Z)^{-1}Z' y\\
&= (Z' X)^{-1}(Z' Z)(Z' Z)^{-1}Z' y\\
&= (Z' X)^{-1}Z' y \\
&= \widehat{\beta}_\mathrm{IV}
\end{align}


となり、\widehat{\beta}_\mathrm{GMM} と \widehat{\beta}_\mathrm{IV} が（適切に識別されている時に）一致することが分かる。参照文献として Davidson and Mackinnnon (1993) を挙げておく。

ここで、m の場合について同値な推定量が存在する。パラメータは線形方程式のシステムの解であるので、方程式  Z'v = 0  を用いた過小識別モデルは一意解を持たない。

## 二段階最小二乗法としての解釈
操作変数法による推定量を計算するために使用できる計算方法の一つが二段階最小二乗法（）である。第一段階において、興味のある方程式において内生性を持つすべての説明変数を、方程式における外生的な共変数と除外されている操作変数を含んだ全てのモデルにおける外生的な変数で回帰する。そしてこの回帰における予測値が得られる。

**第一段階：** X のすべての列を Z で回帰する。（ X = Z \delta + \text{errors} ）

: \widehat{\delta}=(Z' Z)^{-1}Z' X, \,

そして予測値を保存する。

: \widehat{X}= Z\widehat{\delta} = Z(Z' Z)^{-1}Z' X = P_Z X.\, 

第二段階では、内生変数をすべて第一段階での予測値に置き換える以外は、興味のある回帰は通常通り推定される。

**第二段階：** Y を第一段階における予測値で回帰する。

:  Y = \widehat X \beta + \mathrm{noise}.\,

ここから以下が得られる。

:  \beta_{2SLS} = \left(X' P_Z X\right)^{-1} X' P_ZY

通常の最小二乗法推定量は  (\widehat X' \widehat X)^{-1}\widehat X' Y である。 \widehat X = P_Z X と置き換え、P_Z  は対称なべき等行列とする。つまり、P_Z' P_Z=P_Z P_Z = P_Z である。すると
: 
\beta_{2SLS} = (\widehat X' \widehat X)^{-1}\widehat X' Y 
= \left(X' P_Z' P_Z X\right)^{-1} X' P_Z' Y
=\left(X' P_Z X\right)^{-1} X' P_ZY
となる。

\beta の結果的な推定量は上で示した表現と数値的に一致する。\beta の共分散行列を正しく計算するためには、第二段階のモデルにおける残差二乗和について微修正が必要である。

## ノンパラメトリック分析
構造方程式の形式が未知の時、操作変数 Z は以下の方程式を通して依然として定義可能である。

:x = g(z,u) \, 
:y = f(x,u) \, 

ここで f と g は2つの任意な関数で Z は U からは独立である。線形モデルと異なり、Z, X と Y の測定は X の Y に対する平均因果効果（）の識別を可能としない。

:
\text{ACE} = \Pr(y\mid \text{do}(x)) = \operatorname{E}_u[f(x,u)].


Balke and Pearl (1997) はACEの狭い境界を導出し、その境界はACEの符号と大きさについて価値ある情報を提供しうることを示した。

線形分析では、Z が (X,Y) の操作変数であるという仮定が正しいかどうか調べる検定はない。ただし X が離散変数ならばそうではない。Pearl (2000) は全ての f と g について、Z が上に上げた二つの方程式を満たすときはいつでも、以下の制約（"操作不等式"、）を満たさなくてはいけないことを示した。
:\max_x \sum_y [\max_z \Pr(y,x\mid z)]\leq 1.

## 操作変数推定量の解釈において
上での例は興味のある因果効果は観測値の間で不変であることを仮定している、つまり、\beta が定数であるということを仮定している。一般的に、主体が異なれば、"処置" x の変化に対する反応も異なる。この可能性を考慮に入れると、x の y における変化の母集団における平均的な効果は与えられた部分母集団における平均的な効果とは異なるだろう。例えば、職業訓練プログラムの平均効果は、訓練を実際に受講した人々からなるグループと受講しなかった人々からなるグループとで実質的に異なるだろう。これらの理由により、操作変数法は、行動的な反応に対する暗黙的な仮定、もしくはより一般的に処置への反応と処置を受けるかどうかの傾向の間の相関についての仮定を課している。

標準的な操作変数推定量は平均処置効果（）というよりは局所平均処置効果（）を取り出すことができる。Imbens and Angrist (1994) は線形な操作変数推定量は弱い条件の下で局所平均処置効果の加重平均と見なせることを示した。ここでその加重は操作変数の変化に対する内生的な説明変数の弾力性に依存している。端的に言えば、観測された操作変数の変化に反応した部分母集団についてのみ変数の効果は表れ、操作変数の変化に最も反応した部分母集団が操作変数の大きさや程度に最も大きな影響を持つだろうということを意味している。

例えば、研究者がランドグラント大学の存在を所得の回帰における大学教育の操作変数として使う時、研究者は大学が存在していれば学位を取るが、存在していなければ学位を取らなかった部分母集団における大学の所得に対する影響を識別する。この実証研究では、追加的な仮定なしには、大学があろうとなかろうと学位を取る、もしくは取らない人々における大学の効果について研究者は何も言えない。

## 潜在的な問題
操作変数が興味のある方程式の誤差項と相関していると、操作変数推定量は一般的には一致推定量ではない。Bound, Jaeger, and Baker (1995) が記すように、"弱"操作変数の選択によりまた別の問題が起こる。弱操作変数とは、第一段階の方程式において内生的な説明変数の予測値としては不十分であるような操作変数のことである。この場合、弱操作変数による内生変数の予測は不十分で、予測値はほとんど変動しない。結果として第二段階の方程式において内生変数を置き換えるために弱操作変数を使ったとき、最終的な値を予測することが難しくなる。

上で議論した喫煙と健康の例における文脈では、喫煙状態がタバコ税の変動にほとんど反応しないのであれば、タバコ税は喫煙についての弱操作変数となる。もし高い税率によって人々がタバコをやめる（もしくは喫煙を始めない）ことがないのであれば、タバコ税の変動は健康に対する喫煙の影響について何も情報を持たない。もしタバコ税が喫煙に与える影響というより他の経路を通して健康に影響を与えるのであれば、タバコ税は操作変数としては不適格であるし、操作変数法は間違った結果を生みかねない。例えば、相対的に健康志向な人々がいる時と場所においては、高い税率のタバコ税の実施と、たとえタバコ税が低率のままであったとしても健康であろうする行動の両方がなされるだろう。そして、仮にもし喫煙が健康に何の影響も与えないとしても、健康とタバコ税の間に相関が見られるだろう。この場合、タバコ税と健康の間に観測された相関から健康に対するタバコの因果効果を推測するのは過ちであろう。

## 標本性質と仮説検定
共変数が外生的である時、最小二乗推定量の小標本における性質は、X で条件づけた推定量のモーメントを計算するという標準的な方法で得ることができる。共変数のいくつかが内生的で操作変数法を用いている時、推定量のモーメントによる単純な表現は得られなくなる。一般的に操作変数推定量は漸近的に望ましい性質のみを持ち、有限標本における性質は明らかではなく、推量も推定量の標本分布に漸近的に近似することをベースにしている。操作変数が興味のある方程式の誤差項と無相関で、かつ弱操作変数でないとしても、操作変数推定量の有限標本における性質は大したものは得られない。例えば、適切に識別されたモデルからはモーメントを持たない有限標本推定量が得られ、ゆえにこの推定量のバイアスについて言えることは何もないし、検定統計量の名目値も実質的に歪み、推定量は一般的にパラメータの真の値から離れたものになる。

## 操作変数の強さの検定と過剰識別制約
内生的な共変数と操作変数は両方とも観測可能なので、操作変数の強さは直接的に評価できる。内生的な説明変数が一つの場合によく使われる経験的な方法の一つが以下のようなものである。第一段階の回帰において除外されている操作変数が無関係であるという帰無仮説の下でのF検定統計量が10以上かどうかで判別する。

操作変数が興味のある方程式の誤差項と相関していないという仮説は、適切に識別されたモデルでは検定不可能である。もしモデルが過剰識別ならば、この仮説を検定するために用いることのできる利用可能な情報がある。これらの*過剰識別制約*の最も一般的な検定で、と呼ばれるものは、もし操作変数が本当に外生的ならば、回帰残差は外生変数とは無相関だろうという観察に基づいている。サーガン-ハンセン検定統計量は回帰残差を外生的な変数に対し最小二乗法で回帰した際の TR^2（決定係数に観測値の数を掛けたもの）として計算できる。この検定統計量（の分布）は誤差項は操作変数と無相関であるという帰無仮説の下で自由度 m-k のカイ二乗分布に漸近的に収束する。

## 脚注・参照文献

## 発展的参考文献
* 
* 
* 

## 関連項目

*

* 統計的因果推論

## 外部リンク
* http://www.urban.org/toolkit/data-methods/instrumental.cfm Layman による操作変数法の説明
* http://emlab.berkeley.edu/users/mcfadden/e240b_f01/ch4.pdf ダニエル・マクファデンの教科書の第4章
*  Mark Thoma による YouTube ビデオ
*  Mark Thoma による YouTube ビデオ
* ivpack: Instrumental Variable Estimation - CRAN