第19章工具变量法

19.1 它是如何起作用的？

19.1.1 分离变异

如果你想要建造一座雕像，一种方法是找一大块大理石，然后一点点凿，直到露出下面你想要的形状。另一种方法是拿混凝土或者钢水之类的东西，把它倒进一个模具里，这个模具只允许你想要的形状显现出来。基于关闭后门路径和消除不想要的变异的研究设计，有点像前者。工具变量设计则像后者。

工具变量更像是使用模具，因为它不是试图用控制变量来剥离所有不想要的变异，而是找到一种变异来源，这种来源能让你只分离出你感兴趣的前门路径。这和随机对照实验的运作方式类似。在随机实验中，你可以控制谁接受处理，谁不接受。你可以选择以随机的方式分配处理。当你这样做时，你实验中参与者的处理变异是没有后门路径的。通过使用随机分配的处理，而不是现实世界中的处理，你只分离出了处理和结果之间由于处理本身导致的共变，并且完全清除了由于那些讨厌的后门路径导致的变异。

工具变量设计直接利用了随机对照实验的概念，实际上是试图模仿随机实验，但使用的是统计方法和有利的环境，而非真正能够影响或随机化任何事物¹。

在典型情况下，如果我们有一个 “处理$（Treatment）$” 和一个 “结果$（Outcome）$”，我们想知道处理对结果的影响，即 $处理 \rightarrow 结果$。然而，至少可以说，存在一条烦人的后门路径。这里可能有很多情况，但我们就假设还存在 $Treatment \leftarrow Annoyance \rightarrow Outcome$ 这样的关系。此外，由于在社会科学中，$Annoyance$ 代表着很多不同的事物，我们不太可能把所有这些都控制住。随机实验通过添加一个新的变异来源 ——“随机化$（Randomization）$”—— 来改变这个系统，而”$Randomization$“与 $Annoyance$ 完全无关²。所有这些都在 $图19.1$ 中展示了。

工具变量设计的运作方式完全相同。唯一的区别是，我们不是自己对变量进行随机化，而是希望有什么东西已经为我们随机化了它。我们在现实世界中寻找一个没有后门路径的随机化来源³，并用它来模拟随机对照实验。

工具变量设计并没有消除通过关闭后门路径来识别效应的要求。但它确实转移了这一要求，希望能转移到更容易的事情上！我们不再需要关闭 $Treatment$ 和 $Outcome$ 之间的后门路径（这需要我们控制 $Annoyance$ ），而是只需要关闭 $Randomization$ 和 $Outcome$ 之间的后门路径（如果有的话），以及 $Randomization$ 和 $Outcome$ 之间任何不经过 $Treatment$ 的前门路径。

实际上，实施工具变量的机制，其实是一种在你无法完全控制局面时，试图去做实验所做之事的方法 —— 这是有道理的，因为实际情况就是如此。

我们如何模拟随机实验呢？嗯，在随机实验中，我们产生一些随机变异，然后只分离出由该随机变异驱动的那部分处理 —— 我们只使用实验的数据。所以在工具变量中，我们要在统计上分离出仅由工具变量驱动的那部分处理。

换句话说，我们要：

用工具变量来解释处理变量。
去除处理变量中不由工具变量解释的任何部分。
用工具变量来解释结果变量。
去除结果变量中不由工具变量解释的任何部分⁴。
观察结果变量中剩余的、由工具变量解释的部分与处理变量中剩余的、由工具变量解释的部分之间的关系。

实际上，这与控制一个变量的做法相反。当我们在控制变量 $W$ 的情况下探究 $X$ 对 $Y$ 的效应时，我们用 $W$ 来解释 $X$ 和 $Y$，并去除被解释的部分，因为我们想要关闭任何经过 $W$ 的路径。但当我们用 $Z$ 作为工具变量来探究 $X$ 对 $Y$ 的效应时，我们用 $Z$ 来解释 $X$ 和 $Y$，并去除未被解释的部分，这样就只给我们留下来自 $Z$的路径。由 $Z$ 解释的部分是没有后门路径的部分，这正是我们想要的。

整个过程可以在 $图19.2$ 中看到，该图展示了工具变量设计实际上是如何对数据进行操作以得到估计值的。图中使用了一个二元工具变量 —— 一般来说，工具变量并不一定得是二元的，但这样的话图更容易理解。我们从 $图19.2（a）$ 开始，绘制原始数据，它显示了 $X$ 和 $Y$ 之间存在轻微的正相关关系。我们还能在图的不同部分看到 $Z$ 的不同值 ——$0$ 在右下角，$1$ 在左上角。

接下来看 $图19.2（b）$，它展示了我们在不同 $Z$ 值范围内取 $X$ 的均值。工具变量使用的是 $X$ 中由 $Z$ 解释的部分。当 $Z$ 是二元变量时，就是每个组中 $X$ 的均值。我们只想使用 $X$ 中由 $Z$ 解释的部分 —— 那部分没有任何后门路径，所以在 $图19.2（c）$中，我们去除了 $X$ 中所有不由 $Z$ 解释的变异，只留下那些预测值（组内均值的 $X$ 值）。

然后在 $图19.2（d）$和$（e）$中，我们对 $Y$ 重复这个过程。我们再次这样做，是为了关闭 $X$ 和 $Y$ 之间的后门路径。由 $Z$ 驱动的那种良好变异是没有后门路径的，所以我们想要只分离出那部分。

最后，在$图19.2（f）$中，我们观察 $X$ 的预测部分和 $Y$ 的预测部分之间的关系。由于 $Z$ 是二元的，数据中的每个观测值只能得到两种预测中的一种 ——$X$ 和 $Y$ 的$Z = 0$预测值，或者 $X$ 和 $Y$ 的$Z = 1$预测值。只有两种预测意味着图上只剩下两个点，我们只需在这两个点之间画一条线，就能得到 $X$ 和 $Y$ 之间的关系。这条线的斜率就是用 $Z$ 作为工具变量估计出的 $X$ 对 $Y$ 的效应，也就是只分离出了由 $Z$ 预测的$X/Y$关系的那部分。这里的斜率是负的，与$（a）$中的正相关关系不同。而且因为这些数据是我自己生成的，我可以告诉你，这些数据中 $X$ 和 $Y$ 之间的真实关系实际上是负的，就像$（f）$中的斜率一样。

{
  library(tidyverse)
  library(extrafont)
  library(ggpubr)
  library(cowplot)
  library(fixest)
  library(modelsummary)
}

## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.1.4     ✔ readr     2.1.5
## ✔ forcats   1.0.0     ✔ stringr   1.5.1
## ✔ ggplot2   3.5.1     ✔ tibble    3.2.1
## ✔ lubridate 1.9.4     ✔ tidyr     1.3.1
## ✔ purrr     1.0.4     
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
## Registering fonts with R
## 
## 
## Attaching package: 'cowplot'
## 
## 
## The following object is masked from 'package:ggpubr':
## 
##     get_legend
## 
## 
## The following object is masked from 'package:lubridate':
## 
##     stamp

set.seed(1000)
tb <- tibble(W = rnorm(100), Z = sample(0:1, 100, replace = TRUE)) %>%
  mutate(X = 2*W - 2*Z + rnorm(100)) %>%
  mutate(Y = 2*W - .5*X + rnorm(100)) %>%
  group_by(Z) %>% 
  mutate(mean_X = mean(X),mean_Y = mean(Y))
mns <- tb %>%
  group_by(Z) %>%
  summarize(X = mean(X), Y = mean(Y))

p1 <- ggplot(tb, aes(x = X, y = Y, shape = factor(Z))) + 
  geom_point(color = 'black') +
  scale_shape_manual(values = c(1,19)) + 
  labs(shape = 'Instrument\nValues', title = '(a) Raw Data') +
  theme_pubr() + 
  theme(text         = element_text(size = 10, family="Garamond"),
        axis.text.y = element_blank(),
        axis.text.x = element_blank(),
        axis.title.x = element_text(size = 10, family= "Garamond"),
        axis.title.y = element_text(size = 10, family= "Garamond"),
        axis.ticks = element_blank(),
        legend.position = c(.85,.9),
        legend.background = element_rect(size = .5, color = 'black'))

## Warning: The `size` argument of `element_rect()` is deprecated as of ggplot2 3.4.0.
## ℹ Please use the `linewidth` argument instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

## Warning: A numeric `legend.position` argument in `theme()` was deprecated in ggplot2
## 3.5.0.
## ℹ Please use the `legend.position.inside` argument of `theme()` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

p2 <- ggplot(tb, aes(x = X, y = Y, shape = factor(Z))) + 
  geom_point(color = 'black') +
  geom_vline(aes(xintercept = mns$X[1]), linetype = 'dashed', color = 'black') +
  geom_vline(aes(xintercept = mns$X[2]), linetype = 'dashed', color = 'black') +
  annotate(geom = 'label', x = mns$X[1], y = 2, label = 'X mean for\nZ = 0',family = 'Garamond', size = 13/.pt, hjust = 0) +
  annotate(geom = 'label', x = mns$X[2], y = 2, label = 'X mean for\nZ = 1',family = 'Garamond', size = 13/.pt, hjust = 1) +
  scale_shape_manual(values = c(1,19)) + 
  guides(shape = FALSE) +
  labs(title = '(b) Predict X with Z') +
  theme_pubr() + 
  theme(text         = element_text(size = 10, family="Garamond"),
        axis.text.y = element_blank(),
        axis.text.x = element_blank(),
        axis.title.x = element_text(size = 10, family= "Garamond"),
        axis.title.y = element_text(size = 10, family= "Garamond"),
        axis.ticks = element_blank(),
        legend.position = c(.9,.9),
        legend.background = element_rect(size = .5, color = 'black'))

## Warning: The `<scale>` argument of `guides()` cannot be `FALSE`. Use "none" instead as
## of ggplot2 3.3.4.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

p3 <- ggplot(tb, aes(x = mean_X, y = Y, shape = factor(Z))) + 
  geom_point(color = 'black') +
  geom_vline(aes(xintercept = mns$X[1]), linetype = 'dashed', color = 'black') +
  geom_vline(aes(xintercept = mns$X[2]), linetype = 'dashed', color = 'black') +
  annotate(geom = 'label', x = mns$X[1], y = 2, label = 'X mean for\nZ = 0',family = 'Garamond', size = 13/.pt, hjust = 0) +
  annotate(geom = 'label', x = mns$X[2], y = 2, label = 'X mean for\nZ = 1',family = 'Garamond', size = 13/.pt, hjust = 1) +
  annotate(geom = 'segment', x = -4, xend = mns$X[2]-.1, y = 1, yend = 1, arrow = arrow(length = unit(0.03, "npc"))) +
  annotate(geom = 'segment', x = -1, xend = mns$X[2]+.1, y = 1, yend = 1, arrow = arrow(length = unit(0.03, "npc"))) +
  annotate(geom = 'segment', x = -1, xend = mns$X[1]-.1, y = .5, yend = .5, arrow = arrow(length = unit(0.03, "npc"))) +
  annotate(geom = 'segment', x = 2, xend = mns$X[1]+.1, y = .5, yend = .5, arrow = arrow(length = unit(0.03, "npc"))) +
  scale_shape_manual(values = c(1,19)) + 
  scale_x_continuous(limits = c(min(tb$X),max(tb$X))) +
  guides(shape = FALSE) +
  labs(title = '(c) Only Use Predicted Variation',x = 'X') +
  theme_pubr() + 
  theme(text         = element_text(size = 10, family="Garamond"),
        axis.text.y = element_blank(),
        axis.text.x = element_blank(),
        axis.title.x = element_text(size = 10, family= "Garamond"),
        axis.title.y = element_text(size = 10, family= "Garamond"),
        axis.ticks = element_blank(),
        legend.position = c(.9,.9),
        legend.background = element_rect(size = .5, color = 'black'))

p4 <- ggplot(tb, aes(x = mean_X, y = Y, shape = factor(Z))) + 
  geom_point(color = 'black') +
  geom_hline(aes(yintercept = mns$Y[1]), linetype = 'dashed', color = 'black') +
  geom_hline(aes(yintercept = mns$Y[2]), linetype = 'dashed', color = 'black') +
  annotate(geom = 'label', y = mns$Y[1], x = 4, label = 'Y mean for\nZ = 0',family = 'Garamond', size = 10/.pt, vjust = 1) +
  annotate(geom = 'label', y = mns$Y[2], x = 4, label = 'Y mean for\nZ = 1',family = 'Garamond', size = 10/.pt, vjust = 0) +
  scale_shape_manual(values = c(1,19)) + 
  scale_x_continuous(limits = c(min(tb$X),max(tb$X))) +
  guides(shape = FALSE) +
  labs(title = '(d) Predict Y with Z',x = 'X') +
  theme_pubr() + 
  theme(text         = element_text(size = 10, family="Garamond"),
        axis.text.y = element_blank(),
        axis.text.x = element_blank(),
        axis.title.x = element_text(size = 10, family= "Garamond"),
        axis.title.y = element_text(size = 10, family= "Garamond"),
        axis.ticks = element_blank(),
        legend.position = c(.9,.9),
        legend.background = element_rect(size = .5, color = 'black'))

p5 <- ggplot(tb %>% group_by(Z) %>% slice(1), aes(x = mean_X, y = mean_Y, shape = factor(Z))) + 
  geom_point(color = 'black', size = 3) +
  geom_hline(aes(yintercept = mns$Y[1]), linetype = 'dashed', color = 'black') +
  geom_hline(aes(yintercept = mns$Y[2]), linetype = 'dashed', color = 'black') +
  annotate(geom = 'label', y = mns$Y[1], x = 4, label = 'Y mean for\nZ = 0',family = 'Garamond', size = 10/.pt, vjust = 1) +
  annotate(geom = 'label', y = mns$Y[2], x = 4, label = 'Y mean for\nZ = 1',family = 'Garamond', size = 10/.pt, vjust = 0) +
  annotate(geom='segment', x=mns$X[1], xend=mns$X[1], y = mns$Y[1] + 1.1, yend = mns$Y[1]+.3, arrow = arrow(length = unit(0.03, "npc"))) +
  annotate(geom='segment', x=mns$X[1], xend=mns$X[1], y = mns$Y[1] - 1.6, yend = mns$Y[1]-.3, arrow = arrow(length = unit(0.03, "npc"))) +
  annotate(geom='segment', x=mns$X[2], xend=mns$X[2], y = mns$Y[2] + 1.6, yend = mns$Y[2]+.3, arrow = arrow(length = unit(0.03, "npc"))) +
  annotate(geom='segment', x=mns$X[2], xend=mns$X[2], y = mns$Y[2] - 1.1, yend = mns$Y[2]-.3, arrow = arrow(length = unit(0.03, "npc"))) +
  scale_shape_manual(values = c(1,19)) + 
  scale_x_continuous(limits = c(min(tb$X),max(tb$X))) +
  scale_y_continuous(limits = c(min(tb$Y),max(tb$Y))) +
  guides(shape = FALSE) +
  labs(title = '(e) Only Use Predicted Variation',x = 'X',y = 'Y') +
  theme_pubr() + 
  theme(text         = element_text(size = 10, family="Garamond"),
        axis.text.y = element_blank(),
        axis.text.x = element_blank(),
        axis.title.x = element_text(size = 10, family= "Garamond"),
        axis.title.y = element_text(size = 10, family= "Garamond"),
        axis.ticks = element_blank(),
        legend.position = c(.9,.9),
        legend.background = element_rect(size = .5, color = 'black'))

p6 <- ggplot(tb %>% group_by(Z) %>% slice(1), aes(x = mean_X, y = mean_Y)) + 
  geom_line(size = 1) +
  geom_point(aes(shape = factor(Z)), color = 'black', size = 3) +
  annotate(geom = 'text', x = mns$X[1], y = mns$Y[2], label = 'Slope = IV effect',family = 'Garamond', size = 10/.pt, hjust = 0) +
  scale_shape_manual(values = c(1,19)) + 
  scale_x_continuous(limits = c(min(tb$X),max(tb$X))) +
  scale_y_continuous(limits = c(min(tb$Y),max(tb$Y))) +
  guides(shape = FALSE) +
  labs(title = '(f) Relate Predicted Y to Predicted X',x = 'X',y = 'Y') +
  theme_pubr() + 
  theme(text         = element_text(size = 10, family="Garamond"),
        axis.text.y = element_blank(),
        axis.text.x = element_blank(),
        axis.title.x = element_text(size = 10, family= "Garamond"),
        axis.title.y = element_text(size = 10, family= "Garamond"),
        axis.ticks = element_blank(),
        legend.position = c(.9,.9),
        legend.background = element_rect(size = .5, color = 'black'))

## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

