👇 连享会 · 推文导航 | www.lianxh.cn
连享会视频课 · 因果推断实用计量方法
作者： 曹琳君 (南开大学) ；陈波 (深圳大学)
邮箱： linjuncao@yeah.net ；1900123011@email.szu.edu.cn
Source： Practical Guide to Counterfactual Estimators for Causal Inference with Time-Series Cross-Sectional Data." Working Paper, Stanford University, -Link-
温馨提示： 文中链接在微信中无法生效。请点击底部「阅读原文」。或直接长按/扫描如下二维码，直达原文：
目录
1. 续前集
2. 面板数据因果推断的估计量选择
3. 面板数据因果推断的基本原理
3.1 潜在的时变混杂因素的检验原理
3.2 安慰剂测试检验原理
4. Stata 演示
4.1 Model Selection
4.2 Uncertainty Estimates
4.3 Wald Test
4.4 Equivalence Test
4.5 PlaceboTest
5. 小结
6. 致谢
7.参考资料
8. 相关推文
在上一篇推文中 (fect：基于面板数据的因果推断（上）)，我们为大家详细介绍了面板数据进行因果推断的三种有效估计：FE 估计，IFE 估计和 MC 估计。那么我们该如何在这三种估计中进行选择，并对估计结果进行检验呢？
本文将在对 Licheng Liu(Tsinghua) 、Ye Wang (NYU) 和 Yiqing Xu (Stanford) 三位学者的成果进行编译的基础上，为大家介绍面板数据进行因果推断时估计量选择的注意事项和基本原理，并演示具体的 Stata 操作。
基于 fect 命令，我们可以生成三个关于面板数据因果推断的统计量：FE 估计、IFE 估计和 MC 估计。通过对比这三个统计量的估计效果，Liu et al. (2019) 发现：
与 FE 估计相比，IFE 估计和 MC 估计能对反事实进行更好的预测。
当只有少量时变因素存在并且每个时变因素都表现出相对强的信号时，IFE 估计的表现优于 MC 估计，而当存在大量弱时变因素时，MC 的表现优于 IFEct 。
然而，在实际的操作中，我们并没有办法主观的判断影响模型估计效果的因素的多少。所以在 Stata 操作中，我们可以根据不同模型对同一数据的相对预测力进行模型筛选。模型的相对预测能力可以通过 MSPE 测量判断，我们将在下文进行具体的 Stata 演示。
为了检验是否存在潜在的时变混杂因素，Liu et al. (2019) 引入了两种统计检验方法: 联合假设检验 (F 检验) 和双单侧检验 (TOST 检验) 。本环节主要介绍这两种检验的前提假设及统计量构建。
(1) 联合假设检验 (F 检验)
在上列式子中， 和 (1-m) 是预处理周期总数， 。
然而，Hartman and Hidalgo (2018) 认为关于干预前阶段残差和是否为 0 的检验存在两个潜在的问题：
当观察的样本数量有限时，F 检验的结果不能排除残差平均值为 0 并不代表残差平均值为 0 ，即 F 检验不能提供证据不代表没有证据;
当观察的样本数量很大时，一个小的混杂因素 (或几个离群值) 在因果估计中只会造成可忽略的偏差量，进而导致检验结果在大多数情况下拒绝零值假设。这对于 TSCS 数据分析来说尤其成问题，因为这种情况意味着存在一些不可分解但无影响的时变混杂因素的可能性很高。
为了解决上述的问题，Liu et al. (2019) 提出了上述检验的等价检验：双单侧检验 (TOST 检验)。
(2) 等价检验：双单侧检验 (TOST 检验)
其中，，为预先指定的参数。如果拒绝零值假设，则意味着相反的成立概率很高，即对于任何 。换句话说，如果我们收集了足够的数据，并表明干预前的剩余平均值落在预先指定的狭窄范围内，我们就获得了一个证据来支持无时变混杂假设的有效性，所以 可以被称为等价边界。
我们假设干预组中每个单位的干预开始时间比实际开始时间早 S 期，并使用相同的反事实估计量来获得对 ATTs 关于 的估计。同时，我们可以估计 S 干预前的总体 ATT。如果关于混杂时变因素的假设不成立，这个 ATT 估计在统计学上将与零没有区别并利于假设的有效性。
此外，需要注意的是：在实际操作中，S 值的设置要适中。如果 S 设置得太大，估计中剩余的预处理时间就越少。然而，如果 S 和 都太小，测试可能会效力不足。
为确定模型的具体参数，作者引入了交叉验证方法。具体思路为：将样本分为 k 组 (可以根据 kfold 命令设定具体组数，默认组数为 10 )，将每个子集做一次验证集，其余 k-1 组子集作为训练集，得到 k 个模型。这 k 个模型分别在验证集中评估结果，随后所得的误差 MSPE (Mean Squared Prediction Error) 作为交叉验证误差。我们选取 MSPE 最小的模型参数。
交叉验证主要应用与以下三种情形。
情形一
当我们不知道该选择上述三种方法的哪一种时，可以根据设定 method("both") 同时对三种方法的预测情况进行评估。具体命令如下：
//scenario 1
fect Y, treat(D) unit(id) time(time) cov(X1 X2) method("both") r(4) nlambda(15)
其中 r(4) 表示在 ife 方法中，最多引入四个共同因子。nlambda(15) 表示在 mc 方法中，为参数 lambda() 设定 15 个格数 (默认格数为 10 ) 。
下面展示了交叉验证结果，最后一行会提示我们该选择哪种方法。可以发现，根据最小 MSPE 的原则，我们应该选择 r = 2 时的 ife 方法，此时 MSPE 仅为 4.772 。需要注意的是，在 r = 0 时，ife 方法等价于 fe 。
Balanced Panel Data
---------------------------------------------------------------------------------
Cross Validation...
fe r=0 force=two-way mspe=8.538
ife r=1 force=two-way mspe=5.437
ife r=2 force=two-way mspe=4.772
ife r=3 force=two-way mspe=5.339
ife r=4 force=two-way mspe=6.146
mc: lambda=.0178 lambda.norm=1 mspe=8.538
mc: lambda=.0105 lambda.norm=.588 mspe=6.98
mc: lambda=.0062 lambda.norm=.346 mspe=5.693
mc: lambda=.0036 lambda.norm=.203 mspe=5.187
mc: lambda=.0021 lambda.norm=.119 mspe=5.272
mc: lambda=.0013 lambda.norm=.07 mspe=5.309
mc: lambda=.0007 lambda.norm=.041 mspe=5.329
mc: lambda=.0004 lambda.norm=.024 mspe=5.391
mc: lambda=.0003 lambda.norm=.014 mspe=5.605
mc: lambda=.0001 lambda.norm=.008 mspe=6.386
mc: lambda=.0001 lambda.norm=.005 mspe=8.408
mc: lambda=.0001 lambda.norm=.003 mspe=8.461
mc: lambda=0 lambda.norm=.002 mspe=8.493
mc: lambda=0 lambda.norm=.001 mspe=8.511
mc: lambda=0 lambda.norm=.001 mspe=8.523
choose fe/ife model with optimal r=2
下图展示了 r = 2 时的 ife 方法的反事实估计结果：
情形二
当我们确定选择 mc 方法，需要对惩罚项参数 lambda 进行明确设定时，可以在 lambda() 中输入一系列备选参数，通过比较 MSPE 最终选择最后的参数。具体命令如下：
//scenario 2
. fect Y, treat(D) unit(id) time(time) cov(X1 X2) method("mc") ///
lambda(0.001 0.002 0.003 0.004 0.005)
在上述命令中，我们直接确定使用 mc 方法，并引入五个备选参数。交叉验证结果如下所示。根据最小 MSPE 原则，在使用 mc 方法时，我们应该设定 lambda(0.004) 。
Balanced Panel Data
---------------------------------------------------------------------------------
Cross Validation...
mc: lambda=.001 lambda.norm=.056 mspe=5.315
mc: lambda=.002 lambda.norm=.112 mspe=5.283
mc: lambda=.003 lambda.norm=.168 mspe=5.208
mc: lambda=.004 lambda.norm=.224 mspe=5.202
mc: lambda=.005 lambda.norm=.28 mspe=5.374
optimal lambda=.004 in mc model
下图展示了 lambda = 0.004 时的 mc 方法的反事实估计结果：
情形三
当我们明确选择某一方法时，如 ife 或 mc ，直接设定 cv 进行结果。例如当我们选定 ife 方法时，设定最大共同因子数即可。而当选定 mc 方法时，则可以直接运行，因为默认 lambda 的格数为 10 。我们以 ife 方法为例：
//scenario 3
. fect Y, treat(D) unit(id) time(time) cov(X1 X2) method("ife") r(4) cv cvtreat
上述命令中，我们将最大共同因子设定为 4 。同时，使用 cvtreat 将测试集的样本选取限制在处理组。
下面汇报了交叉验证结果。依据 MSPE ，我们应该选取 r = 2 时的 ife 方法。
Balanced Panel Data
---------------------------------------------------------------------------------
Cross Validation...
fe r=0 force=two-way mspe=8.616
ife r=1 force=two-way mspe=5.736
ife r=2 force=two-way mspe=4.905
ife r=3 force=two-way mspe=5.659
ife r=4 force=two-way mspe=6.661
optimal r=2 in fe/ife model
除了估计反事实结果的精确值之外，我们还可以估计其大概出于哪一个范围 (标准差) 。fect 命令引入了两种方法计算标准差： bootstrap (默认方法) 和 jackknife ，可以根据 vartype() 具体设定。
我们先以默认的 bootstrap 方法为例：
. fect Y, treat(D) unit(id) time(time) cov(X1 X2) ///
method("ife") r(2) se nboots(100)
上述命令中，nboots(100) 表示抽样 100 次，且不显示抽样点 (如果显示抽样点的话，屏幕中会渐次出现 100 个点，展示抽样进程)。
下图展示了反事实估计结果，灰色阴影部分是 95% 置信区间 ( 1.96 个标准差) ：
我们还可以对上图进行微调。例如通过 preperiod(-14) 和 offperiod(5) 将绘图的时间窗口限定为处理前的 14 期和处理后的 5 期。还可以根据 xlabel() 和 ylabel() 设定横坐标和纵坐标的标题。最后，我们将标准差的估计方法换为 jackknife 。具体命令如下：
. fect Y, treat(D) unit(id) time(time) cov(X1 X2) ///
method("ife") r(2) se preperiod(-14) offperiod(5) ///
xlabel("s") ylabel("ATTs") vartype("jackknife")
下图即为调整后的图形：
为检验模型的精度与稳健性，fect 命令引入了拟合优度检验 (Wald Test) 来测度处理前的趋势。我们只需加入 wald 即可。Wald Test 的原假设是：在不同时期内，样本处理前的残差均值都等于零。具体命令如下：
. fect Y, treat(D) unit(id) time(time) cov(X1 X2) se ///
method("ife") r(2) preperiod(-14) offperiod(5) ///
wald nboots(100)
preperiod(-14) 和 offperiod(5) 等 options的功能与前文相似。 wald 表示引入 Wald Test 检验样本处理签的时间趋势。
具体检验结果如下所示。Wald Test 的 p 值为 0.91 ，与 0.05 相差甚远，无法拒绝原假设，说明样本在处理前的趋势是相似的。
Wald Testing...
Wald Testing: Already Simulated 100 Times
The p-value in wald test is .91
下图展示了具体估计结果，左上角为 Wald Test 的 p 值 。
为评估模式的识别假设是否有效，fect 命令还引入了等效检验 (Equivalence Test) 。该方法主要用于检验处理前 ATTs 的 90% 置信区间是否超过预定范围。
我们先查看 fe 方法是否通过等效检验。具体命令如下：
. fect Y, treat(D) unit(id) time(time) cov(X1 X2) se ///
method("fe") preperiod(-25) offperiod(0) ///
equiTest nboots(100)
equiTest 即为汇报并图示等效检验结果。
具体检验结果如下所示。可以发现，模型并没有通过等效检验 (Fail) 。
Equivalence Test
Equivalence Test...Fail at s=-25
Equivalence Test...Fail at s=-24
Equivalence Test...Fail at s=-22
Equivalence Test...Fail at s=-21
Equivalence Test...Fail at s=-20
Equivalence Test...Fail at s=-19
Equivalence Test...Fail at s=-16
Equivalence Test...Fail at s=-8
Equivalence Test...Fail at s=-5
Equivalence Test...Fail at s=-3
Equivalence Test...Fail at s=-2
Equivalence Test...Fail at s=-1
Equivalence Test...Fail
图示则提供了更为具象化的信息。图中灰色背景为 90% 置信区间，蓝色虚线 (Equiv.Bond) 为等效边界 (残差标准差的 0.39 倍)，绿色虚线 (Min.Bond) 为最小边界。一般以等效边界为准。可以发现，90% 置信区间多处超过了等效边界。这表明 fe 方法没有通过等效检验。
下面我们对 ife 方法进行等效检验，检验结果显示， ife 方法通过了等效检验 (Pass) 。
. fect Y, treat(D) unit(id) time(time) cov(X1 X2) se ///
method("ife") preperiod(-25) offperiod(0) ///
equiTest nboots(100)
Equivalence Test
Equivalence Test...Pass
图示结果也传递了相似的信息， ATT 的 90% 置信区间都处于等效边界之内，说明 ife 方法通过了等效检验。
最后，我们对 mc 方法进行等效检验。检验结果如下。由于我们未设定具体的 lambda ，所以模型会使用默认方法选择合适的 lambda 。最终结果显示， ife 方法通过了等效检验 (Pass) 。
fect Y, treat(D) unit(id) time(time) cov(X1 X2) se ///
method("mc") preperiod(-25) offperiod(0) ///
equiTest nboots(100)
Balanced Panel Data
------------------------------------------------
Cross Validation...
mc: lambda=.0178 lambda.norm=1 mspe=8.538
mc: lambda=.0075 lambda.norm=.422 mspe=6.118
mc: lambda=.0032 lambda.norm=.178 mspe=5.199
mc: lambda=.0013 lambda.norm=.075 mspe=5.307
mc: lambda=.0006 lambda.norm=.032 mspe=5.351
mc: lambda=.0002 lambda.norm=.013 mspe=5.654
mc: lambda=.0001 lambda.norm=.006 mspe=8.389
mc: lambda=0 lambda.norm=.002 mspe=8.475
mc: lambda=0 lambda.norm=.001 mspe=8.511
mc: lambda=0 lambda.norm=0 mspe=8.527
optimal lambda=.003 in mc model
------------------------------------------------
Bootstrapping...
ATT Estimation: Already Bootstrapped 100 Times
------------------------------------------------
Equivalence Test
Equivalence Test...Pass
图示结果也显示，ATT 的 90% 置信区间都位于等效边界内，mc 方法通过了等效检验。
最后， fect 命令还提供了安慰剂检验来查看是否存在过拟合问题。具体做法是选定处理前的某一时间段作为”安慰剂周期“，删除该周期后对模型进行拟合，然后检验该时间段内的 ATT 是否显著不为零。一般情况下，我们将处理前的 3 期 [-2, 0] 作为默认的安慰剂周期，也可以使用 palceboperiod() 自定义安慰剂周期。
我们先对 fe 方法进行安慰剂检验。 palceboTest 即为汇报安慰剂检验结果：
. fect Y, treat(D) unit(id) time(time) cov(X1 X2) se ///
method("fe") placeboTest nboots(100)
检验结果如下。可以发现 fe 方法未通过安慰剂检验 (Fail) 。
Placebo Test...
Bootstrapping...
ATT Estimation: Already Bootstrapped 100 Times
Placebo Test Fail
下图为安慰剂检验的图示结果，左上角为安慰剂检验的 p 值。可以发现，安慰剂检验的 p 值为 0 ，拒绝安慰剂周期内的 ATT 等于 0 的原假设。图中也显示，安慰剂区域 (Placebo Region) 与零值不存在交集。这都表明 fe 方法未通过安慰剂检验。
接下来我们对 ife 方法进行安慰剂检验。
. fect Y, treat(D) unit(id) time(time) cov(X1 X2) se ///
method("ife") r(2) placeboTest nboots(100)
检验结果如下所示， ife 方法通过了安慰剂检验 (Pass) 。
Placebo Test...
Bootstrapping...
ATT Estimation: Already Bootstrapped 100 Times
Placebo Test Pass
图示结果也显示，安慰剂检验的 p 值为 0.354 ，远大于 0. 05，不拒绝原假设。安慰剂区域 (Placebo Region) 与零值存在明显的交叉，这都表明 ife 方法通过了安慰剂检验。
最后，我们对 mc 方法进行安慰剂检验：
. fect Y, treat(D) unit(id) time(time) cov(X1 X2) se ///
method("mc") lambda(0.004) placeboTest nboots(100)
下面的检验结果显示，mc 方法未通过安慰剂检验 (Fail) 。
Placebo Test...
Bootstrapping...
ATT Estimation: Already Bootstrapped 100 Times
Placebo Test Fail
图示结果也传递了相似的信息。安慰剂检验的 p 值为 0 ，拒绝原假设。安慰剂区域 (Placebo Region) 与零值不存在交集。以上信息都表明 mc 方法未通过安慰剂检验。
以上就是关于使用 fect 命令进行面板数据因果推断所涉及的基本原理和具体 Stata 命令。如果想要进一步了解面板数据的因果推断，也可以深入阅读下方的参考资料和相关推文。
在此感谢命令开发者 Licheng Liu (Tsinghua) 、Ye Wang (NYU) 和 Yiqing Xu (Stanford) 三位学者对面板数据因果推断作出的贡献。特别感谢斯坦福大学徐轶青老师、北大国发院刘子毅同学在我们写作中对相关问题的答疑。
Note：产生如下推文列表的 Stata 命令为：
lianxh 因果推断
安装最新版lianxh命令：
ssc install lianxh, replace
课程推荐：因果推断实用计量方法
主讲老师：邱嘉平教授
🍓 课程主页：https://gitee.com/lianxh/YGqjp
New！ Stata 搜索神器：
lianxh和songblGIF 动图介绍
搜： 推文、数据分享、期刊论文、重现代码 ……
👉 安装：
. ssc install lianxh
. ssc install songbl
👉 使用：
. lianxh DID 倍分法
. songbl all
🍏 关于我们
联系客服
