作者： 邹恬华 (对外经济贸易大学)
邮箱： 201802060@uibe.edu.cn
Source: Cinelli, Carlos, Jeremy Ferwerda, and Chad Hazlett. 2020. “Sensemakr: Sensitivity Analysis Tools for OLS in R and Stata.”, -Link-
目录
1. 我们为什么要做敏感性分析？
1.1 一个小例子
1.2 敏感性分析解决的问题
2. 开始使用 sensemakr 吧！
2.1 命令安装与数据下载
2.2 基础语法结构
2.3 其它进阶选项
3. 反思：敏感性分析应该用来做什么？
4. 参考文献
5. 相关推文
温馨提示： 文中链接在微信中无法生效。请点击底部「阅读原文」。或直接长按/扫描如下二维码，直达原文：
存在可能的遗漏变量问题且数据不可得时，可以使用敏感性分析进行稳健性检验。 即检验：遗漏变量需要达到多强才能使得之前的研究结果改变？ 敏感性分析可以通过 Stata 命令 sensemakr 进行。 开始敏感性分析之旅吧！
在 2003 年和 2004 年，达尔富尔政府策划了一场针对平民的恐怖暴力运动，估计造成 20 万人丧生。我们期望评估，在此事件中直接受伤或致残，是会使人更具有复仇心理，从而不愿与施暴者和平相处，还是更有动力看到这种暴力结束，支持呼吁和平？Darfur 数据集(Hazlett, 2019) 是针对乍得东部达尔富尔难民的一项调查数据，我们建立 OLS 回归模型如下：
PeaceFactor：和平意愿，被解释变量
DirectlyHarmed：是否被直接伤害的虚拟变量，我们关心的核心解释变量
Female：是否为女性的虚拟变量，控制变量之一
OtherControls：除 Female 外的控制变量，包括 Age 年龄， HerderDar 是否为达尔富尔牧民， FarmererDar 是否为达尔富尔农民， PastVoted 是否之前在选举中投票
Village：被调查者居住村庄的固定效应
：我们关心的系数，即被直接伤害如何影响村民的和平意愿
回归结果表明，系数 显著为正，即被直接伤害平均而言更会激发村民和平意愿。但有学者认为以上模型可能存在遗漏变量问题：
问题在于上述可能的遗漏变量数据不可得。在此情况下，想要说明上述结果的稳健性，也可以采取以下思路：
当遗漏变量达到何种强度时，之前的结果才不可信？若这种要求显然过于苛刻，则说明之前的结果大致可信。
此种思路即为 「敏感性分析」。
注：“强度”指与核心解释变量相关性，下同。
结合上述的小例子，我们总结出敏感性分析想要解决的问题：
sensemakr 是 Cinelli et al.(2020) 开发的针对 OLS 估计的敏感性分析的 Stata 命令，同时提供功能类似的 R 语言版本。可以通过向 Stata 命令行输入以下命令下载：
. ssc install sensemakr, replace all前文例子中提到的 Darfur 数据集可以在作者的 Github 主页下载：https://github.com/resonance1/sensemakr-stata/ (Hazlett, 2019)。亦可直接执行如下命令，将数据下载到当前工作路径下：
. net get sensemakr.pkg
checking sensemakr consistency and verifying not already installed...
copying into current directory...
copying darfur.dta
ancillary files successfully copied.打开数据文件：
use darfur.dta, clear数据集中包含 1.1 节小例子中的全部数据。
sensemakr depvar covar [if] [in], treat(varlist) benchmark(varlist)depvar：被解释变量covar：解释变量treat：此处申明最为关心的核心解释变量，这也意味着填写 covar 时的顺序不重要。benchmark：与可能存在的遗漏变量进行对比的变量，即最后得出的结果为遗漏变量需要达到对比变量强度的多少倍才能影响之前的研究结果。在上述的 Darfur 例子中，该语句可以写为：
sensemakr peacefactor directlyharmed age farmer herder ///
pastvoted hhsize female i.village_, ///
treat(directlyharmed) benchmark(female)Stata 汇报结果如下：
. net get sensemakr.pkg // 下载数据
. use darfur.dta, clear // 导入数据
. sensemakr peacefactor directlyharmed age farmer herder ///
pastvoted hhsize female i.village_, ///
treat(directlyharmed) benchmark(female)
DOF = 783
q = 1.00
alpha = 0.05
reduce = TRUE
H0 = 0
Treatment | Coef. S.E. t(H0) R2yd.x RV_q RV_qa
---------------- -----------------------------------------------------------
directlyharmed | 0.0973 0.0233 4.1844 0.0219 0.1388 0.0763
Bound | R2dz.x R2yz.dx Coef. S.E. t(H0) Lower CI Upper CI
--------------- --------------------------------------------------------------
1.00x female | 0.0092 0.1246 0.0752 0.0219 3.4389 0.0323 0.1182
2.00x female | 0.0183 0.2493 0.0529 0.0204 2.6002 0.0130 0.0929
3.00x female | 0.0275 0.3741 0.0304 0.0187 1.6281 -0.0063 0.0670
Extreme Bound | R2dz.x R2yz.dx Coef.
--------------- ----------------------------
1.00x female | 0.0092 1.0000 0.0347
2.00x female | 0.0183 1.0000 0.0084
3.00x female | 0.0275 1.0000 -0.0121部分结果和统计量的含义及解读如下：
R2yd.x：Partial R2 of the treatment with the outcome (R2yd.x):
An extreme confounder (orthogonal to the covariates) that explains 100 percent
of the residual variance of the outcome, would need to explain at least 2.19 percent of the residual variance of the treatment to fully account for the observed estimated effect.
RV_q：Robustness Value, q = 1.00 (RV_q):
Unobserved confounders (orthogonal to the covariates) that explain more than
13.88 percent of the residual variance of both the treatment and the outcome are
strong enough to bring the point estimate to 0 (a bias of 100 percent of the original estimate). Conversely, unobserved confounders that do not explain more than 13.88 percent of the residual variance of both the treatment and the outcome are not strong enough to bring the point estimate to 0.
RV_qa：Robustness Value, q = 1.00, alpha = 0.05 (RV_qa):
Unobserved confounders (orthogonal to the covariates) that explain more than 7.63 percent of the residual variance of both the treatment and the outcome are strong enough to bring the estimate to a range where it is no longer 'statistically different' from 0 (a bias of 100 percent of the original estimate), at the significance level of alpha = 0.05. Conversely, unobserved confounders that do not explain more than 7.63 percent of the residual variance of both the treatment and the outcome are not strong enough to bring the estima te to a range where it is no longer 'statistically different' from 0, at the significance level of alpha = 0.05
Bounds on Omitted Variable Bias:
The table shows the maximum strength of unobserved confounders, bounded by a multiple of the observed explanatory power of the chosen benchmark covariate(s) with the treatment and the outcome.
我们逐个解释汇报指标的经济含义：
R2yd.x：，**在控制其它控制变量情况下，核心解释变量 D 对被解释变量 Y 的偏 **。此处等于 0.0219，意味着在最坏条件下（即遗漏变量解释被解释变量所有剩余方差），遗漏变量需要至少解释核心解释变量 2.19%的剩余方差，才能彻底去除之前观察到的核心解释变量引起的效应。
RV_q：使得估计系数恰好为零的稳健性值。此处等于 0.1388，意味着遗漏变量需要最少同时解释 13.88%的核心解释变量与被解释变量的剩余方差，才能彻底去除之前观察到的核心解释变量引起的效应。
RV_qa：使得估计系数不再在 95%水平上显著的稳健性值。此处等于 0.0763，意味着遗漏变量需要最少同时解释 7.63%的核心解释变量与被解释变量的剩余方差，才能使得 95%置信区间下限为零。
然而，上述指标不够直观，并且我们很难论证为什么遗漏变量确实能够达到上述要求，并以此说明我们模型的稳健性。为此我们引入对比变量，尝试论证当遗漏变量强度只要小于对比变量时，之前的估计结果是有效的。
在 Darfur 这个例子当中，我们**将对比变量设置为 female**，意味着我们尝试论证当遗漏变量强度只要小于 female 时，之前的估计结果是有效的。这在一定程度上是说得通的，因为当地针对特定性别的伤害行为要多于针对城市中心人群/富有人群/特定政治立场人群。
Stata 汇报结果的倒数第二个表格汇报了有关对比变量的信息：
R2yz.dx：，在控制所有控制变量情况下，遗漏变量 Z 对被解释变量 Y 的偏 *。此处等于 0.1246，小于RV_q的值。
R2dz.x：，在控制其它控制变量情况下，遗漏变量 Z 对核心解释变量 D 的偏 *。此处等于 0.0092，小于RV_q的值。
R2yz.dx、R2dz.x均小于RV_q意味着，当遗漏变量强度小于 female 时，之前的估计结果是有效的。R2dz.x 的值甚至低于 R2yd.x ，意味着最坏情况下（即遗漏变量解释被解释变量所有剩余方差），当遗漏变量强度小于 female 时，之前的估计结果是有效的。
为了使得上述对比变量结果更具有说服力，我们继续探讨，若遗漏变量的强度时对比变量 female 的 n 倍时，之前的研究结果是否仍然成立。
我们使用 sensemakr 命令中的 contourplot 与 tcontourplot 绘制系数、t 统计量等值线图：
sensemakr peacefactor directlyharmed age farmer herder pastvoted hhsize ///
female i.village_, treat(directlyharmed) benchmark(female) contourplotsensemakr peacefactor directlyharmed age farmer herder pastvoted hhsize ///
female i.village_, treat(directlyharmed) benchmark(female) tcontourplot绘图结果如下，上图为系数 等值线图，下图为 t 统计量等值线图。
图中横纵轴坐标分别为 与 。
上图等值线为回归系数 值，红线为 。四个数值点分别指不加入遗漏变量、加入 female 相同强度的遗漏变量、加入female2 倍强度的遗漏变量、加入female3 倍强度的遗漏变量的情况。四个数值点的小括号中汇报了其对应计量模型的回归系数 值，可以观察到均大于零，也即全部四个数值点位于红线左侧。这意味着即使加入female3 倍强度的遗漏变量，也不会使得原估计系数由正转负。
下图等值线为 t 统计量，红线为 1.96（95%置信区间临界值）。四个数值点与上图相同，分别指加入female0~3 倍强度的遗漏变量的情况。四个数值点的小括号中汇报了其对应计量模型的 t 统计量。加入female0~2 倍强度的三个数值点 t 统计量小于 1.96，在红线左侧，加入female3 倍强度的数值点 t 统计量大于 1.96，在红线右侧。这意味着即使加入female2 倍强度的遗漏变量，也不会使得原估计系数由显著变为不显著；但当加入female3 倍强度的遗漏变量，原估计系数变得不显著。
此外，sensemakr 命令还允许绘制最坏情况图，代码如下：
sensemakr peacefactor directlyharmed age farmer herder pastvoted hhsize ///
female i.village_, treat(directlyharmed) benchmark(female) extremeplot绘图结果如下：
三条数值线分别指当遗漏变量解释剩余方差的 100%、75%和 50%的情况，即 分别等于 1、0.75 与 0.5 的情况。横轴的小黑点分别指遗漏变量是female1~3 倍强度的情况。
图中可以看到，即使是最坏情况（遗漏变量解释剩余方差的 100%，图中实线），遗漏变量都需要female2 倍以上强度才能推翻之前结论。而对于 75%与 50%的情况，遗漏变量则需要female3 倍以上强度。
sensemakr 还有其它选项可供选择，主要是修改上述分析过程中的一些参数。以下列出这些选项：
alpha(real)：设置显著性水平，默认值 0.05gbenchmark(varlist)：添加单个对比变量使用前文的 benchmark(varlist) ，而当对比变量有一组时，使用gbenchmark(varlist)gname(string)：给一组对比变量取名，默认值为用组内变量名用连字符连接kd(numlist)与ky(numlist)：遗漏变量分别相比核心解释变量与被解释变量的强度，默认值为（1 2 3）。用户单独设定 kd 与 ky 其中一个时，另一个会自动被设置成相同值；若用户两个都设定，则必须设置成相同值。latex(filename)：将结果输出到 filename.texnoreduce：遗漏变量对回归系数 的影响方向，默认值为向 方向影响q(real)：指定估计的哪一部分需要被解释为有问题的。默认值为 1，这意味着将当前效果 100％有问题，核心解释变量的真实效果为 0。r2yz(numlist)：默认值为 1，意味着遗漏变量解释剩余方差的 100%。如果输入(.5 .75)意味着输出遗漏变量分别解释剩余方差的 50%和 75%的结果。surppress：不要汇报详细概念描述与解释等值线图相关：clines：等值线的条数，默认值为 7clim(numlist)：横纵坐标的最小值和最大值，最大范围为(0 1)最坏情况图相关：r2yz：设置线条 取值，默认值为(1 0.75 0.5)，最多可以设置四条线elim(numlist)：设置横坐标的最小值和最大值，最大范围为(0 1)。纵坐标范围自动选择，不可被设定。female 强两倍以上，这个事件不太可能发生。Note：产生如下推文列表的 Stata 命令为：
lianxh 稳健性检验 安慰剂检验
安装最新版lianxh命令：
ssc install lianxh, replace
New！ Stata 搜索神器：
lianxh和songblGIF 动图介绍
搜： 推文、数据分享、期刊论文、重现代码 ……
👉 安装：
. ssc install lianxh
. ssc install songbl
👉 使用：
. lianxh DID 倍分法
. songbl all
🍏 关于我们
联系客服