plot_grid(p1,p2,p3,p4,p5,p6, ncol = 2)

ggsave('InstrumentalVariables/iv_anim.pdf',width = 7, height = 9.5, units = 'in', device = cairo_pdf)

19.1.2 工具变量的假设

为了让工具变量$（IVs）$有效，我们必须满足两个假设：工具变量的相关性和工具变量的有效性⁵。

相关性是相当直接的。工具变量的思路是，我们利用由工具变量 $Z$ 解释的处理变量 $X$ 的那部分。但如果 $X$ 没有任何部分能被 $Z$ 解释呢？如果它们完全不相关呢？在这种情况下，工具变量就不起作用了。

这在直觉上是很容易理解的 —— 如果 $X$ 没有能被 $Z$ 解释的部分，我们就无法分离出 $X$ 中由 $Z$ 解释的部分。而且，如果我们按照本章前面描述的步骤⁶来做，那么实际上工具变量是在问：“对于 $X$ 中每一个由 $Z$ 解释的变动，$Y$ 中由 $Z$ 解释的变动有多少？“换句话说，在其基本形式中，工具变量给出的是 $\text{Cov}(Z,Y)/\text{Cov}(Z,X)$⁷。如果 $Z$ 不能解释 $$，那么 $\text{Cov}(Z,X) = 0$，我们就会得到 $\text{Cov}(Z,Y)/0$。而我们不能除以 0⁸，这样估计根本就行不通。

在现实场景中，几乎没有任何相关性是真正为 $0$ 的。但即使 $\text{Cov}(Z,X)$ 只是很小而非零，我们仍然会有问题。如果协方差 $\text{Cov}(Z,X)$ 很小，我们会称 $Z$ 是 $X$ 的弱工具变量。$\text{Cov}(Z,X)$ 的估计值可能会在不同样本间有较大波动 —— 这是抽样变异的本质。例如，如果 $\text{Cov}(Z,X)$ 平均为 $1$，且在不同样本间有一点变化，我们可能会得到从 $0.95$ 到 $1.05$ 的值，这会使我们的 $\text{Cov}(Z,Y)/\text{Cov}(Z,X)$ 工具变量估计值在不同样本间有大约 $10\%$ 的变化。但如果 $\text{Cov}(Z,X)$ 在一个低得多的基准值下遵循相同的范围，比如从 $0.01$ 到 $0.11$，这会使我们的 $\text{Cov}(Z,Y)/\text{Cov}(Z,X)$ 估计值在不同样本间的变化幅度高达 $1100\%$。大幅且波动的估计值！在生活的某些领域，大幅且波动可能是积极的特质，但在统计学中不是。

所以我们需要确保 $Z$ 实际上与 $X$ 相关。值得庆幸的是，这是本书这一半内容里我们能真正有把握去检验的少数假设之一。只需查看 $X$ 和 $Z$ 之间的关系，看看它有多强。关系越强，你对相关性假设就越有信心，估计值在不同样本间的波动也就越小⁹。

或许这两个假设中更棘手的是有效性假设。实际上，有效性假设是指工具变量 $Z$ 是一个自身没有 “未关闭的后门路径” 的变量¹⁰。

或许更准确地说，工具变量 $Z$ 和结果 $Y$ 之间的任何路径，要么必须经过处理变量 $X$，要么是关闭的。要记住，工具变量并没有免除我们关闭所有存在的后门路径的责任；它只是把这种责任从处理变量转移到了工具变量上，而且希望工具变量的后门路径更容易关闭。

以 $图19.3$ 为例。如果我们想用 $Z$ 作为工具变量来识别 $X$ 对 $Y$ 的影响，为了让有效性假设成立，需要满足什么条件？

有几条从 $Z$ 到 $Y$ 的路径是我们需要考虑的：

$Z \to X \to Y$
$Z \to X \leftarrow D \to Y$
$Z \leftarrow C \to X \to Y$
$Z \leftarrow C \to X \leftarrow D \to Y$
$Z \leftarrow B \to Y$
$Z \to A \to Y$

我们希望所有从 $Z$ 到 $Y$ 的开放路径都包含 $X$。列表中的前四条路径包含 $X$，所以这些是可以的。即使存在 $Z \leftarrow C \to X$ 这样的路径也不是问题。如果那条路径是开放的，我们无法识别 $Z \to X$。但谁在乎呢？我们不必识别 $Z \to X$，我们只需要除了通过 $X$ 之外，没有从 $Z$ 到 $Y$ 的途径¹¹ ¹²。

那 $D$ 呢？这是内生性的恼人来源，我们大概无法对其进行控制，而正是它导致 $X \to Y$ 无法识别，并且首先要求我们使用工具变量$（IV）$。但即使我们有 $Z \to X \leftarrow D \to Y$，这也是没问题的。一旦我们分离出由 $Z$ 驱动的 $X$ 中的变异，任何其他指向 $X$ 的箭头基本上就不再重要了。那些箭头可能指向 $X$，但它们不指向 $Z$，而我们只使用从 $Z$ 开始的变异。

把包含 $X$ 的路径放到一边，就剩下两条：$Z \leftarrow B \to Y$ 和 $Z \to A \to Y$。这两条路径是问题，为了让 $Z$ 有效，必须将它们关闭。$Z \leftarrow B \to Y$ 是很明显的一条。它看起来就像我们会担心的任何普通的后门路径。

但为什么 $Z \to A \to Y$ 也是个问题呢？因为它给了 $Z$ 另一种不通过 $X$ 就与 $Y$ 相关的途径。当我们分离出由 $Z$ 驱动的 $X$ 中的变异时，那种变异也会与 $A$ 紧密相关¹³。所以当我们接着去看 $Y$ 和由 $Z$ 驱动的 $X$ 变异之间的关系时，我们会把 $X$ 的效应和 $A$ 的效应混在一起。

所以，为了让 $Z$ 的有效性假设成立，我们需要关闭 $Z \leftarrow B \to Y$ 路径和 $Z \to A \to Y$ 路径。如果我们能对 $A$ 和 $B$ 都进行控制，那我们就没问题了¹⁴。

我们需要相关性和有效性。不管怎样，有效性有多现实呢？理想情况下，我们希望我们的工具变量表现得就像实验中的随机分配一样。但在现实世界里，这种情况实际发生的可能性有多大呢？或者，如果某个工具变量需要控制变量才能有效，我们又能有多确信这些控制真的做到了我们需要它们做的所有事呢？

在很久以前，研究者们很乐意使用工具变量，而不怎么费心思去思考有效性。如果你回溯到 $20$ 世纪 $70$ 年代或 $80$ 年代，你会发现人们用诸如父母教育程度这样的变量作为个人自身（当然，你父母的教育程度除了通过你自己的教育之外，不可能以其他方式影响你的结果！）的工具变量。那时候简直是”狂野西部”（指缺乏规范约束的状态）

但如今，去参加任何一个展示工具变量相关论文的研讨会，你都会听到没完没了的关于这个工具变量是否有效的担忧和争论。而且随着时间的推移，在有效性方面，似乎人们越来越难被说服了¹⁵。

有充分的理由去担忧。不仅很难证明存在一个与处理变量强相关、且不知怎的与处理变量原本存在的所有难以控制的后门路径的来源完全无关的变量，而且我们也有很多这样的历史：那些我们原本觉得听起来相当不错的工具变量，结果证明效果并不好。

在众多例子中，我们可以先看看阿西莫格鲁和约翰逊$（2007年）$的研究。在这项研究中，他们想要弄清楚一个国家的人口健康状况对其经济增长的影响。我们会预期更健康的人口会带来更好的经济增长，对吧？他们把新医疗技术（比如针对某些疾病的疫苗）引入的时间作为一个国家健康状况（具体来说是某些疾病的死亡率）的工具变量。他们发现，某一年因当年医疗技术变化而带来的健康状况变化，实际上对该国的经济增长有负面影响！

Acemoglu, Daron, and Simon Johnson. 2007. “Disease and Development: The Effect of Life Expectancy on Economic Growth.” Journal of Political Economy 115 (6): 925–85.

布鲁姆、坎宁和芬克$（2014）$对阿西莫格鲁和约翰逊所用的工具变量有不同的看法。布鲁姆、坎宁和芬克指出，原始研究没有考虑到健康状况的变化可能对增长产生长期影响的可能性。为什么这可能很重要呢？因为你原本就越健康，新医疗技术能进一步降低你的死亡率的作用就越小。所以，那种医疗技术的有效性与先前存在的健康状况相关，而先前存在的健康状况又与先前存在的经济因素相关，这些经济因素又与当前的经济增长相关。这样就出现了一个后门路径。布鲁姆、坎宁和芬克发现，当他们加入对先前存在因素的控制时，大概解决了有效性问题，阿西莫格鲁和约翰逊所发现的负面影响实际上变成了正面影响¹⁶。

Bloom, David E., David Canning, and Günther Fink. 2014. “Disease and Development Revisited.” Journal of Political Economy 122 (6): 1355–66.

一个特别著名的、曾经流行但后来声名扫地的工具变量例子是降雨量，它在现代农业社会的研究中曾被广泛使用¹⁷。降雨量水平被用作各种事物的工具变量，比如个人收入或经济活动¹⁸。降雨过少（或过多）会损害农业产出。如果人人都是农民，嗯…… 那可不好！由于农民无法控制降雨量，把它当作外生变量且没有任何后门路径似乎是合理的。

为了看一个使用降雨量作为工具变量的研究例子，米格尔、萨蒂亚纳特和塞尔真蒂$（2004年）$将降雨量用作非洲不同国家经济增长的工具变量，然后研究经济增长如何影响内战冲突。他们发现了一个相当大的影响：收入增长 $5$ 个百分点的负向变化，会使第二年发生内战的概率提高 $50\%$。

Miguel, Edward, Shanker Satyanath, and Ernest Sergenti. 2004. “Economic Shocks and Civil Conflict: An Instrumental Variables Approach.” Journal of Political Economy 112 (4): 725–53.

但降雨量是收入的有效工具变量吗？萨森斯$（2015年）$研究了印度的数据，发现在有大坝水灌溉的地区（在这些地区，降雨量对收入没什么影响）和没有大坝的地区（在这些地区，降雨量对收入影响很大），降雨量和冲突之间的关系是相似的。这意味着降雨量肯定是通过收入之外的某种机制影响冲突的¹⁹。她的一些潜在解释包括降雨量促进移民、降雨量通过附近非大坝灌溉区的溢出效应对大坝灌溉区产生影响，或者降雨量使暴动或抗议变得困难。

Sarsons, Heather. 2015. “Rainfall and Conflict: A Cautionary Tale.” Journal of Development Economics 115: 62–72.

如果有效性如此难以实现，那为什么还要费这个劲呢？你无法用数据来证明你的工具变量是有效的，而且基于你认为自己的因果图是正确的这一点来说服别人它是有效的也很困难。话虽如此，这并不意味着我们需要完全放弃工具变量。我们只是可能需要谨慎选择应用它的时机。

首先，确实存在一些情况，在这些情况下工具变量几乎和随机分配一样好。你只需要留意并足够幸运才能找到这样的工具变量！某种力量恰好起到了分配处理的作用，这种力量真正来自系统之外，或者几乎是完全随机应用的。一个真正好的工具变量通常采取以下两种形式之一。要么它代表真正的随机化，比如实际的随机实验，或者像 “孟德尔随机化” 这样的情况，孟德尔随机化利用父母基因结合产生后代的随机过程作为工具变量²⁰。或者（在社会科学中更常见），一个好的工具变量可能是你永远不会想到要纳入结果变量模型中的变量，而且实际上你可能会惊讶地发现它竟然和处理的分配有任何关系。

为了举一个这类情况的例子，我要自吹一下，提一提我参与过的唯一一项运用了有效性假设的研究（戈德哈伯、格劳特和亨廷顿 - 克莱因，$2017$ 年）²¹。这篇论文是关于一个学区的教师招聘流程的，该学区在招聘教师前会用一个评分系统来评估申请者。当我们查看数据时，意识到不同考核标准的评估分数是通过手工相加来得出总分的。有时候，这种相加会出现错误。加法错误是能否进入招聘流程下一阶段的一个影响因素，而且在我看来，这听起来相当随机²²。

Goldhaber, Dan, Cyrus Grout, and Nick Huntington-Klein. 2017. “Screen Twice, Cut Once: Assessing the Predictive Validity of Applicant Selection Tools.” Education Finance and Policy 12 (2): 197–223.

第二，有一些工具变量分析的形式允许有效性假设有一点不成立。所以，只要你能让有效性假设接近成立，你就能从中得到一些有用的结果。我们会在”专业人士怎么做”部分讨论这个。

第三，在某些情境下，随机性是更可信的。一方面，你可以把工具变量应用到实际的随机实验中，而且事实上，当不是所有人都按照你随机分配的去做（“依从性不佳”）时，这是很常见的做法²³。

你是不是把一些人随机分配去服用一种新药，另一些人不服用，但有些人还是会跳过服药剂量？如果你只是按常规分析这个实验，你会低估这种药物的效果，因为你会把那些实际上从未服用过药物的人混在里面。这会是”意向性治疗”估计，也就是分配人们去服药的效果，而不是药物本身的效果。但如果你把随机分配当作服用药物处理的工具变量，你就能得到处理的效果。而且，那种随机分配听起来非常有效。毕竟，它真的是随机的。没有其他路径！

另一个工具变量以非常可信的有效性假设出现的地方是在回归不连续设计的应用中。我们会在第 $20$ 章讲到这些。但现在，让我们更深入地探讨一下，如果我们确实相信有效性假设，我们自己实际上可以如何使用工具变量。

19.1.3 典型设计

工具变量的使用意味着，你可以借助各种各样巧妙的设计，在那些似乎满是难以处理的后门路径的地方进行因果推断。“巧妙”可能是你听到的用来形容一个好工具变量最常见的形容词²⁴。

不过，超级巧妙的工具变量有一个缺点，那就是它们往往高度依赖具体情境。但有些工具变量设计似乎既相当可信，又能应用于很多不同的情况，并且被反复使用²⁵。其中一些我已经讨论过 —— 第 $9$ 章里的特许学校抽签就是一个例子。当然，还有在依从性不佳的随机实验中工具变量的使用。另一个是 —— 哎呀 —— 把降雨量作为农业生产力的工具变量。但其他一些工具变量设计比降雨量更能经受住多年的批评²⁶。

一个典型工具变量设计的好例子是法官分配（艾泽和小多伊尔，$2015年$）。在许多法院系统中，当你即将受审时，为你的案件分配法官的过程或多或少是随机的。这很重要，因为有些法官更严厉，而另一些则没那么严厉。这意味着”法官严厉程度”可以作为”惩罚”的工具变量。只需利用法官之前的裁决来估计每位法官的严厉程度。然后，你（案件）的法官的严厉程度就是对你所受惩罚的工具变量。这种方法可以用在任何有随机分配法官，且你想了解更严厉的惩罚（甚至是被判有罪）对后续某一结果影响的地方。

另一个例子是将征兵作为服兵役的工具变量（安格里斯特，$1990年$）。征兵通常是半随机进行的，人们被征召的顺序是随机分配的。你是否被征召是你在征召顺序中的位置以及他们想要征召的人数的综合结果。所以，你在随机征兵顺序中的位置是服兵役的一个随机来源。它也可以是其他事情的随机来源，比如上大学，这可能让你得以逃避征兵²⁷。

Angrist, Joshua D. 1990. “Lifetime Earnings and the Vietnam Era Draft Lottery: Evidence from Social Security Administrative Records.” The American Economic Review, 313–36.

或许有史以来最常用的工具变量设计是将义务教育年限作为受教育年限的工具变量（安格里斯特和克鲁格，$1991年$）²⁸。有无数因素会影响你决定继续接受教育的时长以及何时停止。但在很多地方，存在义务教育年龄 —— 在达到一定年龄之前，你不能合法辍学。所以，如果你想在 $15$ 岁辍学，但法律要求你必须待到 $16$ 岁，那这对你来说就是额外接受一点教育的推动力。这些法律随时间在不同地区的变化可以作为你所受教育程度的工具变量。

Angrist, Joshua D., and Alan B. Keueger. 1991. “Does Compulsory School Attendance Affect Schooling and Earnings?” The Quarterly Journal of Economics 106 (4): 979–1014.

最后，还有”巴蒂克转移份额工具变量”（巴蒂克，$1991年$）。在这种设计中，整体经济范围内的共同趋势与不同地区的产业分布相结合，作为经济活动的工具变量。换句话说，如果你住在生产大量冰淇淋的冰淇淋镇，旁边的镇子是生产大量爆米花的爆米花镇，而且在过去十年里冰淇淋在全国范围内变得越来越受欢迎，那这对你的镇子来说将是一个相当随机的利好。很简单！近期的研究更深入地探讨了这种设计，并发现我们对它的解读需要有一些改变（戈德史密斯 - 平卡姆、索金和斯威夫特，$2020年$），但这种设计仍在被使用。

