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用

I. 操作変数法とは

操作変数 (instrumental variable)

  • 人間が操作できないような偶然性の強い状況で生み出される変数

操作変数法 (instrumental varialble estimation)

  • 偶然起こったこと(=操作変数)のために変数間に起こる「連鎖反応」(処置変数 → 結果変数)に着目して 2 段階で分析し、処置群と統制群の比較を行う分析方法
    → 操作変数法は自然実験の一部だと見なすことができる

  • 無視可能な割り付けが成立していない時に、交絡を調整する方法
  • 「内生性」と「外生性」に関する手法として計量経済学では一般的な手法

「傾向スコアマッチング」における因果推論との違い ・「傾向スコアマッチング」では「無視可能な割り付けが成立している」ことが前提
= 処置を割り当てるかどうかのメカニズムが「混同的でない」(unconfounded)
= 「処置を受けるかどうかの意思決定」と「処置を受けることでもたらされる結果」は無関係
= 「内生性 (endogeneity) がない」
・しかし、私たちの現実の生活では「混同的な割当メカニズム」はしばしば存在する

例1:夏休み中に行われる「補講授業」への参加は、前期の試験結果だけで決まるわけではなく、補講授業を受講することで後期の試験結果がどうなるかという予測に基づいて決まる

例2:ある教員を採用するかどうかは、その教員の業績だけで決まるわけではなく、その教員を採用することで大学の業績(学生への人気)がどう変わるかの予想に基づいて決まる

・このような「混合的な割当メカニズム」がある場合
→ バイアスが発生するため「傾向スコアマッチング」などの手法では、因果効果を正しく推定できない
・そこで使われるのが、操作変数法 (instrumental varialble estimation)
・操作変数法を使えば「無視可能な割り付けが成立していない時でも(= 内生性があっても)」因果効果を推定できる

II. 操作変数法を使って分かること

  • ここでは、操作変数法 (instrumental variable estimate) を使う事でどのようなことが明らかになるのかいくつか例を挙げてみよう

1. 補助金配分と犯罪率の関係

  • Harada and Smith (2021) は日本において、政治家が特定の地域へ補助金を配分することと、その地域の犯罪率の関係を分析

背景:

  • 日本では都市と地方では一票の格差が大きく、一票がより「軽い」都市の選挙区と比較すると、一票がより「重い」地方の選挙区により多くの補助金が配分されてきた

理論:

  • 補助金は貧しい住民の雇用や収入を保証する
  • 雇用や収入を保証された地域では犯罪率が低下する  

単純比較の問題点:

  • 補助金受取額が少ない地域と多い地域の犯罪率の単純比較
  • 補助金受取額が少ない地域 → そもそも経済状態が良く犯罪率が低い可能性あり
  • 補助金受取額が多い地域 → そもそも経済状態が悪く犯罪率が高い可能性あり
  • 両地域における犯罪率の単純比較では、補助金受取の因果効果を正しく推定できない

解決策:

  • 1990年代における選挙制度改革における「一票の格差是正」を操作変数として分析  

  • 1994年に衆院選の選挙制度は中選挙区制から小選挙区比例代表並立制に変更
    → 一票の格差が是正され、1996年の総選挙以降、1 票の「重み」が選挙区間でかなり平等になった
    → 1 票の「重み」是正は機械的に実施された
    → 特定の政治家が特定の地域を優遇した可能性は低い

補助金配分と犯罪率の分析結果 (Harada $ Angrist, 2021) ・補助金受取額と犯罪率には負の関係がある
・日本の犯罪率が低いのは、補助金配分が抑制装置として機能していた可能性を示唆

2. 難民の流入と難民に対する態度の関係  

  • Hangartner et al (2021) はトルコからギリシャに難民が流入することで、ギリシャ人の難民への態度をどのように変化させたのか分析した

背景:

  • 世界的に移民や難民が増えている
  • ヨーロッパでは移民や難民の受け入れに対する反発が起こっている

理論:

  • 移民や難民の受け入れた住民は、移民や難民に対して反発するようになる

単純比較の問題点:

  • 難民を多く受け入れた地域と少ない移民しか受け入れない地域の人々の難民・移民への態度を単純比較した場合

  • 難民や移民を容易に受け入れた地域 → そもそも難民に対して友好的な態度だった可能性あり

  • 難民や移民に難色を示した地域 → そもそも難民に対して敵対的な態度だった可能性あり

  • 単純比較では、難民流入の有無が、両地域における移民や難民に対する態度を引き起こす因果効果を正しく推定できない

