library(broom) # For converting models into tables
library(DT)
library(haven)
library(huxtable) # For side-by-side regression table
library(rdrobust) # For robust nonparametric RDD
library(rddensity) # For nonparametric RDD
library(stargazer)
library(tidyverse) # For ggplot
theme_set(theme_gray(base_size = 10, base_family = "HiraginoSans-W3")) # macOS用
theme_set(theme_gray(base_size = 10, base_family = "Meiryo")) # Windows用
instrumental variable
)instrumental varialble estimation
)「傾向スコアマッチング」における因果推論との違い
・「傾向スコアマッチング」では「無視可能な割り付けが成立している」ことが前提
= 処置を割り当てるかどうかのメカニズムが「混同的でない」(unconfounded
)
=
「処置を受けるかどうかの意思決定」と「処置を受けることでもたらされる結果」は無関係
= 「内生性 (endogeneity
) がない」
・しかし、私たちの現実の生活では「混同的な割当メカニズム」はしばしば存在する
・例1:夏休み中に行われる「補講授業」への参加は、前期の試験結果だけで決まるわけではなく、補講授業を受講することで後期の試験結果がどうなるかという予測に基づいて決まる
・例2:ある教員を採用するかどうかは、その教員の業績だけで決まるわけではなく、その教員を採用することで大学の業績(学生への人気)がどう変わるかの予想に基づいて決まる
・このような「混合的な割当メカニズム」がある場合
→ バイアスが発生するため「傾向スコアマッチング」などの手法では、因果効果を正しく推定できない
・そこで使われるのが、操作変数法
(instrumental varialble estimation
)
・操作変数法を使えば「無視可能な割り付けが成立していない時でも(=
内生性があっても)」因果効果を推定できる
instrumental variable estimate
)
を使う事でどのようなことが明らかになるのかいくつか例を挙げてみようHarada and Smith (2021)
は日本において、政治家が特定の地域へ補助金を配分することと、その地域の犯罪率の関係を分析補助金配分と犯罪率の分析結果
(Harada $ Angrist, 2021
)
・補助金受取額と犯罪率には負の関係がある
・日本の犯罪率が低いのは、補助金配分が抑制装置として機能していた可能性を示唆
Hangartner et al (2021)
はトルコからギリシャに難民が流入することで、ギリシャ人の難民への態度をどのように変化させたのか分析した難民を多く受け入れた地域と少ない移民しか受け入れない地域の人々の難民・移民への態度を単純比較した場合
難民や移民を容易に受け入れた地域 → そもそも難民に対して友好的な態度だった可能性あり
難民や移民に難色を示した地域 → そもそも難民に対して敵対的な態度だった可能性あり
単純比較では、難民流入の有無が、両地域における移民や難民に対する態度を引き起こす因果効果を正しく推定できない
難民の流入と難民に対する態度の分析結果
Hangartner et al (2021)
・難民を多く受け入れた島の住民は、移民受け入れを制限する政策を支持する傾向にある
Angrist (1990)
の論文単純な比較の問題点:
・「従軍経験のある人の平均賃金」と「従軍経験のない人の平均賃金」の単純比較してみる
・結果:「従軍経験のある人の平均賃金」<「従軍経験のない人の平均賃金」
→ しかし、これが「従軍経験自体によってもたれされた因果効果」とはいえない
・その理由:
→ 「従軍経験のある人」はもともと稼げる賃金が低いから、平均賃金が低いだけ
・「軍隊に入隊するという意思決定」が「従軍することによって得られる賃金と従軍しないことによって得られる賃金の比較」によってなされる「混合的
(confounded
) な割当メカニズム」だから
= 「処置を受けるかどうか(=
軍隊に入隊するかどうか)の意思決定」と「処置を受けることでもたらされる結果(平均賃金)」は無関係でないから
=「内生性 (endogeneity
) 」があるから
解決策:
・操作変数を用いた分析
→ 内生性を取り除いて因果効果を推定できる
The lottery was to determine the callup order of
nearly 2 million young men turning 19 this year on the basis of pairing
of 365 capsules containing the days of the year and an equal number of
blue capsules containing draft sequence numbers 1 through 365.
Source: https://photos.com/featured/curtis-tarr-spinning-draft-lottery-drum-bettmann.html
従軍兵士と兵士以外の生涯賃金の差の分析結果
(Angrist 1990
) ・従軍経験によって平均賃金は 15%
減少する
\(Y\) | 結果変数 | 賃金 | 観測可 |
\(X_1\) | 介入変数(共変量) | 従軍経験 | 観察可 |
\[\hat{Y} = \hat{β_0} + \hat{β_1}X_1 + \hat{β_2}X_2\]
\(Y\) | 結果変数 | 賃金 | 観測可 |
\(X_1\) | 介入変数(共変量) | 従軍経験 | 観察可 |
\(X_2\) | 共変量 | 交絡因子 | 観察不可 |
Z
(=くじ引き = 操作変数)を利用する\[Y_i = \beta_0 + \beta_1X_{1i}+ \beta_2X_{2i} + ε_i\]
二つの共変量 \(X_{1i}\) と \(X_{2i}\) は相関がある
\[cov[X_1, X_2] ≠ 0\]
誤差項 \(ε_i\) は次の
BLUE
の条件を満たす:
・「誤差項の期待値がゼロ」
・「誤差項の分散均一性」
・「誤差項の条件付き期待値がゼロ」=
誤差項と共変量は独立
しかし、共変量 \(X_{2i}\)
は観測できない
→ 実際に分析する際に使う重回帰式は次のようになる
\[Y_i = \beta_0 + \beta_1X_{1i}+ U_i\]
\[cov[X_1, U] ≠ 0\]
instrument exogeneity
)\[cov[Z, U]=0\]
instrument relevance
)\[cov[Z, X_1] ≠ 0\]
exclusion restriction
)\[\beta_1 = \frac{cov[Z, Y]}{cov[Z, X_1]}\]
操作変数推定量 \(\beta_1\) を求める式の証明 ・操作変数 \(Z_i\) と結果変数 \(Y_i\) の共分散を確認する
・\(Y=\beta_1X_1 + U\) なので
\[cov[Z, Y] = cov[Z, \beta_0 + \beta_1X_1
+ U]\]
・確率変数 \(X\) と \(Y\)
が独立でない場合の共分散は次のように変形できる
\[cov[X_1, Y_1 + Y_2]=cov[X_1, Y_1] +
cov[X_1, Y_2]\]
・従って、\(cov[Z, Y]\)
は次のように表現できる
\[cov[Z, Y] = cov[Z, \beta_0 + \beta_1X_1 + U]\\ cov[Z, \beta_0] + cov[Z, \beta_1X_1] + cov[Z, U]\]
・\(\beta_0\) は定数なので \(cov[Z, \beta_0] = 0\)
・\(k\) が定数なら \(cov[kX, Y] = kcov[X, Y]\)
→ \(\beta_1\) は定数なので \(cov[Z, \beta_1X_1] = \beta_1cov[Z,
X_1]\)
・従って、\(cov[Z, Y]\) は次のように表現できる
\[cov[Z, Y] = cov[Z, \beta_0] + cov[Z, \beta_1X_1] + cov[Z, U]\\ = \beta_1cov[Z, X_1] + cov[Z, U]\]
・「1. 操作変数の外生性」が満たされるなら \(cov[Z, U]=0\) なので
\[cov[Z, Y] = \beta_1cov[Z, X_1] \]
・「2. 操作変数の関連性」が満たされれば、 \(cov[Z, X_1] ≠ 0\)
→ 両辺を \(cov[Z, X_1] ≠
0\)で割り、操作変数推定量 (\(\beta_1\)) を求めることができる
\[\beta_1 = \frac{cov[Z, Y]}{cov[Z, X_1]}\]
\(Y_i\) | 結果変数 | 観測可 |
\(x_{1i}\) | 介入変数(共変量) | 観察可 |
\(x_{2i}\) | 共変量(=交絡因子) | 観察不可 |
\(z_i\) | 操作変数 | 観察可 |
rnorm()
関数を使って、それぞれ個別に生成するset.seed(1)
<- 1000 # N = 1000 と設定
n1 <- rnorm(n1) # N = 1000 の x2 を生成する
x2 <- runif(n1) # N = 1000 の z を生成する z
<- rnorm(n1) e1
BLUE
)
最初の仮定である「誤差項の期待値(=平均値)がゼロ」を確認mean(e1)
[1] -0.01021977
<- 1 + 2 * z - 3 * x2 + e1 x1
<- rnorm(n1) e2
<- x2 + e2 u1
<- 3 + 1.5 * x1 + 1.5 * x2 + e2 y1
\[y_1 = 3 + 1.5x_1 + 1.5x_2 + ε_2\]
plot(x1, y1)
plot(x2, y1)
plot(x2, z)
1. 操作変数の外生性
(instrument exogeneity
) \[cov[Z, U]=0\]
cov(z, u1)
[1] 0.004487261
2. 操作変数の関連性
(instrument relevance
) \[cov[Z, X_1] ≠ 0\]
cov(z, x1)
[1] 0.144279
<- lm(x1 ~ z)
relevance tidy(relevance)
term | estimate | std.error | statistic | p.value |
---|---|---|---|---|
(Intercept) | 1.17 | 0.202 | 5.8 | 8.78e-09 |
z | 1.7 | 0.354 | 4.8 | 1.8e-06 |
p-value = 1.8e-06
で統計的に有意3. 除外制約
(exclusion restriction
) ・操作変数 \(Z_i\) が次の式に含まれていないこと
\[Y_i = \beta_0 + \beta_1X_{1i}+ U_i\\cov[X_1, U] ≠ 0\]
・操作変数 \(Z_i\) と誤差項 \(U_i\) と無相関であること
・\(Z_i\) から \(X_{1i}\) に矢印が出ている
・\(Z_i\) から \(Y_i\) には矢印が出ていない
→ \(Z_i\) は \(Y_i\) に間接的に影響を与えている
→ \(X_i\) の値が決まれば、 \(Z_i\) があってもなくても、\(X_i\) から \(Y_i\) への影響は変わらない
cov(x1, u1)
[1] -3.157631
\(Y_i\) | 結果変数 | 観測可 |
\(x_{1i}\) | 処置変数(共変量) | 観察可 |
\(x_{2i}\) | 共変量(=交絡因子) | 観察不可 |
\(z_i\) | 操作変数 | 観察可 |
<- lm(y1 ~ x1 + x2)
m_1 tidy(m_1)
term | estimate | std.error | statistic | p.value |
---|---|---|---|---|
(Intercept) | 3.07 | 0.0631 | 48.7 | 1.47e-265 |
x1 | 1.49 | 0.0273 | 54.8 | 8.16e-303 |
x2 | 1.47 | 0.0868 | 17 | 7.8e-57 |
<- lm(y1 ~ x1)
m_2 tidy(m_2)
term | estimate | std.error | statistic | p.value |
---|---|---|---|---|
(Intercept) | 3.91 | 0.0439 | 89.1 | 0 |
x1 | 1.06 | 0.0114 | 93.4 | 0 |
<- cov(z, y1)/cov(z, x1)
m_3 tidy(m_3)
x |
---|
1.57 |
モデルの種類 | 推定値 |
真の推定値 | 1.49 |
単回帰モデル | 1.06 |
操作変数モデル | 1.57 |
シミュレーションの結果まとめ
・真の推定値は
1.49
・単回帰モデルの推定値 (1.06) は真の推定値 (1.49)
を過小評価している
・操作変数モデルの推定値
(1.57) は真の推定値 (1.49) に近い値を推定している
→ 未観測の共分散 \(x_2\)
の交絡を除去している
two-stage least squares: 2SLS
) を使っても分析できる\[X_i = \alpha_0 + \alpha_1Z_i + e_i\]
\[Y_i = \beta_0 + \beta_1\hat{X_i} + u_i\]
<- lm(x1 ~ z) # 第 1 段階の最小二乗法
step_1 <- predict(step_1) # x1 の予測値 x1hat を抽出
x1hat
<- lm(y1 ~ x1hat) # 第 2 段階の最小二乗法
step_2 tidy(step_2) # 結果を表示
term | estimate | std.error | statistic | p.value |
---|---|---|---|---|
(Intercept) | 2.9 | 0.477 | 6.08 | 1.68e-09 |
x1hat | 1.57 | 0.231 | 6.78 | 2.03e-11 |
2 段階最小二乗法を使った推定値
・真の推定値は
1.49
・操作変数モデルの推定値
(1.57)
・2 段階最小二乗法を使った推定値
(1.57)の値は操作変数モデルの推定値と同じ
→ 未観測の共分散 \(x_2\)
の交絡を除去している