Bartik, Timothy J. 1991. Who Benefits from State and Local Economic Development Policies? WE Upjohn Institute for Employment Research. Goldsmith-Pinkham, Paul, Isaac Sorkin, and Henry Swift. 2020. “Bartik Instruments: What, When, Why, and How.” American Economic Review 110 (8): 2586–2624.

当然，还有很多其他典型设计。最初一波移民的定居地点作为进一步移民的工具变量，双胞胎的出生作为子女数量的工具变量，住房供应弹性作为住房价格的工具变量，自然灾害作为各种事物的工具变量…… 我可以继续说下去。在本章中，我会讨论另外三种典型设计：风向作为污染的工具变量，你的前两个孩子是否同性作为生育第三个孩子的工具变量，以及”孟德尔随机化”，它使用某人的基因编码作为各种事物的工具变量。

19.2 它是如何实施的？

19.2.1 工具变量估计量

我们可以从对工具变量估计的方法入手，按其字面意思来理解。作为一种研究设计，工具变量就是要分离出由工具变量解释的处理变量的变异。所以，咱们就，呃，这么做吧。

两阶段最小二乘法$（2SLS）$²⁹是一种通过两次回归来估计工具变量模型的方法。“第一阶段”使用工具变量以及其他控制变量来预测处理变量。然后，你从第一阶段获取处理变量的预测（被解释）值，并在第二阶段用它（再次加上控制变量）来预测结果。

给定我们的工具变量 $Z$、处理变量 $X$、结果变量 $Y$ 和控制变量 $W$，我们会估计以下模型：

\[ X = \gamma_0 + \gamma_1 Z + \gamma_2 W + \nu \quad (19.1) \]

\[ Y = \beta_0 + \beta_1 \hat {X} + \beta_2 W + \varepsilon \quad (19.2) \]

其中，$\nu$ 和 $\varepsilon$ 都是误差项，$\hat{X}$ 是 X 的预测值（通过对第一个方程进行普通最小二乘法估计来预测），$\gamma$ 是回归系数，和 $\beta$ 类似，只是用了不同的希腊字母，以免和 $\beta$ 混淆³⁰。

这个过程手工计算相当容易（不过我不建议这么做，原因我马上会说）。只需用 $Z$ 对 $X$ 进行 $OLS$ 回归（在 $R$ 中用 $lm$，在 $Stata$ 中用 $regress$，在 $Python$ 的 $statsmodels$ 中用 $sm.ols().fit()$）。然后，用该回归的结果来预测 $X$（在 $R$ 中用 $predict()$，在 $Stata$ 中用 $predict$，在 $Python$ 中用 $sm.ols().fit().predict()$）。最后，用预测值对 $Y$ 进行回归。不要忘记在第一阶段和第二阶段都纳入所有控制变量。

我们还没完全结束呢 —— 如果我们只是按照我描述的步骤来做，我们的标准误会是错误的。所以需要进行标准误调整，调整它们以考虑到我们是估计了 $\hat{X}$ 而非对其进行测量这一事实，因此这些值的不确定性比普通最小二乘法所捕捉到的要更大³¹。

但掌握了这一点后，我们就创建出了自己的工具变量估计值！

准确来说，它到底在做什么呢？两阶段最小二乘法生成的是效应的比率，用 $Z$ 对 $Y$ 的效应除以 $Z$ 对 $X$ 的效应。它问的是：“通过改变 $Z$ 我们使 $X$ 产生的每一次变动，会导致 $Y$ 产生多少变动？” 由于 $Z$ 没有后门路径，这个问题的答案应该能给我们 $X$ 对 $Y$ 的因果效应。

$2SLS$ 有一些不错的特点 —— 它易于估计，具有灵活性（在第一阶段添加更多工具变量超级容易，尽管添加更多处理变量没那么容易），而且因为它实际上只是使用 $OLS$，所以很容易理解。这些就是 $2SLS$ 到目前为止是实施工具变量最常用方法的原因。

$2SLS$ 也有一些缺点。首先，它在小样本中的表现不是很好。虽然理论上工具变量没有后门路径，但在实际数据集中，$Z$ 和 $Y$ 的非 $X$ 部分之间的关系至少会因为随机因素而有一点非零。样本越小，这种”因随机因素导致的非$0$“就越不只是非零，还可能相当大，这会导致 $Z$ 在给定样本中不完全有效，并给你带来偏差。此外，当误差存在异方差时，$2SLS$ 的表现也不是特别好³²。

尽管两阶段最小二乘法是思考工具变量最直白的方式，但它只是众多估计量中的一种。而且有充分的理由表明，尽管它很流行，但在不同的估计量中，它甚至不是一个那么好的选择。

在这里，我会稍微讲讲用广义矩方法$（GMM）$来估计工具变量的方法。其他一些各有其优势的方法，会在本章中陆续出现。

广义矩方法$（GMM）$是一种比工具变量应用范围广泛得多的估计方法，但在本章中，至少我们只是将它用于工具变量³³。其基本思路是：基于你的假设和理论，构建一些应该具有特定值的统计矩（均值、方差、协方差等）。

例如，如果我们想用 $GMM$ 来估计一个变量 $Y$ 的期望值（基本上就是均值），我们会说，对期望值 $\mu$ 的样本估计与总体期望值 $E(Y)$ 之间的差异应该为零。所以我们会把 $\mu - E(Y) = 0$ 作为我们估计的一个条件。用样本值 $\frac{1}{N}\sum Y$ 代替 $E(Y)$，然后解方程得到 $\hat{\mu} = \frac{1}{N}\sum Y$。$GMM$ 会选出那个从平均意义上使矩条件成立的 $\mu$。或者对于普通最小二乘法来说，我们假设 $X$ 与误差项 $\varepsilon$ 不相关。所以我们使用 $Cov(X, \varepsilon) = 0$ 这个条件。在实际数据中，这相当于 $\sum Xr = 0$，其中 $r$ 是残差。做一点数学运算，你最终会得到和我们之前一样的 $\beta_1$ 的估计值³⁴ ³⁵。

从 $OLS$ 的例子中，你可能能猜到 $GMM$ 在工具变量中是如何运作的。我们假设 $Z$ 与第二阶段的误差项 $\varepsilon$ 不相关。所以我们使用 $Cov(Z, \varepsilon) = 0$ 这个条件。代入协方差的样本数据版本，再做一点数学运算，最终得到 $\beta_1 = \sum (ZY) / \sum (ZX)$，和我们之前得到的结果一样。

所以，答案是一样的。没什么大不了的，对吧？只是当我们存在异方差（我们很可能存在异方差）或者添加更多工具变量时，情况会变得更有意思一些。

$GMM$ 和两阶段最小二乘法$（2SLS）$都能够处理异方差问题。$GMM$ 能自然地做到这一点，因为它并不像构成 $2SLS$ 的 $OLS$ 方程那样对误差项做出假设。不过，我们也可以通过简单地使用异方差稳健的标准误，让 $2SLS$ 也能处理异方差问题。

但当工具变量的数量多于你所拥有的处理 / 内生变量的数量时，情况会变得更有意思，这种情况被称为”过度识别”。

过度识别：工具变量模型中，工具变量的数量多于处理 / 内生变量的数量。恰好识别：工具变量模型中，工具变量的数量与处理 / 内生变量的数量相同。

当你的模型是过度识别的时候，两阶段最小二乘法和广义矩方法就会有分歧³⁶。至少在涉及异方差的情况下，$GMM$ 会更精确。在存在异方差时，即便你添加了异方差稳健的标准误，$GMM$ 的抽样变异也会比 $2SLS$ 更小（因此标准误也更小）。就算你开始应用聚类标准误或者对自相关进行修正，这一点仍然成立。

要记住 ——$GMM$ 不只是调整标准误的另一种方式。当存在过度识别时，估计值本身实际上也会不同。$GMM$ 的标准误更小，不是因为我们选择宣称自己有更多信息，从而更精确，而是因为这种方法本身就会产生更精确的估计值，而更小的标准误反映了这一点。

让我们来编写一些工具变量的代码。在这个练习中，我们将仿照蔡、德扬弗里和萨杜莱$（2015年）$的论文《社会网络与投保决策》进行操作。作者们正在研究农民就是否购买防范天气事件的保险所做的决策。具体来说，他们感兴趣的是关于保险的信息是否通过社会网络传播。

Cai, Jing, Alain De Janvry, and Elisabeth Sadoulet. 2015. “Social Networks and the Decision to Insure.” American

他们研究了中国农村的一项随机实验，在该实验中，家庭被随机分配到两轮不同的关于保险的信息宣讲会中。接下来的问题是：你的朋友们了解到的关于保险的情况，在多大程度上会影响你自己的保险购买行为？他们观察了参加第二轮宣讲会的人，以及这些人在第一轮宣讲会中的朋友的所作所为和所见所闻，以此来看看前者是如何受到后者影响的³⁷。

蔡、德扬弗里和萨杜莱大体上确实发现，农民的决策受到了他们朋友的所见以及朋友所接收信息的影响。我将关注他们的一项特定分析，在这项分析中，他们探究农民的决策是否受到其朋友行为的影响：你朋友实际购买保险的行为会让你更有可能购买保险吗？

我们想要识别在参加第二轮信息宣讲会的人群中，朋友的购买行为$（FriendsPurchaseBehavior）$对自身购买行为$（YourPurchaseBehavior）$的影响，观察的是他们那些参加了第一轮信息宣讲会的朋友的平均购买行为。这种影响存在一些明显的后门路径。不同地区对保险的偏好可能有高有低，或者你可能只是更有可能拥有与你偏好相似的朋友，包括在保险这类话题上的偏好。

作为朋友购买行为$（FriendsPurchaseBehavior）$的工具变量，他们使用了”第一轮默认设置”$（FirstRoundDefault）$这个变量，它是一个二元指标，用于表示你的朋友是否被随机分配到”默认购买”的信息宣讲会（在该宣讲会中，参与者被默认分配为购买保险，若不买则需明确表明不买的偏好），或者”默认不购买”的宣讲会（在该宣讲会中，参与者被默认分配为不购买保险，若要买则需明确表明购买的偏好）。每个人都有相同的选项，也得到了相同的信息，但默认设置不同。人们会遵循默认设置！参加”默认购买”宣讲会的人比参加”默认不购买”宣讲会的人购买保险的可能性高出 $12$ 个百分点。“第一轮默认设置”是随机分配的这一事实，使得它是一个有效工具变量的说法相当可信。此外，$12$ 个百分点的增幅似乎足以满足相关性假设³⁸。

好了，现在，终于要编写代码了。我们将从两阶段最小二乘法$（2SLS）$开始，然后使用相同的数据，进行广义矩方法$（GMM）$、有限信息极大似然法$（LIML）$以及带有固定效应的工具变量（“面板工具变量”）分析。等等，我把 $LIML$ 加进来了 —— 那是什么？那是有限信息极大似然法。我会在这里展示如何编写它的代码，因为切换方法很容易，但我实际上会在本章后面讨论这种方法。

# There are many ways to run 2SLS; 
# the most common is ivreg from the AER package. But we'll use feols from fixest for speed and ease of fixed-effects additions later
library(tidyverse); library(modelsummary); library(fixest)

d <- causaldata::social_insure
# Include just the outcome and controls first, then endogenous ~ instrument in the second part, and for this study we cluster on address
m <- feols(takeup_survey ~ male + age + agpop + ricearea_2010 + literacy + intensive + risk_averse + disaster_prob + 
             factor(village) | pre_takeup_rate ~ default, cluster = ~address, data = d)

## NOTE: 32 observations removed because of NA values (RHS: 32).

# Show the first and second stage, omitting all the controls for ease of visibility
msummary(list('First Stage' = m$iv_first_stage$pre_takeup_rate,'Second Stage' = m), coef_map = c(default = 'First Round Default',
                fit_pre_takeup_rate = "Friends' Purchase Behavior"), 
         stars = c('*' = .1, '**' = .05, '***' = .01))

	First Stage	Second Stage
* p < 0.1, p < 0.05, * p < 0.01
First Round Default	0.118***
	(0.034)
Friends' Purchase Behavior		0.791***
		(0.273)
Num.Obs.	1378	1378
R2	0.469	0.127
R2 Adj.	0.448	0.092
AIC	-775.7	1911.4
BIC	-498.6	2188.5
RMSE	0.18	0.47
Std.Errors	by: address	by: address

import pandas as pd
from linearmodels.iv import IV2SLS
from causaldata import social_insure

d = social_insure.load_pandas().data
# Add a [endogenous ~ instrument] segment to the formula
m = IV2SLS.from_formula('''takeup_survey ~ male + age + agpop + ricearea_2010 + literacy + intensive + risk_averse + disaster_prob + 
C(village) + [pre_takeup_rate ~ default]''', data = d)

## C:\Users\capen\AppData\Roaming\Python\Python313\site-packages\linearmodels\iv\model.py:557: MissingValueWarning: 
## Inputs contain missing values. Dropping rows with missing observations.
##   super().__init__(

# since we want to cluster, and will use m.notnull to see which observations to drop
second_stage = m.fit(cov_type = 'clustered', clusters = d['address'][m.notnull])
# If we want the first stage we must do it ourselves! 
# move the endogenous variable to the dependent position and make the instrument a predictor, removing the []
first_stage = IV2SLS.from_formula('''
pre_takeup_rate ~ male + age + agpop + ricearea_2010 + literacy + intensive + risk_averse + disaster_prob + C(village) + default''', 
data = d).fit(cov_type = 'clustered', clusters = d['address'][m.notnull])

## C:\Users\capen\AppData\Roaming\Python\Python313\site-packages\linearmodels\iv\model.py:557: MissingValueWarning: 
## Inputs contain missing values. Dropping rows with missing observations.
##   super().__init__(

first_stage

OLS Estimation Summary
Dep. Variable:	pre_takeup_rate	R-squared:	0.4685
Estimator:	OLS	Adj. R-squared:	0.4476
No. Observations:	1378	F-statistic:	1.163e+16
Date:	周日, 9月 07 2025	P-value (F-stat)	0.0000
Time:	18:16:35	Distribution:	chi2(52)
Cov. Estimator:	clustered

Parameter Estimates
	Parameter	Std. Err.	T-stat	P-value	Lower CI	Upper CI