解決策:

  • 2015年から2016年にかけてのヨーロッパでの難民危機に注目
  • 調査対象:ギリシャ国内のエーゲ海の島の住人たち
  • これらの島々はトルコから海を越えてギリシャにやってくる難民たちの到着地
  • とること距離の近い島々ほど多くの難民がやってきたことに注目
  • 「トルコと各島の距離」を操作変数として使用

  • 1994年に衆院選の選挙制度は中選挙区制から小選挙区比例代表並立制に変更
    → 一票の格差が是正され、1996年の総選挙以降、1 票の「重み」が選挙区間でかなり平等になった
    → 1 票の「重み」是正は機械的に実施された
    → 特定の政治家が特定の地域を優遇した可能性は低い

難民の流入と難民に対する態度の分析結果 Hangartner et al (2021) ・難民を多く受け入れた島の住民は、移民受け入れを制限する政策を支持する傾向にある

3. 従軍兵士と兵士以外の生涯賃金の差

  • Angrist (1990)の論文
  • ベトナム戦争で従軍経験が平均賃金への影響を測定した
  • Angrist の問題意識:
    ・アメリカでは退役軍人には様々な社会保障が与えられていると言われているが、本当に十分な社会保障が与えられているのか?
    ・従軍を経験することが、その後の人生に悪影響があるのではないか?

一般的な見解:「軍隊以外では十分な稼ぎが得られないから、軍隊に行く」

  • 単純な比較の問題点:
    ・「従軍経験のある人の平均賃金」と「従軍経験のない人の平均賃金」の単純比較してみる
    結果:「従軍経験のある人の平均賃金」<「従軍経験のない人の平均賃金」
    → しかし、これが「従軍経験自体によってもたれされた因果効果」とはいえない
    その理由
    → 「従軍経験のある人」はもともと稼げる賃金が低いから、平均賃金が低いだけ
    ・「軍隊に入隊するという意思決定」が「従軍することによって得られる賃金と従軍しないことによって得られる賃金の比較」によってなされる「混合的 (confounded) な割当メカニズム」だから
    = 「処置を受けるかどうか(= 軍隊に入隊するかどうか)の意思決定」と「処置を受けることでもたらされる結果(平均賃金)」は無関係でないから
    =「内生性 (endogeneity) 」があるから

  • 解決策
    ・操作変数を用いた分析
    → 内生性を取り除いて因果効果を推定できる

Angrist の分析

  • ベトナム戦争中にアメリカで一時的に行われていた「くじによる徴兵」に注目
  • 「くじ当選の有無」を操作変数として分析

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

  • くじの実施方法:
    ・「1」から「366」までの番号を書いたくじを作り、壺に入れる
    ・対象者がくじを一枚引く
    ・「1」 = 1月1日、「2」 = 1月2日、・・・、「366」 = 12月31日
    ・引いたくじが自分の誕生日に該当した場合「当選」
    ・予定人数に達するまで、このくじ引きを繰り返す
  • 「当選者」全員が徴兵されるとは限らない
  • 「当選者」は身体検査を受けて、合格すれば徴兵される
  • 「落選者」全員が徴兵されないとは限らない
  • 「落選者」が自発的に志願することも可能
  • 「当選者」の方が「落選者」よりも従軍する確率は高いことは確か
  • 「従軍するか否かの決定」が対象者の意思でなく「くじ」というランダムな割当によって発生したことがポイント
    → くじの当選者の平均賃金と落選者の平均賃金を比較すれば、従軍経験という「処置」がランダムな割当メカニズムによって割り当てられたことになる
    → 因果効果を推定できる(当選者と落選者の従軍確率が100%の場合)
  • 例えば、当選者の従軍確率が 80%、落選者の従軍確率が 20% の場合
  • 因果効果 = 「平均賃金の差」/ (80 − 20)

従軍兵士と兵士以外の生涯賃金の差の分析結果 (Angrist 1990) ・従軍経験によって平均賃金は 15% 減少する

4. 操作変数が交絡を取り去るメカニズム

4.1 ベン図による解説

  • 上で紹介した Angrist の研究成果を使って、操作変数法を可視化して解説する
  • ここで求めたいのは「従軍経験 (\(X_1\))」の「賃金 (\(Y\))」への因果効果 (\(γ_1\))
    \[\hat{Y} = \hat{γ_0} + \hat{γ_1}X_1\]
