Stata连享会 主页 || 视频 || 推文 || 知乎 || Bilibili 站
温馨提示: 定期 清理浏览器缓存,可以获得最佳浏览体验。
New!
lianxh
命令发布了:
随时搜索推文、Stata 资源。安装:
. ssc install lianxh
详情参见帮助文件 (有惊喜):
. help lianxh
连享会新命令:cnssc
,ihelp
,rdbalance
,gitee
,installpkg
⛳ Stata 系列推文:
作者:张子楠 (浙江财经大学)
邮箱:zinanzh@gmail.com
编者按:本文主要摘译自下文,特此致谢!
Source:Callaway B, Sant’Anna P H C. Difference-in-differences with multiple time periods[J]. Journal of Econometrics, 2021, 225(2): 200-230. -PDF-
目录
在连享会推文「DID偏误问题:两时期DID的双重稳健估计量(上)-drdid」中,我们介绍了对于两时期 DID 模型,双重固定效应 (TWFE) 估计 DID 处理效应可能存在的偏误问题,以及如何通过双重稳健估计方法来缓解。此外,我们还介绍了如何通过 Sant 和 Zhao(2020) 提供的命令 drdid
来实现上述过程。
但 drdid
命令有两个不足:一是只适用两期 DID 情形,无法直接拓展到多期 DID 中;二是要求面板数据为强平衡数据,不能有缺失值。而这两点要求,在实践应用中都往往难以满足。特别是第一点,多时期 DID 要远比两时期 DID 更为常见。为了缓解上述问题,Callaway 和 Sant (2021) 提供了 drdid
命令的多时期版本——csdid
。
与其它多时期 DID 估计的纠偏估计命令类似,csdid
命令也是将样本分为不同的子组 (group),分别估计不同组别的 ATT(g),再通过特定策略将不同组别的 ATT(g) 加总算出样本期的 ATT。加总策略原则同样也是对那些可能存在偏误组的 ATT(g),降低它们的加总权重。
与其它纠偏估计量不同的是,csdid
是基于 drdid
的双重稳健估计思路,避免了 TWFE 估计量的问题。同时只需要 Outcome Regression 或 Inverse Propesnity Weight Regression 两种估计方法的依赖条件至少有一个满足,就可以显著缓解 DID 估计偏误问题。
csdid
命令还有一个优势在于,提供了多种加总不同子组 DID 回归结果的选项。除了常见的利用事件发生法去分析政策的动态效应以外,命令还允许我们计算不同子组处理效应的异质性,不同年份处理效应的异质性,以及处理效应的在一段时期的累计值。这不仅丰富了 DID 回归结果,还可以更为深入地了解 DID 总平均处理效应的底层来源。
注:Callaway 和 Sant Anna(2021) 同时分析了面板数据和混合截面数据两种数据结构下的双重稳健估计问题。本推文只介绍了面板数据部分,对混合截面数据相关内容感兴趣的读者可以自行阅读原文。
* 命令的安装
ssc install csdid, all replace
* 命令的语法
csdid depvar [indepvars] [if] [in] [weight], [ivar(varname)] time(varname) gvar(varname) [ options ]
对于面板数据,模型设定的选项包括以下四类:
depvar
为被解释变量,indepvars
为解释变量和控制变量;ivar(varname)
和 time (varname)
为面板设置选项。其中 ivar(varname)
设置个体标识变量,time (varname)
设置时间标识变量;gvar(varname)
为多时期 DID 标识变量,和通常的多时期 DID 赋值不太一样。如果在样本期受到过政策冲击,赋值等于其受到冲击的期数序号。例如,在第 1 期受到冲击,赋值为 1。需要注意的是,如果某个个体在样本期开始前就已经受到冲击,即属于所谓的 always-treated 组,默认是不进入回归中;notyet
要求对照组里不仅要有 not-yet 组里尚未受到冲击的样本,同时也要包括 always-treated 组样本。默认是不包括后者;long
和 long2
对于冲击发生前的样本期的样本,使用 long gaps 而非 shortgaps;drimp
)、基于 ipw 的双重稳健估计量 (dripw
)、outcome regression 估计量 (reg
)、标准化 ipw 估计量 (stdipw
)、ipw 估计量 (ipw
);robust and asymptotic standard errors
,同时也可以选择 wboot
、cluster
等;agg(aggtype)
。加总方式 (aggtype
) 包括:加总所有组所有时期 (simple
)、每一组或者一群在所有时期加总 (group
)、所有组在每一个时期分别加总 (calendar
)、采用事件发生法在每一个时期分别加总 (event
) 四大类。
我们采用作者在帮助文档里提供的数据 mpdta.dta 来演示 csdid
的用法。
. use https://friosavila.github.io/playingwithstata/drdid/mpdta.dta, clear
. des
Contains data from https://friosavila.github.io/playingwithstata/drdid/mpdta.dta
Observations: 2,500 Written by R.
Variables: 6 17 May 2021 11:45
-------------------------------------------------------------------------------------
Variable Storage Display Value
name type format label Variable label
-------------------------------------------------------------------------------------
year int %9.0g year
countyreal long %9.0g countyreal
lpop double %9.0g lpop
lemp double %9.0g lemp
first_treat int %9.0g first.treat
treat byte %9.0g treat
-------------------------------------------------------------------------------------
. tab year
year | Freq. Percent Cum.
------------+-----------------------------------
2003 | 500 20.00 20.00
2004 | 500 20.00 40.00
2005 | 500 20.00 60.00
2006 | 500 20.00 80.00
2007 | 500 20.00 100.00
------------+-----------------------------------
Total | 2,500 100.00
. tab first_treat treat
first.trea | treat
t | 0 1 | Total
-----------+----------------------+----------
0 | 1,545 0 | 1,545
2004 | 0 100 | 100
2006 | 0 200 | 200
2007 | 0 655 | 655
-----------+----------------------+----------
Total | 1,545 955 | 2,500
数据共包含 6 个变量,其中 year 表示年份、countyreal 表示个体、lpop 为被解释变量、lemp 为协变量。treat 为处理组指示变量,如果样本期内受到冲击,则为 1,否则为 0。first_treat 为个体第一次受到冲击的年份,取值包括 0,2004,2006 和 2007,0 表示样本期内没有受到冲击。观察上面结果可知,样本量为 2500,样本时期为 2003-2007,共五年。总共有 2500/5 个个体,其中分别有 100/5、200/5 和 655/5 数量的个体在 2004、2006 和 2007 年受到冲击。
接下来,我们用 csdid
命令来识别处理组平均处理效应 (ATT)。
. * 需要使用 agg(simple) 选项, 以估计总 ATT
. * method(dripw) 表示使用增进型双重稳健估计
. csdid lemp lpop , ivar(countyreal) time(year) gvar(first_treat) method(dripw) agg(simple)
Outcome model : least squares Number of obs = 2,500
Treatment model: inverse probability
------------------------------------------------------------------------------
| Coefficient Std. err. z P>|z| [95% conf. interval]
-------------+----------------------------------------------------------------
ATT | -0.042 0.012 -3.63 0.000 -0.064 -0.019
------------------------------------------------------------------------------
Control: Never Treated
See Callaway and Sant'Anna (2021) for details
此外,我们也用 TWFE 方法识别 ATT 来作为对比。
. * 生成 did 选项: 受到冲击的当期及以后为 1, 其它情况为 0
. gen did=0
. replace did=1 if first_treat!=0 & year>= first_treat
. * TWFE 回归
. reghdfe lemp did lpop, absorb(countyreal year) cluster(countyreal)
HDFE Linear regression Number of obs = 2,500
Absorbing 2 HDFE groups F( 1, 499) = 7.59
Statistics robust to heteroskedasticity Prob > F = 0.0061
R-squared = 0.9932
Adj R-squared = 0.9915
Within R-sq. = 0.0042
Number of clusters (countyreal) = 500 Root MSE = 0.1391
(Std. err. adjusted for 500 clusters in countyreal)
------------------------------------------------------------------------------
| Robust
lemp | Coefficient std. err. t P>|t| [95% conf. interval]
-------------+----------------------------------------------------------------
did | -0.037 0.013 -2.76 0.006 -0.063 -0.010
lpop | 0.000 (omitted)
_cons | 5.777 0.002 3741.28 0.000 5.774 5.780
------------------------------------------------------------------------------
Absorbed degrees of freedom:
-----------------------------------------------------+
Absorbed FE | Categories - Redundant = Num. Coefs |
-------------+---------------------------------------|
countyreal | 500 500 0 *|
year | 5 1 4 |
-----------------------------------------------------+
* = FE nested within cluster; treated as redundant for DoF computation
对比可以发现,两种回归结果的回归系数和标准误都比较接近,置信区间基本重合。这说明对于这个例子,TWFE 估计偏误问题并不严重。注:有兴趣的读者可以用 Bacon 分解方法来理解背后原因。
除了更加稳健,以及能有更大概率获得无偏估计结果,csdid
命令的另一大优势在于,可以更为灵活地划分子组,并加总获得不同特征的 ATT,这在异质性分析等方面极大丰富 DID 估计结果。对这部分细节感兴趣的读者,可以阅读 Callaway 和 Sant Anna (2021) 的第 3 章,对此有详细的叙述。
接下来我们分别演示常见的几种加总获得样本总 ATT 的思路。第一、利用事件发生法的思路,以接受冲击时期的长短来来区分不同组别,估计事件冲击的动态效应。
. * 需使用 agg(event) 选项来实现
. csdid lemp lpop, ivar(countyreal) time(year) gvar(first_treat) method(dripw) agg(event)
Outcome model : least squares Number of obs = 2,500
Treatment model: inverse probability
------------------------------------------------------------------------------
| Coefficient Std. err. z P>|z| [95% conf. interval]
-------------+----------------------------------------------------------------
T-3 | 0.027 0.014 1.90 0.057 -0.001 0.054
T-2 | -0.004 0.013 -0.28 0.780 -0.029 0.022
T-1 | -0.023 0.014 -1.60 0.109 -0.052 0.005
T+0 | -0.021 0.011 -1.83 0.067 -0.044 0.001
T+1 | -0.053 0.016 -3.24 0.001 -0.085 -0.021
T+2 | -0.140 0.035 -3.97 0.000 -0.210 -0.071
T+3 | -0.107 0.033 -3.25 0.001 -0.171 -0.042
------------------------------------------------------------------------------
Control: Never Treated
See Callaway and Sant'Anna (2021) for details
类似地,我们也用 TWFE 估计方法进行了对比,如下所示。由结果可知,对于我们这个例子,两种估计量的估计结果基本相一致。
. * 生成事件发生发所需要的变量
. gen didlength=year-first_treat
. replace didlength=0 if first_treat==0
. tab didlength,gen(event)
. drop event1 // 删除第一个 event 变量, 以和 csdid 命令处理策略一致
. * TWFE 回归
. reghdfe lemp event* lpop, absorb(countyreal year) cluster(countyreal)
HDFE Linear regression Number of obs = 2,500
Absorbing 2 HDFE groups F( 7, 499) = 3.60
Statistics robust to heteroskedasticity Prob > F = 0.0009
R-squared = 0.9933
Adj R-squared = 0.9915
Within R-sq. = 0.0103
Number of clusters (countyreal) = 500 Root MSE = 0.1388
(Std. err. adjusted for 500 clusters in countyreal)
------------------------------------------------------------------------------
| Robust
lemp | Coefficient std. err. t P>|t| [95% conf. interval]
-------------+----------------------------------------------------------------
event2 | 0.021 0.015 1.39 0.166 -0.009 0.051
event3 | 0.020 0.019 1.02 0.308 -0.018 0.058
event4 | -0.004 0.023 -0.16 0.877 -0.048 0.041
event5 | -0.022 0.025 -0.85 0.393 -0.072 0.028
event6 | -0.047 0.029 -1.63 0.103 -0.104 0.010
event7 | -0.135 0.039 -3.44 0.001 -0.213 -0.058
event8 | -0.096 0.042 -2.27 0.024 -0.179 -0.013
lpop | 0.000 (omitted)
_cons | 5.788 0.022 263.18 0.000 5.745 5.831
------------------------------------------------------------------------------
Absorbed degrees of freedom:
-----------------------------------------------------+
Absorbed FE | Categories - Redundant = Num. Coefs |
-------------+---------------------------------------|
countyreal | 500 500 0 *|
year | 5 1 4 |
-----------------------------------------------------+
* = FE nested within cluster; treated as redundant for DoF computation
第二、计算不同年份的所有组别平均效应。即对于某一年份的所有样本,其年平均处理效应是多少。需要注意的是,事件发生法是以受冲击的时期长度分组,这里是以某年份来区分不同的组别。
. * 需使用 agg(calendar) 选项来实现
. csdid lemp lpop, ivar(countyreal) time(year) gvar(first_treat) method(dripw) agg(calendar)
Outcome model : least squares Number of obs = 2,500
Treatment model: inverse probability
------------------------------------------------------------------------------
| Coefficient Std. err. z P>|z| [95% conf. interval]
-------------+----------------------------------------------------------------
T2004 | -0.015 0.022 -0.66 0.511 -0.058 0.029
T2005 | -0.076 0.029 -2.67 0.008 -0.133 -0.020
T2006 | -0.046 0.021 -2.18 0.029 -0.088 -0.005
T2007 | -0.040 0.013 -3.06 0.002 -0.065 -0.014
------------------------------------------------------------------------------
Control: Never Treated
See Callaway and Sant'Anna (2021) for details
第三、计算不同组别的同一年份平均效应。在本示例里,事件冲击先后发生在 2004、2006 和 2007 这三年,因而有三个子组,分别用 G2004、G2006和 G2007 来标识。
. * 需使用 agg(group) 选项来实现
. csdid lemp lpop, ivar(countyreal) time(year) gvar(first_treat) method(dripw) agg(group)
Outcome model : least squares Number of obs = 2,500
Treatment model: inverse probability
------------------------------------------------------------------------------
| Coefficient Std. err. z P>|z| [95% conf. interval]
-------------+----------------------------------------------------------------
G2004 | -0.085 0.025 -3.44 0.001 -0.133 -0.036
G2006 | -0.020 0.017 -1.15 0.248 -0.054 0.014
G2007 | -0.029 0.016 -1.77 0.076 -0.061 0.003
------------------------------------------------------------------------------
Control: Never Treated
See Callaway and Sant'Anna (2021) for details
第四、不仅使用 Never Treated,也使用 Not yet Treated 作为控制组。需要注意的是,虽然回归结果表格下显示的是 "Control: Not yet Treated"。但根据帮助文件里选项设定的描述来看,这里面应该是同时使用了 Never Treated 和 Not yet Treated 来作为控制组的。
. * 需使用 notyet 选项来实现
. csdid lemp lpop, ivar(countyreal) time(year) gvar(first_treat) method(dripw) agg(simple) notyet
Outcome model : least squares Number of obs = 2,500
Treatment model: inverse probability
------------------------------------------------------------------------------
| Coefficient Std. err. z P>|z| [95% conf. interval]
-------------+----------------------------------------------------------------
ATT | -0.041 0.011 -3.62 0.000 -0.064 -0.019
------------------------------------------------------------------------------
Control: Not yet Treated
See Callaway and Sant'Anna (2021) for details
第五、对于非平衡面板数据进行回归。在 drdid
命令中,是无法使用非平衡面板数据回归的。csdid
则在这方面拓展了双重稳健估计量在 DID 分析中的应用。
. * 生成随机数 sample, 赋值为 0 或者 1, 均值为 0.9
. set seed 1
. gen sample = runiform()<.9
. * 需使用 if sample==1 选项来挑选部分样本回归, 以实现类似存在缺失值的效果
. csdid lemp lpop if sample==1, cluster(countyreal) time(year) gvar(first_treat) method(dripw) agg(simple)
Outcome model : least squares Number of obs = 2,245
Treatment model: inverse probability
------------------------------------------------------------------------------
| Coefficient Std. err. z P>|z| [95% conf. interval]
-------------+----------------------------------------------------------------
ATT | -0.052 0.018 -2.86 0.004 -0.087 -0.016
------------------------------------------------------------------------------
Control: Never Treated
See Callaway and Sant'Anna (2021) for details
在帮助文件里,作者还给出了一些注意事项,主要有以下三点:
drdid
命令类似,csdid
假定协变量不随时间变化。当使用面板数据回归时,即使协变量是时变的,csdid
命令也只使用协变量基期值来估计;csdid
命令来获得 DID 无偏估计量,一般只能使用 never-treated 和 (或) not-yet-treated 的数据作为控制组,而不能使用 always-treated 组。否则由于异质性处理效应,平行趋势假设并不能满足,从而估计结果仍旧有偏;drdid
命令不同,csdid
命令并没有要求数据为强平衡面板数据,即并没有要求没有缺失值。但是,在估计不同组的 ATT 时,存在缺失值的样本并不会进入回归。这和 TWFE 里对于缺失值的处理方式是一致的。
Note:产生如下推文列表的 Stata 命令为:
lianxh did, m
安装最新版lianxh
命令:
ssc install lianxh, replace
免费公开课
最新课程-直播课
专题 | 嘉宾 | 直播/回看视频 |
---|---|---|
⭐ 最新专题 | 文本分析、机器学习、效率专题、生存分析等 | |
研究设计 | 连玉君 | 我的特斯拉-实证研究设计,-幻灯片- |
面板模型 | 连玉君 | 动态面板模型,-幻灯片- |
面板模型 | 连玉君 | 直击面板数据模型 [免费公开课,2小时] |
⛳ 课程主页
⛳ 课程主页
关于我们
课程, 直播, 视频, 客服, 模型设定, 研究设计, stata, plus, 绘图, 编程, 面板, 论文重现, 可视化, RDD, DID, PSM, 合成控制法
等
连享会小程序:扫一扫,看推文,看视频……
扫码加入连享会微信群,提问交流更方便
✏ 连享会-常见问题解答:
✨ https://gitee.com/lianxh/Course/wikis
New!
lianxh
和songbl
命令发布了:
随时搜索连享会推文、Stata 资源,安装命令如下:
. ssc install lianxh
使用详情参见帮助文件 (有惊喜):
. help lianxh