Intercept	0.1248	0.1581	0.7897	0.4297	-0.1850	0.4347
male	-0.0062	0.0165	-0.3743	0.7082	-0.0386	0.0262
age	0.0004	0.0004	1.0166	0.3094	-0.0004	0.0012
agpop	-0.0016	0.0030	-0.5518	0.5811	-0.0075	0.0042
ricearea_2010	-0.0004	0.0002	-1.8801	0.0601	-0.0008	1.672e-05
literacy	0.0121	0.0137	0.8844	0.3765	-0.0148	0.0390
intensive	0.0146	0.0050	2.9028	0.0037	0.0048	0.0245
risk_averse	-0.0213	0.0173	-1.2313	0.2182	-0.0551	0.0126
disaster_prob	4.229e-05	0.0002	0.1790	0.8579	-0.0004	0.0005
C(village)[T.beixing]	0.1815	0.1540	1.1784	0.2386	-0.1204	0.4834
C(village)[T.caijia]	0.1482	0.1665	0.8904	0.3733	-0.1781	0.4746
C(village)[T.daqiao]	0.5772	0.1544	3.7381	0.0002	0.2746	0.8799
C(village)[T.daxi]	0.1425	0.1594	0.8935	0.3716	-0.1700	0.4550
C(village)[T.dayu]	0.2003	0.1834	1.0922	0.2748	-0.1592	0.5598
C(village)[T.dazhou]	0.1794	0.1641	1.0932	0.2743	-0.1423	0.5011
C(village)[T.dongan]	0.2576	0.2025	1.2719	0.2034	-0.1393	0.6545
C(village)[T.dukou]	0.2478	0.1779	1.3932	0.1636	-0.1008	0.5964
C(village)[T.fusheng]	0.3864	0.2032	1.9017	0.0572	-0.0118	0.7846
C(village)[T.fuzhou]	0.4337	0.1858	2.3344	0.0196	0.0696	0.7979
C(village)[T.gangtou]	0.2977	0.1857	1.6032	0.1089	-0.0663	0.6616
C(village)[T.gangxia]	0.5351	0.1585	3.3768	0.0007	0.2245	0.8457
C(village)[T.guojia]	-0.0256	0.1552	-0.1652	0.8688	-0.3298	0.2785
C(village)[T.hefeng]	0.0693	0.1648	0.4206	0.6740	-0.2537	0.3923
C(village)[T.helin]	0.3530	0.1771	1.9934	0.0462	0.0059	0.7000
C(village)[T.hongxing]	0.1347	0.1649	0.8168	0.4140	-0.1885	0.4580
C(village)[T.huangshan]	0.3362	0.1790	1.8783	0.0603	-0.0146	0.6870
C(village)[T.jingang]	0.2112	0.1682	1.2559	0.2092	-0.1184	0.5408
C(village)[T.jinggang]	0.0998	0.2066	0.4833	0.6289	-0.3050	0.5047
C(village)[T.lianqian]	0.2425	0.1615	1.5012	0.1333	-0.0741	0.5590
C(village)[T.lianxing]	0.4722	0.1669	2.8293	0.0047	0.1451	0.7994
C(village)[T.longqing]	0.5360	0.1746	3.0705	0.0021	0.1939	0.8781
C(village)[T.lusikou]	0.1626	0.1551	1.0479	0.2947	-0.1415	0.4666
C(village)[T.mazhou]	0.4843	0.1743	2.7784	0.0055	0.1427	0.8260
C(village)[T.minzhu]	0.0254	0.1646	0.1545	0.8772	-0.2972	0.3481
C(village)[T.shigang]	0.1683	0.1682	1.0012	0.3168	-0.1612	0.4979
C(village)[T.wanshang]	0.4769	0.1542	3.0930	0.0020	0.1747	0.7791
C(village)[T.xiabao]	0.6779	0.1604	4.2261	0.0000	0.3635	0.9924
C(village)[T.xianghu]	0.4876	0.1819	2.6809	0.0073	0.1311	0.8441
C(village)[T.xiaofang]	0.2308	0.2422	0.9529	0.3406	-0.2439	0.7054
C(village)[T.xiashatou]	0.2166	0.1542	1.4051	0.1600	-0.0855	0.5188
C(village)[T.xiecheng]	0.2749	0.1577	1.7432	0.0813	-0.0342	0.5841
C(village)[T.xihe]	0.1040	0.1698	0.6122	0.5404	-0.2289	0.4368
C(village)[T.xihu]	0.4205	0.1550	2.7130	0.0067	0.1167	0.7243
C(village)[T.xilian]	0.0053	0.1877	0.0284	0.9773	-0.3626	0.3733
C(village)[T.xinguang]	0.3204	0.1748	1.8327	0.0668	-0.0222	0.6630
C(village)[T.xingzeng]	0.2808	0.1544	1.8179	0.0691	-0.0219	0.5835
C(village)[T.xinlian]	0.0570	0.1577	0.3612	0.7180	-0.2522	0.3661
C(village)[T.yanjiang]	0.3405	0.1694	2.0096	0.0445	0.0084	0.6726
C(village)[T.yazhou]	-0.0276	0.1655	-0.1669	0.8675	-0.3521	0.2968
C(village)[T.yongfeng]	0.4131	0.1877	2.2013	0.0277	0.0453	0.7809
C(village)[T.zhangxi]	0.1491	0.1619	0.9207	0.3572	-0.1683	0.4664
C(village)[T.zixi]	0.1562	0.1717	0.9093	0.3632	-0.1804	0.4927
default	0.1180	0.0337	3.5049	0.0005	0.0520	0.1840

id: 0x259ff871fd0

second_stage

IV-2SLS Estimation Summary
Dep. Variable:	takeup_survey	R-squared:	0.1265
Estimator:	IV-2SLS	Adj. R-squared:	0.0923
No. Observations:	1378	F-statistic:	3.484e+15
Date:	周日, 9月 07 2025	P-value (F-stat)	0.0000
Time:	18:16:35	Distribution:	chi2(53)
Cov. Estimator:	clustered

Parameter Estimates
	Parameter	Std. Err.	T-stat	P-value	Lower CI	Upper CI
male	0.0375	0.0526	0.7143	0.4751	-0.0655	0.1406
age	0.0052	0.0013	4.1541	0.0000	0.0028	0.0077
agpop	-0.0024	0.0066	-0.3582	0.7202	-0.0152	0.0105
ricearea_2010	0.0028	0.0011	2.4895	0.0128	0.0006	0.0049
literacy	0.0821	0.0350	2.3457	0.0190	0.0135	0.1507
intensive	-0.0090	0.0261	-0.3449	0.7302	-0.0602	0.0422
risk_averse	0.2367	0.0442	5.3549	0.0000	0.1501	0.3234
disaster_prob	0.0013	0.0007	1.7254	0.0844	-0.0002	0.0027
C(village)[beilian]	0.0205	0.1255	0.1632	0.8703	-0.2255	0.2664
C(village)[beixing]	-0.4723	0.1447	-3.2636	0.0011	-0.7559	-0.1886
C(village)[caijia]	-0.3722	0.1416	-2.6276	0.0086	-0.6498	-0.0946
C(village)[daqiao]	-0.4732	0.2398	-1.9736	0.0484	-0.9431	-0.0033
C(village)[daxi]	-0.2504	0.1447	-1.7304	0.0836	-0.5341	0.0332
C(village)[dayu]	-0.3084	0.1473	-2.0931	0.0363	-0.5972	-0.0196
C(village)[dazhou]	-0.4851	0.1411	-3.4376	0.0006	-0.7616	-0.2085
C(village)[dongan]	-0.4074	0.1624	-2.5087	0.0121	-0.7257	-0.0891
C(village)[dukou]	-0.4782	0.1754	-2.7268	0.0064	-0.8220	-0.1345
C(village)[fusheng]	-0.4183	0.2067	-2.0236	0.0430	-0.8235	-0.0132
C(village)[fuzhou]	-0.4174	0.2206	-1.8918	0.0585	-0.8499	0.0150
C(village)[gangtou]	-0.3962	0.1705	-2.3233	0.0202	-0.7303	-0.0620
C(village)[gangxia]	-0.1921	0.2214	-0.8679	0.3854	-0.6260	0.2418
C(village)[guojia]	-0.0886	0.0893	-0.9922	0.3211	-0.2637	0.0864
C(village)[hefeng]	-0.3479	0.1316	-2.6428	0.0082	-0.6059	-0.0899
C(village)[helin]	-0.3470	0.2512	-1.3810	0.1673	-0.8394	0.1455
C(village)[hongxing]	-0.3217	0.1410	-2.2808	0.0226	-0.5981	-0.0453
C(village)[huangshan]	-0.3733	0.2030	-1.8389	0.0659	-0.7711	0.0246
C(village)[jingang]	-0.4446	0.1502	-2.9598	0.0031	-0.7390	-0.1502
C(village)[jinggang]	-0.1799	0.1257	-1.4312	0.1524	-0.4263	0.0665
C(village)[lianqian]	-0.3318	0.1614	-2.0551	0.0399	-0.6482	-0.0154
C(village)[lianxing]	-0.3106	0.2187	-1.4204	0.1555	-0.7392	0.1180
C(village)[longqing]	-0.3650	0.2225	-1.6405	0.1009	-0.8010	0.0711
C(village)[lusikou]	0.0362	0.1245	0.2909	0.7712	-0.2078	0.2802
C(village)[mazhou]	-0.3578	0.2170	-1.6487	0.0992	-0.7831	0.0675
C(village)[minzhu]	-0.2737	0.1163	-2.3536	0.0186	-0.5017	-0.0458
C(village)[shigang]	-0.3159	0.2205	-1.4327	0.1519	-0.7482	0.1163
C(village)[wanshang]	-0.5445	0.2154	-2.5273	0.0115	-0.9667	-0.1222
C(village)[xiabao]	-0.4267	0.2747	-1.5535	0.1203	-0.9650	0.1116
C(village)[xianghu]	-0.1475	0.2043	-0.7221	0.4702	-0.5479	0.2529
C(village)[xiaofang]	-0.2530	0.1420	-1.7816	0.0748	-0.5312	0.0253
C(village)[xiashatou]	-0.1648	0.1518	-1.0854	0.2777	-0.4624	0.1328
C(village)[xiecheng]	-0.3672	0.1708	-2.1504	0.0315	-0.7020	-0.0325
C(village)[xihe]	-0.4123	0.1374	-3.0013	0.0027	-0.6816	-0.1431
C(village)[xihu]	-0.4864	0.2062	-2.3587	0.0183	-0.8905	-0.0822
C(village)[xilian]	-0.1305	0.1102	-1.1845	0.2362	-0.3465	0.0855
C(village)[xinguang]	-0.3958	0.1665	-2.3779	0.0174	-0.7221	-0.0696
C(village)[xingzeng]	-0.2510	0.1786	-1.4054	0.1599	-0.6011	0.0990
C(village)[xinlian]	-0.4622	0.1214	-3.8055	0.0001	-0.7002	-0.2241
C(village)[yanjiang]	-0.4406	0.2054	-2.1445	0.0320	-0.8432	-0.0379
C(village)[yazhou]	-0.2693	0.0940	-2.8649	0.0042	-0.4535	-0.0851
C(village)[yongfeng]	-0.4882	0.1907	-2.5599	0.0105	-0.8619	-0.1144
C(village)[zhangxi]	-0.2801	0.1476	-1.8978	0.0577	-0.5693	0.0092
C(village)[zixi]	-0.3247	0.1529	-2.1235	0.0337	-0.6245	-0.0250
pre_takeup_rate	0.7911	0.2671	2.9617	0.0031	0.2676	1.3146

Endogenous: pre_takeup_rate
Instruments: default
Clustered Covariance (One-Way)
Debiased: False
Num Clusters: 166
id: 0x259ff8716a0

dotext = "
  causaldata social_insure.dta, use clear download
  
  * We want village fixed effects, but that's currently a string
  encode village, g(villid)
  * The order doesn't matter, but we need controls here as well as (endogenous = instrument)
  * don't forget to specify the estimator 2sls! and we cluster on address and also show the first stage with the first option
  ivregress 2sls takeup_survey (pre_takeup_rate = default) male age agpop ///
      ricearea_2010 literacy intensive risk_averse disaster_prob i.villid, cluster(address) first
"
library(RStata)

stata(dotext, data.in=NULL, data.out=TRUE)

## . 
## .   causaldata social_insure.dta, use clear download
## .   
## .   * We want village fixed effects, but that's currently a string
## .   encode village, g(villid)
## .   * The order doesn't matter, but we need controls here as well as (endogenou
## > s = instrument)
## .   * don't forget to specify the estimator 2sls! and we cluster on address and
## >  also show the first stage with the first option
## .   ivregress 2sls takeup_survey (pre_takeup_rate = default) male age agpop ///
## >       ricearea_2010 literacy intensive risk_averse disaster_prob i.villid, cl
## > uster(address) first
## 
## First-stage regressions
## -----------------------
## 
##                                                 Number of obs     =      1,378
##                                                 N. of clusters    =        166
##                                                 F(  49,   1325)   =    4925.91
##                                                 Prob > F          =     0.0000
##                                                 R-squared         =     0.4685
##                                                 Adj R-squared     =     0.4476
##                                                 Root MSE          =     0.1792
## 
## ------------------------------------------------------------------------------
##              |               Robust
## pre_takeup~e |      Coef.   Std. Err.      t    P>|t|     [95% Conf. Interval]
## -------------+----------------------------------------------------------------
##         male |  -.0061916   .0169212    -0.37   0.714    -.0393869    .0270037
##          age |   .0004062   .0004087     0.99   0.320    -.0003956    .0012079
##        agpop |  -.0016431   .0030458    -0.54   0.590    -.0076183    .0043321
## riceare~2010 |  -.0003936   .0002142    -1.84   0.066    -.0008138    .0000265
##     literacy |     .01214   .0140415     0.86   0.387     -.015406    .0396861
##    intensive |   .0146403    .005159     2.84   0.005     .0045196     .024761
##  risk_averse |  -.0212538    .017656    -1.20   0.229    -.0558905    .0133829
## disaster_p~b |   .0000423   .0002417     0.17   0.861    -.0004318    .0005164
##              |
##       villid |
##     beixing  |   .1815029   .1575562     1.15   0.250     -.127584    .4905897
##      caijia  |   .1482456   .1703072     0.87   0.384    -.1858555    .4823467
##      daqiao  |   .5772212   .1579522     3.65   0.000     .2673575    .8870848
##        daxi  |    .142465   .1631002     0.87   0.383    -.1774978    .4624279
##        dayu  |    .200324   .1876202     1.07   0.286    -.1677411    .5683892
##      dazhou  |   .1794321   .1678957     1.07   0.285    -.1499382    .5088025
##      dongan  |   .2575805   .2071531     1.24   0.214    -.1488033    .6639644
##       dukou  |   .2477883    .181935     1.36   0.173    -.1091239    .6047004
##     fusheng  |   .3863708   .2078294     1.86   0.063    -.0213398    .7940814
##      fuzhou  |   .4337347   .1900602     2.28   0.023      .060883    .8065863
##     gangtou  |   .2976852    .189939     1.57   0.117    -.0749288    .6702991
##     gangxia  |   .5351047   .1620959     3.30   0.001     .2171122    .8530973
##      guojia  |  -.0256415   .1587544    -0.16   0.872    -.3370789    .2857959
##      hefeng  |   .0693112   .1685605     0.41   0.681    -.2613633    .3999857
##       helin  |   .3529597   .1811225     1.95   0.052    -.0023584    .7082778
##    hongxing  |   .1347103   .1687035     0.80   0.425    -.1962448    .4656654
##   huangshan  |   .3361827   .1830855     1.84   0.067    -.0229864    .6953518
##     jingang  |   .2112013   .1720256     1.23   0.220    -.1262709    .5486735
##    jinggang  |   .0998288   .2112916     0.47   0.637    -.3146737    .5143314
##    lianqian  |   .2424598   .1652056     1.47   0.142    -.0816333    .5665528
##    lianxing  |   .4722216   .1707286     2.77   0.006     .1372938    .8071495
##    longqing  |    .536004    .178564     3.00   0.003      .185705     .886303
##     lusikou  |   .1625567   .1586862     1.02   0.306     -.148747    .4738603
##      mazhou  |   .4843098   .1783063     2.72   0.007     .1345164    .8341031
##      minzhu  |   .0254326   .1683929     0.15   0.880    -.3049132    .3557784
##     shigang  |   .1683445   .1720019     0.98   0.328    -.1690812    .5057702
##    wanshang  |   .4769198   .1577239     3.02   0.003     .1675039    .7863357
##      xiabao  |    .677943   .1640924     4.13   0.000     .3560337    .9998523
##     xianghu  |   .4876279   .1860538     2.62   0.009     .1226357    .8526201
##    xiaofang  |   .2307609   .2477087     0.93   0.352    -.2551832    .7167049
##   xiashatou  |   .2166193   .1576941     1.37   0.170    -.0927381    .5259766
##    xiecheng  |   .2749492    .161342     1.70   0.089    -.0415644    .5914627
##        xihe  |   .1039802   .1737241     0.60   0.550     -.236824    .4447845
##        xihu  |   .4205349   .1585568     2.65   0.008     .1094852    .7315847
##      xilian  |   .0053379   .1920303     0.03   0.978    -.3713786    .3820544
##    xinguang  |   .3203665   .1788064     1.79   0.073    -.0304081     .671141
##    xingzeng  |   .2807815   .1579875     1.78   0.076    -.0291514    .5907144
##     xinlian  |    .056969   .1613345     0.35   0.724      -.25953     .373468
##    yanjiang  |    .340501   .1733142     1.96   0.050     .0005007    .6805012
##      yazhou  |  -.0276229   .1693311    -0.16   0.870    -.3598093    .3045634
##    yongfeng  |   .4131141   .1919647     2.15   0.032     .0365261     .789702
##     zhangxi  |   .1490763   .1656292     0.90   0.368    -.1758478    .4740004
##        zixi  |   .1561528   .1756664     0.89   0.374    -.1884618    .5007674
##              |
##      default |   .1180257   .0344461     3.43   0.001      .050451    .1856005
##        _cons |   .1248373   .1616991     0.77   0.440    -.1923769    .4420516
## ------------------------------------------------------------------------------
## 
## 
## Instrumental variables (2SLS) regression          Number of obs   =      1,378
##                                                   Wald chi2(52)   =    5218.44
##                                                   Prob > chi2     =     0.0000
##                                                   R-squared       =     0.1265
##                                                   Root MSE        =     .46586
## 
##                               (Std. Err. adjusted for 166 clusters in address)
## ------------------------------------------------------------------------------
##              |               Robust
## takeup_sur~y |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
## -------------+----------------------------------------------------------------
## pre_takeup~e |    .791097    .267112     2.96   0.003      .267567    1.314627
##         male |     .03754   .0525577     0.71   0.475    -.0654712    .1405511
##          age |   .0052447   .0012625     4.15   0.000     .0027702    .0077192
##        agpop |   -.002354   .0065718    -0.36   0.720    -.0152344    .0105265
## riceare~2010 |   .0027682   .0011119     2.49   0.013     .0005888    .0049475
##     literacy |   .0821087   .0350045     2.35   0.019     .0135011    .1507163
##    intensive |  -.0090163     .02614    -0.34   0.730    -.0602497    .0422171
##  risk_averse |   .2367229    .044207     5.35   0.000     .1500788    .3233671
## disaster_p~b |   .0012724   .0007374     1.73   0.084    -.0001729    .0027177
##              |
##       villid |
##     beixing  |  -.4927345   .0744112    -6.62   0.000    -.6385779   -.3468912
##      caijia  |  -.3926604   .0701404    -5.60   0.000     -.530133   -.2551877
##      daqiao  |  -.4936617   .1748929    -2.82   0.005    -.8364455   -.1508778
##        daxi  |  -.2709225   .0851532    -3.18   0.001    -.4378196   -.1040254
##        dayu  |   -.328883   .1065975    -3.09   0.002    -.5378103   -.1199557
##      dazhou  |  -.5055345   .0976986    -5.17   0.000    -.6970202   -.3140487
##      dongan  |   -.427896   .1012469    -4.23   0.000    -.6263363   -.2294557
##       dukou  |  -.4987223   .1087694    -4.59   0.000    -.7119065   -.2855381
##     fusheng  |  -.4388294   .1428507    -3.07   0.002    -.7188116   -.1588472
##      fuzhou  |  -.4378898   .1651094    -2.65   0.008    -.7614982   -.1142813
##     gangtou  |  -.4166364   .0998372    -4.17   0.000    -.6123136   -.2209592
##     gangxia  |   -.212618   .1514221    -1.40   0.160    -.5093998    .0841638
##      guojia  |  -.1091161   .0540893    -2.02   0.044    -.2151291    -.003103
##      hefeng  |  -.3683648   .0832482    -4.42   0.000    -.5315282   -.2052014
##       helin  |  -.3674614   .2044574    -1.80   0.072    -.7681905    .0332676
##    hongxing  |  -.3421713   .0864178    -3.96   0.000     -.511547   -.1727956
##   huangshan  |  -.3937542   .1429747    -2.75   0.006    -.6739795   -.1135288
##     jingang  |    -.46507    .093309    -4.98   0.000    -.6479523   -.2821876
##    jinggang  |  -.2004033   .0617144    -3.25   0.001    -.3213613   -.0794453
##    lianqian  |  -.3522646    .092291    -3.82   0.000    -.5331518   -.1713775
##    lianxing  |  -.3310802   .1642526    -2.02   0.044    -.6530094    -.009151
##    longqing  |  -.3854579   .1559286    -2.47   0.013    -.6910725   -.0798434
##     lusikou  |   .0157262   .0545588     0.29   0.773    -.0912071    .1226595
##      mazhou  |  -.3782782   .1674892    -2.26   0.024     -.706551   -.0500054
##      minzhu  |  -.2942304    .064151    -4.59   0.000     -.419964   -.1684967
##     shigang  |  -.3364245   .1861184    -1.81   0.071    -.7012098    .0283607
##    wanshang  |  -.5649449   .1488353    -3.80   0.000    -.8566568    -.273233
##      xiabao  |  -.4471856    .211903    -2.11   0.035    -.8625078   -.0318634
##     xianghu  |  -.1679932   .1460592    -1.15   0.250    -.4542639    .1182776
##    xiaofang  |  -.2734452   .0752336    -3.63   0.000    -.4209004     -.12599
##   xiashatou  |  -.1853048   .0843979    -2.20   0.028    -.3507216   -.0198879
##    xiecheng  |  -.3877247   .1037921    -3.74   0.000    -.5911536   -.1842959
##        xihe  |  -.4328249    .076155    -5.68   0.000     -.582086   -.2835638
##        xihu  |  -.5068664   .1357542    -3.73   0.000    -.7729398    -.240793
##      xilian  |  -.1510135    .046111    -3.27   0.001    -.2413895   -.0606375
##    xinguang  |  -.4163167   .1033951    -4.03   0.000    -.6189674    -.213666
##    xingzeng  |  -.2715014   .1069633    -2.54   0.011    -.4811456   -.0618572
##     xinlian  |  -.4826528    .070666    -6.83   0.000    -.6211556     -.34415
##    yanjiang  |  -.4610531   .1743666    -2.64   0.008    -.8028055   -.1193008
##      yazhou  |  -.2897794   .0528273    -5.49   0.000     -.393319   -.1862397
##    yongfeng  |  -.5086532   .1275676    -3.99   0.000    -.7586811   -.2586253
##     zhangxi  |  -.3005388   .0938223    -3.20   0.001    -.4844272   -.1166505
##        zixi  |  -.3452332   .0904892    -3.82   0.000    -.5225889   -.1678776
##              |
##        _cons |   .0204836   .1254827     0.16   0.870     -.225458    .2664253
## ------------------------------------------------------------------------------
## Instrumented:  pre_takeup_rate
## Instruments:   male age agpop ricearea_2010 literacy intensive risk_averse
##                disaster_prob 2.villid 3.villid 4.villid 5.villid 6.villid
##                7.villid 8.villid 9.villid 10.villid 11.villid 12.villid
##                13.villid 14.villid 15.villid 16.villid 17.villid 18.villid
##                19.villid 20.villid 21.villid 22.villid 23.villid 24.villid
##                25.villid 26.villid 27.villid 28.villid 29.villid 30.villid
##                31.villid 32.villid 33.villid 34.villid 35.villid 36.villid
##                37.villid 38.villid 39.villid 40.villid 41.villid 42.villid
##                43.villid 44.villid default