\(Y\) 結果変数 賃金 観測可
\(X_1\) 介入変数(共変量) 従軍経験 観察可
  • 図の青い部分が因果効果 (\(γ_1\))

  • しかし、通常、従軍経験 (\(X_{1}\))」だけで「賃金 (\(Y\))」が決まるわけではない  
  • 観察できない他の共変量 (\(X_{2}\)) を想定するべき
  • 観察できない他の共変量 (\(X_{2}\)) をモデルに加えて、交絡因子を取り除く

\[\hat{Y} = \hat{β_0} + \hat{β_1}X_1 + \hat{β_2}X_2\]

\(Y\) 結果変数 賃金 観測可
\(X_1\) 介入変数(共変量) 従軍経験 観察可
\(X_2\) 共変量 交絡因子 観察不可
  • 交絡因子を取り除いた因果効果 (\(\hat{β_1}\)) が推定できる

  • しかし、共変量 \(X_2\) は観測できない
    = 未観測の共変量 \(X_2\) をモデルに加えることができない

解決策:

  • 観測できる共変量 Z(=くじ引き = 操作変数)を利用する
  • 操作変数 \(Z\)\(X_1\) に対して直接影響を及ぼす
  • 操作変数 \(Z\)\(Y\) に対して直接影響を及ぼさない
  • 操作変数 \(Z\)\(Y\) に対して直接影響を及ぼさず、\(X_1\) に対してのみ直接影響を及ぼすため、緑色の点線が単独で結果変数 \(Y\) とオーバーラップするところはない
  • 操作変数 \(Z\) を使う
    → 未観測の共変量 \(X_2\) の影響を取り除き、処置効果(= 従軍経験 \(X_1\))を適切に推定できる

4.2 フォーマルな表現による解説

  • 操作変数が交絡を取り去るメカニズムをフォーマルに表現してみる
  • 上で考えたモデルの重回帰式は次のとおり

\[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\]

  • この回帰式には誤差項 \(ε_i\) の代わりに誤差項 \(U_i\) が含まれている
    → 誤差項 \(U_i\) の中に削除された \(X_{2i}\) が含まれている
    → \(X_1\) と 誤差項 \(U_i\) の共分散がゼロではない

\[cov[X_1, U] ≠ 0\]

  • 最小二乗法で推定するためには \(cov[X_1, U] = 0\) である必要がある
    = 「誤差項の条件付き期待値がゼロ」= 誤差項と共変量は独立
    「誤差項の条件付き期待値がゼロ」の仮定を満たさない
    → 最小二乗法によって \(\beta_1\) を推定できない

操作変数 \(Z\) を使った解決策

  • ここで操作変数 \(Z_i\) をモデルに加えることで、\(\beta_1\) を適切に推定できる
「誤差項の条件付き期待値がゼロ」の意味
・誤差項の期待値がゼロ \(E[ε_i|X] = 0\) とは「誤差項と共変量は独立している」という意味であり「統制すべき交絡因子が十分にモデルに含まれている状態」のこと
・ここでは「\(X_2\) は観測できない」と想定しているため、統制すべき交絡因子が含まれていないことになる
→ 偏りのない推定ができない

操作変数 \(Z\) を使った解決策

  • ここで操作変数 \(Z_i\) をモデルに加えることで、\(\beta_1\) を適切に推定できる

操作変数に関する仮定のまとめ

1. 操作変数の外生性 (instrument exogeneity)

\[cov[Z, U]=0\]

  • \(cov[Z, U]=0\)\(Z\)\(U\) の共分散がゼロであるという意味
  • 「共分散がゼロである」= 「相関係数はゼロ」
  • \(Z_i\) に向けて \(X_{2i}\) からも \(Y_i\) からも矢印が出ていない
    = \(Z_i\)\(X_{2i}\) は関係がない

2. 操作変数の関連性 (instrument relevance)

\[cov[Z, X_1] ≠ 0\]

  • \(cov[Z, X_1] ≠ 0\)\(Z\)\(U\) の共分散がゼロでないという意味
  • 「共分散がゼロでない」= 「相関係数はゼロでない」
  • \(Z_i\) から \(X_{1i}\) に矢印が向いている
    = \(Z_i\)\(X_{1i}\) は関係がある