这就得到了$表19.1$中的结果（要记住，这个表没有展示所有控制变量的大量系数）³⁹。第一阶段回归将内生变量（你的朋友是否购买了保险）作为因变量，并且有一个我们工具变量的系数。该系数为 $0.118$，且在统计上显著。它表明，你的朋友被分配到”默认购买”的实验条件下，会使他们购买保险的概率提高 $11.8$ 个百分点。我们利用这 $0.118$ 的增幅（以及表格中未显示的其他预测变量）来预测你的朋友是否购买了保险，并在第二阶段使用这些预测值。

第二阶段将你实际购买保险的结果作为因变量。仅利用随机实验分配和控制变量所驱动的随机变异，你朋友购买保险的比例每增加 $1$ 个单位，你购买保险的概率就会提高 $0.791$。这是一个相当强的溢出效应！

这就是两阶段最小二乘法。那广义矩方法$（GMM）$和有限信息极大似然法$（LIML）$呢？在 $Stata$ 和 $Python$ 中，转换很容易。在 $Stata$ 中，运行 $ivregress \space gmm$或 $ivregress \space liml$，而不是 $ivregress \space 2sls$。不幸的是，有一些重要的 $LIML$ 参数 $\alpha$ 和 $\kappa$，我们之后会讨论，$ivregress \space liml$ 不允许你自行设置这些参数。在 $Python$ 中，你可以使用 $IVGMM$ 或 $IVLIML$ 来代替 $IV2SLS$，而且 $IVLIML$ 确实允许你控制 $\alpha$ 和 $\kappa$。

那在 $R$ 中呢？在 $R$ 中，不幸的是，至少目前我们得切换包。$ivmodel()$ 函数能够进行 $LIML$ 分析，并且有设置 $\alpha$ 和 $\kappa$ 的选项。那 $GMM$ 呢？目前，你得自己使用 $gmm$ 包中的 $gmm()$ 函数来搭建整个双方程 $GMM$ 模型，尽管你对标准误调整的选项会有点受限。还有 $momentfit$ 包，它有一些不错的改进，但同样也很难搭建。这个例子的语法是：

library(modelsummary); library(gmm)

## Loading required package: sandwich

d <- causaldata::social_insure
# Remove all missing observations ourselves
d <- d %>%
    select(takeup_survey, male, age, agpop, ricearea_2010,literacy, intensive, risk_averse, disaster_prob,
           village, address, pre_takeup_rate, default) %>%
    na.omit()
m <- gmm(takeup_survey ~ male+age+agpop+ricearea_2010+literacy+intensive+risk_averse+disaster_prob+factor(village)+pre_takeup_rate,
         ~ male+age+agpop+ricearea_2010+literacy+intensive+risk_averse+disaster_prob+factor(village)+default, data = d)

# We can apply the address clustering most easily in msummary
msummary(m, vcov = ~address, stars = c('*' = .1, '**' = .05, '***' = .01))

	(1)
* p < 0.1, p < 0.05, * p < 0.01
(Intercept)	0.020
	(0.129)
male	0.038
	(0.047)
age	0.005***
	(0.001)
agpop	-0.002
	(0.006)
ricearea_2010	0.003***
	(0.001)
literacy	0.082**
	(0.036)
intensive	-0.009
	(0.026)
risk_averse	0.237***
	(0.044)
disaster_prob	0.001*
	(0.001)
factor(village)beixing	-0.493***
	(0.112)
factor(village)caijia	-0.393***
	(0.128)
factor(village)daqiao	-0.494**
	(0.219)
factor(village)daxi	-0.271***
	(0.102)
factor(village)dayu	-0.329***
	(0.111)
factor(village)dazhou	-0.506***
	(0.100)
factor(village)dongan	-0.428***
	(0.152)
factor(village)dukou	-0.499***
	(0.111)
factor(village)fusheng	-0.439***
	(0.139)
factor(village)fuzhou	-0.438***
	(0.168)
factor(village)gangtou	-0.417***
	(0.117)
factor(village)gangxia	-0.213
	(0.166)
factor(village)guojia	-0.109
	(0.124)
factor(village)hefeng	-0.368***
	(0.091)
factor(village)helin	-0.367***
	(0.137)
factor(village)hongxing	-0.342***
	(0.109)
factor(village)huangshan	-0.394***
	(0.148)
factor(village)jingang	-0.465***
	(0.101)
factor(village)jinggang	-0.200**
	(0.092)
factor(village)lianqian	-0.352***
	(0.125)
factor(village)lianxing	-0.331**
	(0.168)
factor(village)longqing	-0.385**
	(0.166)
factor(village)lusikou	0.016
	(0.173)
factor(village)mazhou	-0.378**
	(0.176)
factor(village)minzhu	-0.294***
	(0.095)
factor(village)shigang	-0.336**
	(0.142)
factor(village)wanshang	-0.565***
	(0.185)
factor(village)xiabao	-0.447
	(0.278)
factor(village)xianghu	-0.168
	(0.173)
factor(village)xiaofang	-0.273**
	(0.138)
factor(village)xiashatou	-0.185
	(0.139)
factor(village)xiecheng	-0.388***
	(0.131)
factor(village)xihe	-0.433***
	(0.094)
factor(village)xihu	-0.507***
	(0.177)
factor(village)xilian	-0.151
	(0.114)
factor(village)xinguang	-0.416***
	(0.111)
factor(village)xingzeng	-0.272**
	(0.138)
factor(village)xinlian	-0.483***
	(0.084)
factor(village)yanjiang	-0.461**
	(0.185)
factor(village)yazhou	-0.290***
	(0.107)
factor(village)yongfeng	-0.509***
	(0.140)
factor(village)zhangxi	-0.301***
	(0.099)
factor(village)zixi	-0.345***
	(0.107)
pre_takeup_rate	0.791***
	(0.251)
Num.Obs.	1378
Std.Errors	by: address

最后，如果我们的工具变量模型中有很多固定效应，那该怎么办呢？在这些情况下，必须进行一些技术调整。这一次，$R$ 的转换是最容易的。我们已经使用过的$feols()$函数可以轻松地纳入固定效应 ——$+factor(village)$只需变成 $village$ 就行。不过，$Stata$ 和 $Python$ 也不算太难。在 $Stata$ 中，你要从 $ivregress$ 切换到 $xtivreg$。这会带来一些其他变化 —— 你必须用 $xtset$ 来告诉 $Stata$ 面板结构是什么样的，它只进行两阶段最小二乘法分析，而不是广义矩方法或有限信息极大似然法，而且你必须告诉它你想要固定效应$（fe）$、随机效应$（re）$还是其他什么（查看帮助文件）。目前，还没有与之对应的 $Python$ 方法，不过我们常用的 $Python$ 语法可以通过添加二进制指示变量集，让你估计包含固定效应的 $IV$ 模型。

19.2.2 工具变量与处理效应

工具变量估计的是什么？当思考平均工具变量会产生何种处理效应时，我们可以回顾第 $10$ 章的内容。

总的来说，我们知道工具变量是一种专门用于分离由工具变量解释的处理变量变异的方法。这意味着我们关注的是局部平均处理效应，其中个体处理效应的权重取决于该个体观测对工具变量的响应程度。

在像两阶段最小二乘法或广义矩方法这类标准估计量的情况下，当只有一个处理 / 内生变量和一个工具变量时，权重就是在第一阶段中工具变量对你产生的个体效应大小。

例如，假设最近有一系列电视广告鼓励人们更多地锻炼。你想把对这些广告的接触情况作为你锻炼量的工具变量，然后考察锻炼对血压的影响。

考虑样本中的三个人：贾凯拉、凯尔和李。这些广告会让贾凯拉每周额外锻炼半小时，而每周额外锻炼一小时会使她的血压降低 $2$ 个点。你可以在 $表19.2$ 中看到这些数值，以及凯尔和李的数值。

要记住 —— 广告对锻炼时长的那些影响，是理论上对那些个体而言会产生的效应。显然，我们无法同时看到贾凯拉既接触广告又不接触广告的情况。但这是在说，接触广告的贾凯拉每周比不接触广告的贾凯拉多锻炼半小时。

两阶段最小二乘法会告诉我们锻炼时长对血压的效应是多少呢？嗯，贾凯拉对广告的反应最强烈，所以她所获得的锻炼带来的 $-2$ 的效应会被赋予更大的权重。具体来说，它采用了她在广告效应上的 $0.5$ 的权重。类似地，凯尔得到的权重是 $0.25$。另一方面，李对广告完全没有反应 —— 广告对他没有任何影响。所以结果是，他对 $2SLS$ 的估计没有任何影响，他的权重为 $0$。

当我们进行 $2SLS$ 估计时，它无法看到这些理论上的个体效应。尽管如此，它给我们的局部平均处理效应$（LATE）$是基于这些效应的！具体来说，$LATE$ 是：

\[ \frac{0.5 \times (-2) + 0.25 \times (-8) + 0 \times (-10)}{0.5 + 0.25 + 0} = \frac{-1 + -2}{0.75} = -4 \quad (19.3) \]

这就是 $2SLS$ 会给出的结果。这与平均处理效应$(-2) + (-8) + (-10))/3 = -6.67$形成了对比。

这立刻让我们得出一个非常重要的结论：$2SLS$ 会根据所使用的工具变量的不同而给出不同的结果。如果我们选择了一个对让李锻炼非常有效，但对贾凯拉不那么有效的工具变量，那么 $2SLS$ 会估计出锻炼对血压更强的效应。

如果有不止一个工具变量呢？同样的逻辑适用 —— 工具变量对你的效应越强，你的处理效应被赋予的权重就越大。但在这一点上，数学变得有点复杂，因为它会是你从不同工具变量那里得到的不同权重的混合。所以，如果一个工具变量给贾凯拉 $0.5$ 的权重，而另一个给她 $0.8$ 的权重，那么同时使用这两个工具变量时她所得到的权重会是 $0.5$ 和 $0.8$ 的某种混合（而且不一定只是两者的平均值）。

这里要强调的一点是，我们得到的特定权重取决于估计量。不同的工具变量估计方法会产生不同的加权平均处理效应。如果你没有使用两阶段最小二乘法，一定要去弄清楚你到底得到了什么。尽管在一些特殊情况之外，一般来说，你受到工具变量的影响越大，被赋予的权重通常也越大。我们在估计局部平均处理效应的想法在大多数情况下是成立的，尽管我们得到的具体是哪一种局部平均处理效应会因估计量的不同而改变。

在思考这些权重，尤其是工具变量中的”$O$“时，有一个常用的术语体系。我们可以把样本分成三组：

依从者：对于依从者来说，工具变量对处理的效应符合预期方向。在 $表19.2$ 中，贾凯拉和凯尔是依从者，因为告诉他们要多锻炼的广告让他们更多地锻炼了。
总是接受者 / 从不接受者：总是接受者 / 从不接受者完全不受工具变量的影响。在 $表19.2$ 中，李是总是接受者 / 从不接受者，因为广告对他的锻炼量的效应为 $0$⁴⁰。
违背者：违背者受到工具变量的影响，但其方向与预期方向相反⁴¹。

从这个术语体系中，我们可以得到一个结果和一个需要做出的假设：

首先，结果：如果所有依从者受到工具变量的影响程度相同，那么 $2SLS$ 给出的是依从者中的平均处理效应……⁴²不错！尽管工具变量不太可能以完全相同的方式影响每一个人。

其次，假设：为了让这一切成立，我们需要假设不存在违背者。想象一下，我们在 $表19.2$ 中添加第四个人：梅西，她非常讨厌这些广告，以至于如果看到广告，她决定少锻炼 $0.25$ 小时。每周锻炼一小时会使她的血压降低 $8$ 个点。

现在局部平均处理效应的权重是怎么计算的呢？梅西的权重是 $-0.25$。现在加权平均的计算结果为：

\[ \frac{0.5 \times (-2) + 0.25 \times (-8) + 0 \times (-10) + (-0.25) \times (-4)}{0.5 + 0.25 + 0 + (-0.25)} = \frac{-1 + (-2) + 2}{0.5} = -2 \quad (19.4) \]

实际上，锻炼对梅西的效果比我们之前估计的效果更显著，但加入她之后，整体效果反而变小了！这是因为她的权重为负，这使得加权平均的计算变得很奇怪，所以它们不再是真正的加权平均了⁴³。

不存在违背者的假设也被称为单调性假设。和有效性、相关性假设一样，这是我们在使用工具变量时需要做出的另一个关键假设，尽管这个假设往往受到的关注要少得多。

单调性：在工具变量的背景下，假设数据中的每一个人，要么工具变量对他们的影响方向相同（正向或负向），要么完全没有影响（效应为零）。

所以，如果你有一个平均来看有效果的工具变量，要仔细思考这种效果是否可能对每个人都朝着相同的方向产生影响。有很多情况并非如此，比如，有些人可能会因为一项干预措施而非常恼火，以至于它会产生与预期相反的效果。

或者，更广泛地说，人们是不同的，会有不同的反应。安格里斯特和埃文斯$（1998年）$的研究就是思考这个问题的一个很好的例子。他们观察到，家庭似乎更倾向于既有男孩又有女孩。一个恰好头两个孩子都是男孩，或者都是女孩的家庭，更有可能要第三个孩子，试图凑成一男一女。所以，“你的头两个孩子性别相同” 在很多研究中被用作”要第三个孩子”的工具变量，这一做法是受到安格里斯特和埃文斯研究的启发。

Angrist, Joshua D., and William N. Evans. 1998. “Children and Their Parents’ Labor Supply: Evidence from Exogenous Variation in Family Size.” American Economic Review 88 (3): 450–77.

但要让这种工具变量的应用奏效，必须不存在违背者。即便大多数人在头两个孩子性别相同时更有可能要第三个孩子（或者完全不根据头两个孩子的性别来决定是否要第三个孩子），但如果有一些人因为头两个孩子性别相同而更不可能要第三个孩子，那么单调性假设就会被违背。也许有些父母会因为有可能有三个同性别的孩子而感到恐惧？不管我们编造出什么理由，人们选择要更多孩子有很多复杂的原因，其中一些可能与我们需要的假设情况相矛盾。所以，使用这种工具变量的研究需要仔细思考单调性假设是否可能成立，以及他们能对此做些什么。

19.2.3 检验工具变量假设

为了让工具变量发挥我们期望的作用，我们依赖于相关性假设。

记住，相关性假设是指工具变量 $Z$ 和处理 / 内生变量 $X$ 之间存在关联。我们可以进一步说，我们需要假设 $X$ 和 $Z$ 之间的关联足够强，这样我们就不会遇到”弱工具变量问题”。

弱工具变量：仅能对处理 / 内生变量进行弱预测的工具变量。

弱工具变量是指有效且确实能预测处理变量，但它对处理变量的预测作用很小，预测能力很弱。记住我们的一般直觉，即 $IV$ 给出的是 $Cov(Z,Y)/Cov(Z,X)$，那么如果 $Cov(Z,X)$ 很小，我们就接近除以 $0$ 的问题了！整个估计值会变得非常大（因为你在除以一个很小的数），而且抽样变异会变得极大。

值得庆幸的是，与我们涉及的很多假设不同，弱工具变量（进而相关性）是相当容易检验的⁴⁴。毕竟，我们只需要确保 $Z$ 和 $X$ 是相关的，而且不是弱相关。我们知道如何测量和检验关联！

到目前为止，检验相关性最常见的方法是第一阶段$F-统计量检验$⁴⁵。方便的是，它也是最容易操作的。你所要做的就是：

估计模型的第一阶段（将处理 / 内生变量对控制变量和工具变量进行回归）
对工具变量进行$联合F检验$
从该检验中得到$F统计量$，并利用它来判断工具变量是否具有相关性

就是这样！如果有不止一个处理 / 内生变量，计算会变得有点复杂（因为不存在真正意义上的单一第一阶段），但核心思想是一样的。

我们在第 $13$ 章关于回归的内容里讲过做联合$F$检验的代码，但因为这是个很常见的检验，你的工具变量命令通常会自动为你做这件事。在 $R$ 中，如果你的 $feols$回归被保存为 $m$，$summary(m)$ 会报告”第一阶段$F$统计量”⁴⁶。在 $Stata$ 中，你可以在 $ivregress$ 命令之后使用 $estat \space firststage$ ⁴⁷。在 $Python$ 中，用 $statsmodels$ 的话，你可以做 $IV2SLS().fit().first\_stage$ ，然后查看”$partial \space F-stat$“。—— 我们不只是担心这种关系为 $0$，我们还担心它很小。

因为小没有一个单一精确的定义，所以也没有一个单一正确的$F$统计量临界值可供参考。相反，我们面临一种权衡。你的$F$统计量越大，你得到的偏差就越小。所以你想要的$F$统计量会基于你愿意接受多少偏差。

等等 —— 什么偏差？谁提到偏差了？弱工具变量会导致偏差，因为即使工具变量确实有效，在实际的数据样本中，工具变量也会仅仅由于随机因素而与误差项存在非$0$关系，这会损害有效性并给你带来偏差。工具变量越弱，这种情况就越严重。

我们可以围绕这种权衡来确定$F$统计量临界值的选择。斯托克和约戈$（2005年）$计算了使用工具变量时相对于单独使用普通最小二乘法所产生的偏差。工具变量越强，$IV$ 偏差相对于 $OLS$ 偏差就越小。对于特别弱的工具变量，你还不如回到 $OLS$（或者，如果 $OLS$ 也不行，那你可能就没什么办法了）。他们论文末尾的表格会告诉你，例如，如果你有一个处理 / 内生变量和三个工具变量，你需要 $F$ 统计量超过 $13.91$ 才能将 $2SLS$ 中的 $IV$ 偏差降低到小于 $OLS$ 偏差的 $5\%$，但只需要 $F$ 统计量超过 $9.08$ 就能将 $IV$ 偏差降低到小于 $OLS$ 偏差的 $10\%$。如果你有四个工具变量呢？那你需要 $16.85$ 或 $10.27$。一些第一阶段检验命令（比如 $Stata$ 的）会自动告诉你相关的临界值⁴⁸。

Stock, James H., and Motohiro Yogo. 2005. “Testing for Weak Instruments in Linear IV Regression.” In Identification and Inference for Econometric Models: Essays in Honor of Thomas Rothenberg, 80–108. Cambridge: Cambridge University Press.

如果你的$F$统计量不达标，该怎么办呢？在这一点上，你可以直接放弃 —— 不少研究者会建议放弃那些第一阶段$F$统计量小的项目。但如果每个人都这么做，这种做法就不好了，尤其是当你谈论的是可能最终会被发表的项目时。如果基于$F$统计量临界值舍弃了低$F$统计量的项目，那么最终被发表的项目就会是真正强的工具变量项目和研究者恰好得到一个$F$统计量大的样本的实际弱工具变量项目的混合。当实际弱的工具变量碰巧产生了大的$F$统计量时，这就会带来我们之前担心的那种弱工具变量偏差。我们是特意基于随机得到工具变量无效的样本来挑选要完成的项目。所以，$F$统计量临界值会导致已发表文献中存在偏差的结果比例更高。

那该怎么办呢？在弱工具变量的情况下，只要你确实非常确定工具变量是有效的，与其放弃项目，你可能想要尝试一种对弱工具变量更稳健的工具变量方法。其中一些方法，比如安德森 - 鲁宾置信区间，会在”专业人士怎么做”部分展示。

事实上，无论如何，你可能只想直接使用”专业人士怎么做”部分中那些对弱工具变量稳健的方法。弱工具变量的问题在于，它们会引入大量的抽样变异，并干扰你的标准误。预检验的思路是，你只在知道问题不太大的情况下进行分析。但这并不意味着问题就消失了。李等人$（2022年）$表明，即使第一阶段$F$统计量达到 $104.7$，你的置信区间至少也会有一点点错误。基恩和尼尔$（2022年）$表明，典型的临界值忽略了统计功效之类的因素，而且在非常现实的场景中，你可能至少需要$F$统计量达到 $50$。哎呀！这并不意味着要舍弃任何$F$值低于 $104.7$ 或 $50$ 之类的项目。相反，这只是意味着预检验可能无法解决真正的问题。

ee, David S., Justin McCrary, Marcelo J Moreira, and Jack Porter. 2022. “Valid t-Ratio Inference for IV.” American Economic Review 112 (10): 3260–90. Keane, Michael P., and Timothy Neal. 2022. “A Practical Guide to Weak Instruments.” Annual Review of Economics 16.

所以我们可以检验相关性和弱工具变量。那有效性呢？有几个检验我或许应该提一下，尽管这些检验让我很不爽。特别是，有几个非常知名的检验是专门用来检验有效性假设的。虽然这些检验有其用途，但它们的应用方式我觉得完全没有说服力，原因我之后会提到。不过，它们非常常见，所以至少要足够了解它们，知道是怎么回事。

首先，检验有效性到底确切意味着什么呢？我们想要检验工具变量 $Z$ 和结果 $Y$ 之间是否存在任何开放的后门路径。我们可以将其重新表述为检查 $Z$ 是否与第二阶段的误差项 $\varepsilon$ 相关。我们希望 $Z$ 和 $\varepsilon$ 不相关。如果它们相关，有效性就会被违背。

检验有效性存在一些明显需要克服的障碍。首先，我们实际上无法观测到误差项 $\varepsilon$，所以不能直接去看 $Z$ 和 $\varepsilon$ 之间的关系。

那为什么不直接看 $Z$ 和残差 $r$ 之间的关系呢？我们一直用 $r$ 来代替 $\varepsilon$ 计算诸如标准误之类的东西。然而，如果 $Z$ 是无效的，那么第二阶段的估计就会有偏差，这会使得残差不能很好地代表误差。因为误差是真实模型预测值和结果之间的差异，而我们知道模型是有偏差的，所以估计出的模型平均来说甚至都不是真实模型，因此残差不能很好地代表误差。

那我们能做什么呢？

有一种你经常能看到，但实际上毫无意义的方法，就是运行第二阶段模型，但把工具变量作为控制变量包含进去。

\[ Y = \beta_0 + \beta_1 X + \beta_2 Z + \varepsilon \quad (19.5) \]

如果工具变量 $Z$ 的系数非零，这表明有效性存在违背。为什么呢？因为所有从 $Z$ 到 $Y$ 的开放路径都应该通过 $X$，对吧？所以如果我们在控制 $X$ 的情况下观察 $Z$ 对 $Y$ 的效应，这应该关闭所有路径。如果在 $\hat{\beta}_2$ 中仍然存在明显的关系，那么似乎还有其他破坏有效性的路径仍然是开放的！我一直说”这表明似乎”，是因为它实际上什么都没告诉你！为什么呢？因为存在对撞变量偏差！如果 $Z$ 和 $\varepsilon$ 都导致 $X$，那么即使 $Z$ 和 $\varepsilon$ 不相关，当你控制 $X$ 时（我们在这里就是这么做的），它们会突然又变得相关了。所以我们得到的 $\beta_2$ 的估计是有偏差的，而且这个估计并没有告诉我们 $Z$ 对 $Y$ 的直接效应的任何信息。

另一种方法是使用 $Durbin-Wu-Hausman \space test$。$Durbin-Wu-Hausman \space test$ 会比较两个模型的结果，其中一个模型如果某些假设不成立，结果可能不一致，而另一个模型不依赖那些假设，但精确性较低。如果结果相似，这表明假设至少没有错到会搞乱结果的程度，也意味着可以使用更精确的模型。

在工具变量的背景下，$Durbin-Wu-Hausman \space test$ 有两种用法。首先，它可以用来比较普通最小二乘法（若 $X$ 与 $\varepsilon$ 相关则不一致）和工具变量法（精确性较低）。如果结果不同，这意味着 $X$ 确实存在开放的后门路径，我们可能应该使用工具变量法。

其次，$Durbin-Wu-Hausman \space test$ 可以在过度识别检验中用来比较两个不同的工具变量模型。如果我们拥有的工具变量数量超过所需，处于过度识别状态，那么只要我们非常确定拥有的有效工具变量数量达到所需，就可以将使用所有工具变量的工具变量模型（若其中一些工具变量无效则不一致）与仅使用我们非常确定的工具变量的工具变量模型（精确性较低）进行比较。如果结果不同，这会告诉我们额外的工具变量可能是无效的。

实际上，对于两阶段最小二乘法，更常使用 $Sargan \space test$ 来进行过度识别检验，该检验会从 $2SLS$ 模型的第二阶段获取残差，并查看这些残差与工具变量之间的关系。对于广义矩方法模型，通常会使用 $Hansen \space test$ ，这里我就不深入讲解了。对于所有这些检验，都有很多不同的操作方式。

在 $R$ 中，当你对从 $feols$ 得到的工具变量模型进行 $summary()$ 时，会自动得到 $Durbin-Wu-Hausman \space test$ 和 $Sargan \space test$。如果你使用的是其他命令，可能需要使用 $diagnostics=TRUE$ 选项。在 $Stata$ 中，你可以在回归后使用 $estat \space endog$ 来检验 $OLS$ 与 $IV$，或者使用 $estat \space overid$ 来进行萨根检验。在 $Python$ 中，你可以通过模型的 $wu\_hausman()$ 得到 $Durbin-Wu-Hausman \space test$，或者通过 $sargan()$ 得到 $Sargen \space test$。

所以这些就是我们用来检验工具变量有效性的主要工具。为什么这些检验让我不爽呢？因为人们试图用它们来为有问题的工具变量辩护。想象一下，有人向你展示一项研究，其中的工具变量你认为是无效的，或者一个普通最小二乘法模型，其中的处理变量你认为仍存在一些开放的后门路径。但是啊 —— 别担心！他们做了内生性检验，并且没能拒绝有问题的处理 / 工具变量是有效的这一原假设。没问题，对吧？

还是有问题！原假设是否被拒绝，除了取决于你所检验的根本情况是否为真之外，还涉及很多其他因素 —— 比如统计功效、抽样变异等等。没能发现有效性违背的证据，可能意味着工具变量是有效的，也可能意味着是其他那些因素在起作用。

这些检验更合理的用途是，当你有一个你非常确定的工具变量，但担心自己可能只是运气不好抽到了不好的样本时。用这些检验来对一个你原本非常确定的工具变量感到失望，比用它们来安抚自己对一个不确定的工具变量的担忧，要更有道理得多⁴⁹。

Mogstad, Magne, Alexander Torgovitsky, and Christopher R. Walters. 2021. “The Causal Interpretation of Two-Stage Least Squares with Multiple Instrumental Variables.” American Economic Review 111 (11): 3663–98.

这些检验，尤其是过度识别检验，让我不爽还有另一个原因。这和工具变量模型所产生的局部平均处理效应估计有关。当你有不止一个工具变量时，这些工具变量会产生不同的局部平均处理效应估计，因为它们解释了处理变量变异的不同部分。当你把两个工具变量结合起来时，你会得到另一个局部平均处理效应⁵⁰。所以，如果你的过度识别检验发现模型在使用不同工具变量时产生了不同的结果，嗯…… 这不一定意味着工具变量是无效的。这只是意味着这些工具变量产生的结果不一样（帕伦特和席尔瓦， $2012年$）。

Parente, Paulo M. D. C., and J. M. C. Santos Silva. 2012. “A Cautionary Note on Tests of Overidentifying Restrictions.” Economics Letters 115 (2): 314–17.

19.3 专业人士怎么做