3. 除外制約 (exclusion restriction)

  • 操作変数 \(Z_i\) と誤差項 \(U_i\) と無相関である
  • \(Z_i\) から \(X_{1i}\) に矢印が出ている
  • \(Z_i\) から \(Y_i\) には矢印が出ていない
    \(Z_i\)\(Y_i\) に間接的に影響を与えている
    \(X_i\) の値が決まれば、 \(Z_i\) があってもなくても、\(X_i\) から \(Y_i\) への影響は変わらない

操作変数推定量 \(\beta_1\) を求める式

\[\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]}\]

5. 操作変数のシミュレーション

  • 上で説明した変数間の関係を実際のデータを使ってシミュレーションして確認する

5.1 シミュレーションのための変数の生成

  • ここで生成する変数は次の 4 つ:
\(Y_i\) 結果変数 観測可
\(x_{1i}\) 介入変数(共変量) 観察可
\(x_{2i}\) 共変量(=交絡因子) 観察不可
\(z_i\) 操作変数 観察可
  • ここで推定したいのは共変量 \(x_{1i}\) が結果変数 \(Y_i\) に与える影響
  • 共変量 \(x_{2i}\) は観測できないが、ここではシミュレーションをするため、観測できると想定する

\(x_{2i}\)\(z_i\) の生成

  • まず N = 1000 で 交絡因子 \(x2\) と操作変数 \(z\) を生成する
  • 交絡因子 \(x2\) と操作変数 \(z\) は無関係なので、rnorm()関数を使って、それぞれ個別に生成する
set.seed(1)
n1 <- 1000    # N = 1000 と設定
x2 <- rnorm(n1)  # N = 1000 の x2 を生成する   
z <- runif(n1)   # N = 1000 の z を生成する 
  • 誤差項 \(e1\) を生成する
e1 <- rnorm(n1)
  • 回帰分析における最良線形不偏推定量 (BLUE) 最初の仮定である「誤差項の期待値(=平均値)がゼロ」を確認
mean(e1)
[1] -0.01021977

\(x_{1i}\) の生成

  • 誤差項 \(e1\) を加味して、\(x_{2i}\)\(z_i\) から影響を受けている共変量 \(x_{1i}\) を生成する
x1 <- 1 + 2 * z - 3 * x2 + e1
  • ここでは次のような関係を想定している
    \[x_1 = 1 + 2z_1 - 3x_2 + ε_1\]

\(y_{i}\) の生成

  • 誤差項 \(e2\) を生成する
e2 <- rnorm(n1)
  • 観察されない共変量 \(x_{2i}\) をモデルに含めない時の誤差項 \(u_1\) を作成する
u1 <- x2 + e2
  • \(x_{2i}\)\(x_{1i}\) から影響を受けている共変量 \(y_i\) を生成する
y1 <- 3 + 1.5 * x1 + 1.5 * x2 + e2
  • ここでは次のような関係を想定している

\[y_1 = 3 + 1.5x_1 + 1.5x_2 + ε_2\]

  • 生成した \(y_i\) と共変量 \(x_{1i}\) の散布図を描いてみる
plot(x1, y1)

  • 生成した \(y_i\) と未観測の共変量 \(x_{2i}\) の散布図を描いてみる
plot(x2, y1)

  • 生成した \(z_i\) と未観測の共変量 \(x_{2i}\) の散布図を描いてみる
plot(x2, z)

  • これで、下のような関係の変数が完成したことがわかる

5.2 操作変数に関する 3 つの仮定の確認

  • 次に上で生成したデータが操作変数を使って適切に因果推論するために必要な 3 つの仮定を満たしているかどうかチェックする

1. 操作変数の外生性 (instrument exogeneity) \[cov[Z, U]=0\]

  • \(Z_i\) に向けて \(X_{2i}\) からも \(Y_i\) からも矢印が出ていないこと
  • \(Z_i\)\(X_{2i}\) は関係がないこと
  • 操作変数 \(Z_i\) と観察されない共変量 \(x_{2i}\) をモデルに含めない時の誤差項 \(u_1\) の共分散がゼロであることを確認する
cov(z, u1)
[1] 0.004487261

2. 操作変数の関連性 (instrument relevance) \[cov[Z, X_1] ≠ 0\]

  • \(Z_i\) から \(X_{1i}\) に矢印が向いていること
  • \(Z_i\)\(X_{1i}\) は関係があること
  • \(Z_i\)\(X_{1i}\) の共分散が関係がある(=ゼロでない)ことを確認する