19.3.1 别只检验工具变量是否弱，要解决它！

弱工具变量是个实实在在的问题。我已经讨论过用$F检验$来检测弱工具变量的方法，但这种检验本身也存在问题。相反，我们可以直接采用一种受弱工具变量影响没那么大的工具变量估计方法，这类方法确实是存在的！

第一种解决办法也是最容易的。在只有一个处理 / 内生变量和一个工具变量的情况下⁵¹，我们可以调整标准误，以考虑工具变量可能较弱的情况。

这个解决办法很古老，早在 $1949$ 年就有了。$Anderson-Rubin \space confidence \space intervals$ 即使在工具变量较弱的情况下，也能为我们对效应的估计提供有效的不确定性度量（安德森和鲁宾， $1949$）⁵²。这其实并没有真正解决工具变量弱的问题。也就是说，你仍然面临弱工具变量以及随之而来的抽样变异问题。但它确实能确保你的结果反映出这种变异。换句话说，在你面临弱工具变量问题时，它会让结果更”诚实”。

Anderson, Theodore W., and Herman Rubin. 1949. “Estimation of the Parameters of a Single Equation in a Complete System of Stochastic Equations.” The Annals of Mathematical Statistics 20 (1): 46–63.

就是这么回事。如果你有一个工具变量和一个处理 / 内生变量，不管你的工具变量是否弱，你都不妨跳过检验弱工具变量的整个过程，直接报告$Anderson-Rubin \space confidence \space intervals$ 。虽然这有点令人意外，但这并不是一种完全通用的做法。正如你在下一段会看到的，在常见的工具变量命令中，它甚至都没有被作为一个选项。我们已经知道这个方法很多年了，唉，好吧。

在 $R$ 中，你可以先使用 $AER$ 包中的 $ivreg$ 函数估计你的工具变量模型，然后将其传递给 $ivpack$ 包中的 $anderson.rubin.ci$ 函数，以此得到$Anderson-Rubin \space confidence \space intervals$。在 $Stata$ 中，你可以在 $ivregress$ 模型之后，使用 $weakiv$ 包中的 $weakiv$ 函数来得到$Anderson-Rubin \space confidence \space intervals$。据我所知，在 $Python$ 中没有现成的方法可以得到$Anderson-Rubin \space confidence \space intervals$。

一种较新的方法，“有效 $t$ 比率推断”，出自李等人（$2022年$）的研究。这种方法沿用了更标准的两阶段最小二乘法统计量类型，但会调整标准误和置信区间，以考虑第一阶段估计以及弱工具变量可能带来的额外噪声。有效 $t$ 比率推断可以在 $R$ 语言的 $ivDiag$ 包中实现。在 $Stata$ 中，可以使用$net \space install \space tf, force \space from (http://www.princeton.edu/~davidlee/wp)$ 命令安装 $tf$包。

Lee, David S., Justin McCrary, Marcelo J Moreira, and Jack Porter. 2022. “Valid t-Ratio Inference for IV.” American Economic Review 112 (10): 3260–90.

我们可以更诚实地报告弱工具变量问题，但我们真的能解决它吗？也许可以！

当工具变量较弱时，一种试图表现更优的常用估计方法是有限信息极大似然法$（LIML）$⁵³。

简而言之，工具变量法会利用处理 / 内生变量 $X$ 和结果 $Y$ 中由工具变量 $Z$ 预测的部分。在工具变量的背景下，$LIML$ 做的是同样的事，只是它会用一个参数 $\kappa$ 来缩小预测值，使其变弱，而 $\kappa$ 通常是根据数据在模型中估计得到的。如果 $\kappa = 1$，我们将预测值按 $1$ 的比例缩放，不做任何改变，最终又回到 $2SLS$。但如果 $\kappa < 1$，我们会找回一点原本在普通最小二乘法中会用到的内生变量，从而减少对那个弱工具变量的依赖。

我们也可以使用被称为 $Fuller$ 的 $\alpha$ 参数对 $\kappa$ 进行调整。不再使用基于数据估计得到的 $\hat{\kappa}$，$\alpha$ 会让 $\kappa$ 变得更小，改用 $\kappa = \hat{\kappa} - \alpha/(N - N_I)$，其中 $N$ 是观测数量，$N_I$ 是工具变量数量。例如，$\alpha = 4$ 能使估计在不同样本中的均方误差最小化。它也可以根据数据进行估计。$LIML$ 过程减少了弱工具变量带来的偏差，但这可能 —— 不出所料 —— 以降低精度为代价（布卢姆奎斯特和达尔伯格，$1999年$）。

Blomquist, Sören, and Matz Dahlberg. 1999. “Small Sample Properties of LIML and Jackknife IV Estimators: Experiments with Weak Instruments.” Journal of Applied Econometrics 14 (1): 69–88.

还有别的方法吗？最有前景的研究方向之一是关注大量弱工具变量的使用。当然，单独一个弱工具变量对处理变量的影响不大。但一大群弱工具变量呢？这就像是电影的高潮部分。我不可能在这里把整个相关文献都讲清楚，但一个起始点可以看赵和斯旺森（$2005年$）的研究。

Chao, John C., and Norman R. Swanson. 2005. “Consistent Estimation with a Large Number of Weak Instruments.” Econometrica 73 (5): 1673–92.

19.3.2 超越局部平均处理效应

我们是不是只能局限于局部平均处理效应呢？完全不是！局部平均处理效应是两阶段最小二乘法默认给出的结果，但还有其他方法。关于工具变量估计的各种方法，有大量文献旨在得出平均处理效应或其他平均值，而非局部平均处理效应。当然，这些方法都依赖不同的假设，或者适用于不同的情境。

这里我仅简要提及两种。这只是代表两种不同的获取平均处理效应通用方法的两篇论文，而且这两种方法远非仅有的两种⁵⁴。它们只是起点。

首先，赫克曼和维特拉奇$（1999年）$研究了取大量不同值的工具变量。他们意识到，通过比较不同的值可以识别出一系列效应 —— 那些被工具变量从最低值推到略高于该值的群体，其处理效应是多少？从略高于该值推到再略高一点的值的群体，效应又是多少？依此类推。这样一来，针对工具变量的取值范围，你就有了一大堆效应。在工具变量的某些取值下处理概率为 $0$，而在另一些取值下处理概率为 $1$ 的条件下，你可以对这整个效应分布取平均，从而得到平均处理效应。

Heckman, James J., and Edward J. Vytlacil. 1999. “Local Instrumental Variables and Latent Variable Models for Identifying and Bounding Treatment Effects.” Proceedings of the National Academy of Sciences 96 (8): 4730–34.

其次是结合相关随机效应（还记得第 $16$ 章讲的相关随机效应吗？）与 $2SLS$ 的文献。另一种思考异质处理效应的方式是，处理变量与一些我们无法观测到的个体差异存在交互作用。相关随机效应试图对未观测到的个体差异进行建模，这似乎是个不错的起点。

通过让处理变量与一组控制变量交互，然后用相关随机效应估计模型，你可以得出平均处理效应。或者，甚至可以得到你喜欢的任何群体的条件平均处理效应，该群体由你给模型的控制变量来定义。伍德里奇$（2008年）$对这种方法以及它在何时效果好、何时效果不好进行了通俗易懂的概述。一个重要的要点是，这种方法往往在处理变量连续性越强时效果越好，而对于二元处理变量，情况就变得相当棘手了。

Wooldridge, Jeffrey M. 2008. “Instrumental Variables Estimation of the Average Treatment Effect in the Correlated Random Coefficient Model.” In Modelling and Evaluating Treatment Effects in Econometrics. Emerald Group Publishing Limited.

19.3.3 非线性工具变量

我们可以用工具变量和非线性模型做些什么呢？到目前为止，我们所讨论的工具变量，尤其是两阶段最小二乘法，都假定我们处于线性模型的范畴内。但如果我们有二元结果变量或处理 / 内生变量呢？我们难道不想用概率单位$（probit）$或对数单位$（logit）$模型吗？

在常见的存在二元处理 / 内生变量的情况下，一种明显的诱惑是，用概率单位或对数单位模型来运行$2SLS$的”第一阶段”，获取预测值，然后用这些预测值来估计第二阶段。如果结果变量是二元的，或许在第二阶段也用概率单位或对数单位模型。计量经济学家 $Jreey \space Hausman$ 将这种方法称为”禁忌回归”，因为它根本无法如你所愿地发挥作用（伍德里奇，$2010年$）。它不起作用的原因与非线性回归中每个变量的效应依赖于其他变量的值这一事实有关。所以，尽管工具变量 $Z$ 和第二阶段的误差项 $\varepsilon$ 不相关，但拟合值 $\hat{X}$ 不再能从 $Z$ 那里继承那种良好的不相关性了。后门路径仍然存在。

Wooldridge, Jeffrey M. 2010. Econometric Analysis of Cross Section and Panel Data. Cambridge, MA: MIT Press.

另一种诱惑是，在了解到不能随意将概率单位 / 对数单位模型和 $2SLS$ 混用时，干脆运行一个线性模型，也就是线性概率模型，忽略非线性问题⁵⁵。这并非最糟糕的想法 ——$2SLS$确实可以这样运作，不像”禁忌回归”那样不行。但另一方面，与概率单位 / 对数单位模型相比，你确实会遇到线性概率模型的所有缺点。别忘了，当二元变量的均值接近 $0$ 或 $1$ 时，普通最小二乘法相较于概率单位 / 对数单位模型，表现很差。我们可能感兴趣的不少二元处理变量，要么非常罕见，要么非常普遍。而且，这种方法得出的估计值，精度不如那些适当考虑了非线性的模型⁵⁶。

那接下来该怎么办呢？此时我们有很多选择，这取决于你想在哪个阶段（第一阶段或第二阶段）对非线性进行建模。在本节中，我们会保持简单，只关注二元处理 / 结果变量，而不是我们可能担心的无数种非线性建模的替代形式。我也不会深入讲解每一种估计这些模型的方法，因为方法实在太多了。

咱们从二元处理变量开始讲。有一个便捷且易于实施的方法，能把二元处理变量纳入两阶段最小二乘法，这是由 $Jeffery \space Wooldridge（2010年）$推广开的。

Wooldridge, Jeffrey M. 2010. Econometric Analysis of Cross Section and Panel Data. Cambridge, MA: MIT Press.

首先，通过对工具变量和控制变量进行处理 / 内生变量的非线性回归，来估计第一阶段。运行那个概率单位模型！然后，获取预测值。

接着，不像”禁忌回归”那样把这些预测值直接代入第二阶段，而是在 $2SLS$ 中用这些预测值来替代工具变量。

就这么简单！在这个过程中，非线性不再会对估计产生偏差，而且和运行线性概率模型相比，我们能得到更精确的估计。

你或许更倾向于”处理效应回归”方法，它完全避开$2SLS$，直接对数据实际的二元结构进行建模。这些方法基本上是同时估计概率单位第$1$阶段和线性第$2$阶段⁵⁷，使得工具变量能像往常一样对处理 / 内生变量产生影响。要估计处理效应回归，你可以去看看 $R$ 语言中 $sampleSelection$ 包的 $treatReg$，或者 $Stata$ 里的 $etregress$ 函数。这些函数也有工具能帮着估计除局部平均处理效应之外的处理效应平均值。至少目前，在 $Python$ 里你得自己去做。

要是结果变量是二元的呢？在这种情况下，一种常见的解决办法是采用工具变量的控制函数法。

控制函数法和 $2SLS$ 很像，不同之处在于，它不是分离出 $X$ 中被解释的部分并在第二阶段使用，而是使用常规的 $X$，但同时控制 $X$ 中未被解释的部分。在常规的线性工具变量中，$2SLS$ 和控制函数法给出的点估计是一样的⁵⁸。

然而，当涉及非线性工具变量时，这两种方法就不再相同了。而且，当 $2SLS$ 方法不再奏效时，在概率单位第二阶段加入线性第一阶段残差作为控制变量的控制函数法，并不会失效。

话虽如此，和两阶段最小二乘法一样，要对二元结果变量正确使用控制函数法，必须进行一些调整。所以一般来说，你会想要使用专为这项工作设计的估计函数。在 $R$ 语言中，这可以在 $ivprobit$ 包中实现。在 $Stata$ 中，这是带有 $twostep$ 选项的 $ivprobit$⁵⁹。据我所知，如果你想在 $Python$ 中做这件事，就得自己动手。研究 $R$ 语言 $ivprobit$ 包的代码是开展这个项目的一个好起点。

哦，不好。要是两者都是二元的呢？在这些情况下，通常会应用二元概率单位模型$（bivariate \space probit \space model）$。还记得几段之前讲的处理效应回归吗？在那个模型中，我们用极大似然法同时估计概率单位第一阶段和线性第二阶段。

在二元概率单位模型中，我们用极大似然法同时估计两个概率单位模型，这两个模型的因变量是相关的。所以，如果我们把其中一个因变量作为另一个的预测变量，同时给它一个从另一个模型中排除的预测变量，那就意味着我们在同时估计概率单位第一阶段和概率单位第二阶段。这就是工具变量！处理效应回归和二元概率单位模型的技术细节不同，但概念是一样的⁶⁰。

二元概率单位模型作为 $2SLS$ 的替代方法很重要，因为当两个阶段都是二元的时候，$2SLS$ 可能会特别不精确。另外，二元概率单位模型能给你平均处理效应，而不是局部平均处理效应，这可能会很有用（奇布里斯、达斯和洛克申，$2012$）。

Chiburis, Richard C., Jishnu Das, and Michael Lokshin. 2012. “A Practical Comparison of the Bivariate Probit and Linear IV Estimators.” Economics Letters 117 (3): 762–66.

在 $R$ 语言中，可以使用 $GJRM$ 包中的 $gjrm$，通过选项 $Model="B"$ 和 $margins = c('probit', 'probit')$ 来估计二元概率单位模型。在 $Stata$ 中，你可以使用 $biprobit$ 命令。遗憾的是，在 $Python$ 中你又得自己想办法了。

19.3.4 好吧，你可以稍微违背一下有效性

我们再次扪心自问：这值得吗？我们在所有这些细节上折腾 —— 估计方法、对假设合理性的检验、面板数据下的工具变量、含非线性变量的工具变量、弱工具变量、单调性。但如果有效性不成立，这一切都毫无意义。而且，在实际随机化处理之外，我们完全有理由怀疑有效性是否成立。

幸好，有不断发展的文献探讨即使有效性并非完美成立时仍能使用工具变量的方法。

需要明确的是，这里的”有效”指的是你能从分析中得出有用的推论。你可能得提供更多种类的其他信息，或者或许你得出的推论可能不像有完全有效工具变量时那么精确。

这里有两种常用的方法。第一种是在工具变量基本有效的情况下，尽可能利用它；第二种是通过使用大量工具变量来弥补有效性的不足。

只要工具变量对处理变量有影响，我们肯定能利用它做点什么，对吧？对的！

我们能做的一件事是，不把工具变量无效当作放弃分析的理由，而是思考这种有效性违背的程度有多严重 —— 工具变量与第二阶段误差项的关联有多强 —— 并基于合理的违背程度范围，构建一个合理的估计值范围。

这种方法为我们带来了所谓的”部分识别”，因为它不是给你一个单一的估计值（精确识别出一个估计值），而是基于你可能做出的一系列不同假设，给出一个范围。我会在第 $21$ 章更详细地总体介绍这种方法，但就工具变量中的部分识别而言，一篇关键论文是康利、汉森和罗西（$2012年$）的研究。他们展示了，无需假设工具变量完全有效，你可以构建一个合理的无效性取值范围，并生成与这些取值对应的合理估计值范围。主要的有效性违背情况如何？那你的范围就会很大。但工具变量越强，你的取值范围也会越好。这种方法可以在 $R$ 语言中使用 $ivDiag$ 包实现，或者在 $Stata$ 中使用 $plausexog$ 包实现。据我所知，在 $Python$ 中没有相关实现。

Conley, Timothy G., Christian B. Hansen, and Peter E. Rossi. 2012. “Plausibly Exogenous.” Review of Economics and Statistics 94 (1): 260–72.

另一种常用方法是，用关于工具变量的其他可能比有效性更合理的假设来替代有效性假设，然后看看这个假设能让你走多远。有很多这类论文，每篇都依赖于对工具变量的不同假设。也许有一篇论文使用的假设对你的工具变量来说非常合理！

内沃和罗森（$2012年$）的论文就是这样一篇。在这里，我们必须对工具变量做出的关键额外假设是：工具变量与第二阶段误差项的相关性，和内生 / 处理变量与误差项的相关性符号相同，但程度更弱。这会基于工具变量相关性具体处于 $0$ 到内生 / 处理变量相关性之间的哪个位置，生成一个部分识别的估计值范围。

Nevo, Aviv, and Adam M. Rosen. 2012. “Identification with Imperfect Instruments.” Review of Economics and Statistics 94 (3): 659–71.

另一篇是弗洛雷斯和弗洛雷斯 - 拉古内斯（$2013年$）的论文。他们在这里依赖的假设是单调性 —— 具体来说，工具变量的效应在组内符号相同。基于此，他们再次能够利用工具变量总是朝着同一方向产生影响这一点，得出估计值的部分识别边界。

Flores, Carlos A., and Alfonso Flores-Lagunes. 2013. “Partial Identification of Local Average Treatment Effects with an Invalid Instrument.” Journal of Business & Economic Statistics 31 (4): 534–45.

当然，这些方法要求我们从理论上检验我们的工具变量是否满足这些替代假设。希望它们是满足的！但或许我们还能再添点火力。

另一大且不断发展的关于无效工具变量的文献，确实很强调”工具变量”中的”$s$“（即强调工具变量的复数形式，也就是多个工具变量）。换句话说，它聚焦于有多个工具变量可用的情况。而且，确实存在有多个工具变量可用的情况。这种情况常出现在工具变量有许多不同子成分，且每个子成分都能被单独测量并作为单独工具变量的地方。例如，或许工具变量是分配到某个特定组，这会产生一大堆二元变量的固定效应工具变量，就像本章前面典型设计里的一些设计：随机分配的法官设计，或者孟德尔随机化（用基因作为工具变量）。基因的每个部分（有很多个）都可以是它自己的工具变量！

如果你回想一下第 $13$ 章关于正则化回归，特别是 $LASSO$ 的讨论，它在这里也被应用可能就不足为奇了。$LASSO$ 是一种修改普通最小二乘法的方法，它要求 $OLS$ 在最小化残差平方和的同时，也最小化系数绝对值的和，实际上是鼓励它把变量完全从模型中剔除。

$LASSO$ 在这里可以像平常一样使用，在控制变量和工具变量中挑选最重要的预测因子（切尔诺朱可夫、汉森和斯宾德勒，$2015$）。不过，$LASSO$ 也能用来帮助识别无效工具变量。通过在第二阶段纳入所有工具变量，并用 $LASSO$ 估计模型，你会预期，在控制了处理变量后，有效工具变量将没有预测能力，因此会被从第二阶段剔除（回到它们所属的第一阶段）⁶¹。温德梅杰等人（$2019$）改进了这种方法，使得只要无效工具变量的数量少于一半，就有望奏效。

Chernozhukov, Victor, Christian Hansen, and Martin Spindler. 2015. “Post-Selection and Post-Regularization Inference in Linear Models with Many Controls and Instruments.” American Economic Review 105 (5): 486–90. Windmeijer, Frank, Helmut Farbmacher, Neil Davies, and George Davey Smith. 2019. “On the Use of the LASSO for Instrumental Variables Estimation with Some Invalid Instruments.” Journal of the American Statistical Association 114 (527): 1339–50.

我们有大量工具变量（其中很多我们认为是有效的）这一事实，或许也能让我们的标准方法本身就奏效。科莱萨尔等人（$2015$）表明，在存在大量工具变量的情况下，在某些条件下，即便存在一些无效工具变量，两阶段最小二乘法也能正常发挥作用⁶²。由于 $2SLS$ 的第一阶段会在一定程度上综合所有工具变量的效应，有效工具变量能抵消无效工具变量的影响。这种情况并非总能发生，但他们阐述了其发生的条件，这些条件依赖于有效工具变量和无效工具变量效应之间的关系。他们的结论也并非来自常规的 $2SLS$；他们使用的是带有偏差校正的版本。

Kolesár, Michal, Raj Chetty, John Friedman, Edward Glaeser, and Guido W. Imbens. 2015. “Identification and Inference with Many Invalid Instruments.” Journal of Business & Economic Statistics 33 (4): 474–84.

你会从第 $9$ 章 “寻找前门路径” 以及分离前门路径的内容中，认出这里的很多逻辑。↩︎
在所有这些内容中，我们可以回顾第 $13$ 章的术语：某个变量是内生的，这是由于它是由另一个因素决定的。但另一个变量是外生的，或者我们可以说它是”一个外生的变异来源”，因为它不是由系统中的其他任何事物引起的。↩︎
或者至少（存在的）后门路径是可以通过控制（变量等方式）轻易关闭的。↩︎
实际上，在典型情况下，步骤 3 和步骤 4 对你的估计其实并没有什么影响 —— 如果你用工具变量解释的处理（变量）去解释常规的结果（变量），最终还是会得到结果中被解释的部分。但我还是把它们保留下来，因为我觉得这样能让概念更清晰。↩︎
更准确地说，我们必须满足至少两个假设。不过这两个是真正具有明星效应、人人都在谈论的假设。常规的普通最小二乘法假设同样适用，还有单调性假设 —— 即工具变量与内生变量之间的关系总是同号（或为零）。我们会在 “它是如何实施的？” 部分更多地讨论单调性。↩︎
并且我们使用的是工具变量的标准线性单处理、单工具变量版本。↩︎
回顾上面的 $图19.2$，这就是在分离出两个变量中由 $Z$ 解释的部分后，你得到的斜率。毕竟，斜率是 “上升量”（$Z$ 增加时 $Y$ 增加了多少）除以 “运行量”（给定 $Z$ 的增加量时 $X$ 增加了多少），这是思考斜率的另一种方式。↩︎
这和用连词开头写句子不同，用连词开头写句子我们完全可以做到↩︎
它需要有多强呢？好问题。我们会在下一部分讨论这个↩︎
有效性有时也被称为 “排除限制”，因为它是这样一个假设：在纳入 $Z \to X$ 这条路径后，Z 可以合理地从 Y 的模型中排除。↩︎
我们可能不想对 $C$ 进行控制（尽管这么做不会让 $Z$ 失效）。$C$ 的存在会增强 $Z$ 对 $X$ 的预测效应，这取决于那些箭头是正还是负。所以对 $C$ 进行控制可能会让 $Z$ 的相关性变弱。这挺让人意外的 —— 在这种情况下，我们对 $X \to Y$ 的识别，实际上因为没能识别 $Z \to X$ 而得到了帮助。↩︎
在这个图表里， $C$ 本身也会是一个很棒的工具变量。↩︎
事实上，只要我们讨论的是线性模型，由 $Z$ 驱动的 $X$ 的变异就会与由 $Z$ 驱动的 $A$ 的变异完全相关。↩︎
因为在这种情况下，我们需要对一些因素进行控制才能让有效性成立，所以我们会说，有效性是在 A 和 B 的条件下成立的。↩︎
这种对有效性的关注是好的，但有时会以忽略其他工具变量方面的考虑为代价，比如单调性（我们会讲到），甚至是像数据质量如何这样的基础问题。↩︎
就我个人而言，当数据是全国层面的时候，我基本上从不相信有效性假设 —— 对于如此庞大且相互关联的事物来说，关闭所有必要的后门路径似乎是不可能的！“又一篇宏观工具变量论文？”我笑着，翻了个白眼，对着世界啧啧不已。这些傻瓜！然后他们获得成千上万的引用，可能还得了诺贝尔奖之类的，而我在凌晨 $3$ 点独自坐在黑暗中，琢磨着如果我向一家掠夺性的付费出版期刊投稿，是否真的有人会注意到。↩︎
如今，在某一群因果推断爱好者中，说你要用降雨量作为工具变量，算是一种近乎好笑的内部笑话。这个笑话对所有普通人来说完全无法理解，就像因果推断爱好者们说的那些不是笑话的话一样。↩︎
一般来说，同一个工具变量被用作多个不同处理变量的工具变量，应该引起怀疑。毕竟，这意味着你的工具变量会导致多种同样会导致你所关注结果的事物。所以，为了消除所有不通过你所研究的处理变量的路径，你将不得不对工具变量影响的所有其他处理变量进行控制。这往往是不可能的。↩︎
其他学者还有其他怀疑降雨量有效性的理由；这只是其中之一。如果你感兴趣的话，萨森斯在她的论文中列举了其他一些研究。↩︎
即使是像孟德尔随机化这样真正随机的情况，有效性也并非有保障。如果某个特定的基因影响了你所关注的处理，但也影响了另一个影响结果的变量（比如同一个基因影响身高和体重，而身高和体重都影响运动能力），那这就是一条糟糕的路径！↩︎
我们在这篇论文中运行的是一个选择模型，而非使用工具变量，但这种方法同样需要一个有效性假设。↩︎
当然，也许我们遗漏了一些东西。也许评估者如果因为某种原因不喜欢某个人，就更有可能出错。思考是什么驱动了你的工具变量非常重要，就像在常规的非工具变量情境下，你会思考是什么驱动了你的处理变量一样。↩︎
如果你还记得的话，我们早在第 $10$ 章关于处理效应的内容里就已经讨论过这个了。↩︎
不过，这个词是暗示对设计的兴奋，还是对工具变量这一整体概念的指责，取决于说话的人是谁。↩︎
这些通常是会被重复使用的工具变量和处理变量的组合。如果你有一个工具变量被反复用于不同的处理变量，那就意味着这个工具变量不是很好，因为对于任何一个处理变量来说，你可能会有一些从工具变量到结果的开放前门路径，这些路径不经过该处理变量 —— 它们也经过其他处理变量。↩︎
你可能会注意到，这些留存下来的设计大多都有某种明确的随机化因素；你可以自行体会这一点。Aizer, Anna, and Joseph J. Doyle Jr. 2015. “Juvenile Incarceration, Human Capital, and Future Crime: Evidence from Randomly Assigned Judges.” The Quarterly Journal of Economics 130 (2): 759–803.↩︎
哎呀，一个用于两种处理的工具变量！不过，与降雨量相比，这个情况的问题要小一点，因为这些活动在某种程度上是相互替代的。用于任何一种处理的工具变量，也会是用于”不接受该处理”的工具变量。你可以反过来把征兵顺序看作是 “高中毕业后做什么” 的工具变量。↩︎
好家伙，那个叫安格里斯特的家伙还真是到处搞研究啊，是吧？↩︎
或者说 TSLS，但要是你问我的话，那个看起来就是不对劲。↩︎
如果目标是关闭与 $X$ 中未被 $Z$ 解释的部分相关的后门路径，那我们为什么不采用残差而非预测值，然后在第二阶段将残差与原来的 $X$ 一起进行控制呢？嗯…… 你是可以这么做的！这被称为控制函数法。在标准的线性工具变量中，这会产生和两阶段最小二乘法基本相同的结果。但它在非线性工具变量方面有一些重要的应用，我会在本章后面讲到。↩︎
出于这个原因，一般来说，你会想要使用专门为工具变量设计的软件命令来运行工具变量分析。↩︎
异方差稳健标准误在修正这个问题上只能说做得还行。↩︎
很多优秀的教材会更深入地探讨（广义矩方法等内容）。令人颇为惊讶的是，其中很多教材的标题就是”广义矩方法”，再加上几个额外的词。↩︎
怎么推导呢？让我们通过假设所有变量都是均值为零的来简化，这样就不需要截距项了。从 $\sum Xr = 0$ 开始。然后代入 $r = Y - \beta_1 X$，得到 $\sum X(Y - \beta_1 X) = \sum (XY - \beta_1 X^2) = 0$。解 $\beta_1$ 可得 $\beta_1 = \sum (XY) / \sum (X^2)$，这和我们之前得到的结果一致。↩︎
就均值和普通最小二乘法而言，它得出的解和我们已经得到的是一样的 —— 取样本数据的均值。但我们得到这个结果的方式是不同的。↩︎
在过度识别的情况下，广义矩方法（GMM）实际上无法精确满足所有的矩条件，因此它需要为这些条件设定权重，以决定哪些条件更为重要。大致来说，它会根据满足这些条件的难易程度来赋予相应的权重。↩︎
通过以这种方式构建研究，即只观察第一轮对第二轮的影响，他们避免了 “反射” 问题，而这是研究朋友如何影响你的行为时的一个常见问题；这些研究也被称为 “同伴效应” 研究。“反射问题” 由常年与诺贝尔奖无缘的查尔斯・曼斯基提出，它指出了如果你想了解朋友如何影响你时所面临的识别问题…… 问题在于你也会影响你的朋友！存在一个反馈循环，而我们知道因果图最讨厌反馈循环。通过设置第一轮和第二轮，我们知道因果箭头指向的方向。↩︎
我们将使用他们对第二轮参与者子样本所做的分析，这些参与者被告知了其第一轮朋友所做的购买决策。他们对未被告知的第二轮参与者还做了另一项分析，发现没有效果。↩︎
这里的术语源自处理变量为二元的情况。“总是接受者”无论工具变量如何，总是会接受处理；“从不接受者”无论工具变量如何，永远不会接受处理。这些术语应用于像我们这里（锻炼时长这种）连续的处理时有点奇怪，但核心思想是一样的。↩︎
如果样本中存在受工具变量影响方向相反的人，把其中一类称为”依从者”、另一类称为”违背者”有点随意，但真正重要的是样本中存在这两种不同影响方向的人。↩︎
会有很多人告诉你，局部平均处理效应和依从者中的平均处理效应是一回事。但这其实并不对，它依赖于工具变量对所有依从者影响程度相同这一假设。↩︎
重要的是，这里的问题不是梅西受到了负面的影响，而是她受到影响的方向与其他人不同。一致性在这里很有帮助。如果每个人都有负面效应，那就不会有问题。那样的话，权重上的所有负值都会被抵消掉。这就好比将原来的局部平均处理效应计算式乘以$-1 / -1 = 1$，不会有任何差别。↩︎
不过，这可能不是解决弱工具变量问题的最佳方法…… 读到本节末尾就知道了。↩︎
这有时被称为识别不足检验。↩︎
如果你使用的是其他工具变量命令，这可能不适用。很多其他命令，比如 $AER$ 包中的 $ivreg$，需要用 $summary(m,diagnostics=TRUE)$ 来实现。↩︎
如果你做了两阶段最小二乘法，这还会给你不同的相关临界值，我们稍后会讨论。↩︎
还有一个经验法则，笼统地认为你的$F$统计量总体上必须达到 $10$ 或以上。这无疑比在表格中查找数值更容易记住，但也非常粗略。就像统计学中许多一刀切的数值一样，这一传统或许最好被摒弃。↩︎
经典的统计学啊，总是有办法让人失望。↩︎
若想了解它们具体是如何组合的更详细内容，请参见莫格斯塔德、托尔戈维茨基和沃尔特斯（$2021年$）的研究，该研究还表明，要让这种组合奏效，在其他工具变量保持不变的情况下，单调性对于你的每个工具变量都必须成立。↩︎
老实说，这是目前为止最常见的情况。我们当初为什么又要费神去做那些过度识别检验呢？↩︎
这可别和检验弱工具变量的 $Anderson-Rubin \space test$ 弄混了。我知道，这很容易混淆。↩︎
有限信息极大似然法通常与支持两阶段最小二乘法（2SLS）和广义矩方法（GMM）的同一款工具变量（IV）估计软件兼容，而且本章前面已经讨论过使用有限信息极大似然法进行估计的代码。↩︎
而事实上，第三种方法会在下一节中出现。↩︎
经济学家往往会这么做。总的来说，经济学家很可能会审视一种让非线性模型难以处理的研究设计，比如固定效应或工具变量，然后认为，解决非线性问题所带来的额外困难，会引发比解决的问题更多的新问题，于是就采用线性模型。我承认自己对这种观点抱有同感。↩︎
而且工具变量本来就比常规回归的精度低。我们真的需要更低的精度吗？↩︎
或者至少极大似然法是这样的。这些估计量的两阶段方法与下一节描述的控制函数法更为相似，但会对第一阶段的预测值进行一些调整，以对其进行适当的缩放。↩︎
无论你是分离出被解释的变异，还是控制掉未被解释的变异，得到相同的结果，这在直觉上是很说得通的。而它实际上确实有效，这也算是一种慰藉吧。↩︎
对于此处列出的 $R$ 语言和 $Stata$ 版本，它们默认采用的是结合极大似然法的非控制函数法。这种做法并无问题。↩︎
从技术角度来说，该系统的非线性特性会使得即便没有实际工具变量，你也能对模型进行估计 —— 这一点与线性工具变量不同；但通常认为，省去工具变量是个相当糟糕的做法。↩︎
在进行此项操作时，我们必须格外小心：要纳入必要的适当控制变量，以解释我们这么做必然会引入的对撞偏差 —— 正如之前在介绍类似方法时所提及的那样。↩︎
不过，极大似然法并非如此。↩︎

第19章 工具变量法