cov(z, x1)
[1] 0.144279
  • 操作変数の関連性のモデルを診断してみる
relevance <- lm(x1 ~ z)
tidy(relevance)
termestimatestd.errorstatisticp.value
(Intercept)1.170.2025.88.78e-09
z1.7 0.3544.81.8e-06 
  • \(z_1\)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\) への影響は変わらない

  • 共変量 \(x_1\) と観察されない共変量 \(x_{2i}\) をモデルに含めない時の誤差項 \(u_1\) の共分散がゼロでないことの確認
cov(x1, u1)
[1] -3.157631

5.3 シミュレーション結果

  • ここでは次の 4 つの変数を生成した
\(Y_i\) 結果変数 観測可
\(x_{1i}\) 処置変数(共変量) 観察可
\(x_{2i}\) 共変量(=交絡因子) 観察不可
\(z_i\) 操作変数 観察可
  • これら 4 つの変数を使って 3 種類のモデルを推定してみる

重回帰モデル(真の値)

m_1 <- lm(y1 ~ x1 + x2)
tidy(m_1)
termestimatestd.errorstatisticp.value
(Intercept)3.070.063148.71.47e-265
x11.490.027354.88.16e-303
x21.470.086817  7.8e-57  

単回帰モデル

m_2 <- lm(y1 ~ x1)
tidy(m_2)
termestimatestd.errorstatisticp.value
(Intercept)3.910.043989.10
x11.060.011493.40

操作変数モデル

m_3 <- cov(z, y1)/cov(z, x1)
tidy(m_3)
x
1.57
  • 推定した 3 つのモデルを推定した結果は次の通り
モデルの種類 推定値
真の推定値 1.49
単回帰モデル 1.06
操作変数モデル 1.57

シミュレーションの結果まとめ真の推定値は 1.49
・単回帰モデルの推定値 (1.06) は真の推定値 (1.49) を過小評価している
操作変数モデルの推定値 (1.57) は真の推定値 (1.49) に近い値を推定している
→ 未観測の共分散 \(x_2\) の交絡を除去している

6. 2 段階最小二乗法

  • 上記と同じ分析を 2 段階最小二乗法 (two-stage least squares: 2SLS) を使っても分析できる

第 1 段階

\[X_i = \alpha_0 + \alpha_1Z_i + e_i\]

  • \(e_i\) は誤差項
  • \(\alpha_0\)\(\alpha_1\) を最小二乗法を使って推定
  • さらにそれらの推定値と \(Z_i\) のデータを使えば \(X_i\) の予測値 \(\hat{X_i}\) を計算できる

第 2 段階

  • \(X_i\) の予測値 \(\hat{X_i}\) を使って第 2 段階の関係を推定する

\[Y_i = \beta_0 + \beta_1\hat{X_i} + u_i\]

  • \(u_i\) は誤差項
step_1 <- lm(x1 ~ z)     # 第 1 段階の最小二乗法
x1hat <- predict(step_1) # x1 の予測値 x1hat を抽出  

step_2 <- lm(y1 ~ x1hat) # 第 2 段階の最小二乗法
tidy(step_2)       # 結果を表示  
termestimatestd.errorstatisticp.value
(Intercept)2.9 0.4776.081.68e-09
x1hat1.570.2316.782.03e-11

2 段階最小二乗法を使った推定値真の推定値は 1.49
操作変数モデルの推定値 (1.57)
2 段階最小二乗法を使った推定値 (1.57)の値は操作変数モデルの推定値と同じ
→ 未観測の共分散 \(x_2\) の交絡を除去している

Reference
  • Angrist, Joshua D. 1990. Lifetime Earnings and the Vietman era draft lottery: Evidence social security administrative records, American Economic Review, 80, 313–336.
  • 矢内勇生「KUT計量経済学応用」
  • 高橋将宣『統計的因果推論の理論と実装』、共立出版、2022年
  • 宋財泫 (Jaehyun Song)・矢内勇生 (Yuki Yanai)「私たちのR: ベストプラクティスの探究」
  • 浅野正彦・ 矢内勇生『Rによる計量政治学』オーム社、2018年
  • 松林哲也 『政治学と因果推論』、2021年
  • 森田果 『実証分析入門』、2014年
  • Kosuke Imai, Quantitative Social Science: An Introduction, Princeton University Press, 2017