查看原文
其他

DID偏误问题:多时期DID的双重稳健估计量(下)-csdid

连享会 连享会 2023-10-24

👇 连享会 · 推文导航 | www.lianxh.cn

连享会 · 因果推断实用计量方法

作者:张子楠 (浙江财经大学)
邮箱: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-


目录

  • 1. 绪论

  • 2. csdid 命令介绍

  • 3. csdid 基本用法

  • 4. csdid 进阶用法

  • 5. 进一步讨论

  • 6. 相关推文



1. 绪论

在连享会推文「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) 同时分析了面板数据和混合截面数据两种数据结构下的双重稳健估计问题。本推文只介绍了面板数据部分,对混合截面数据相关内容感兴趣的读者可以自行阅读原文。

2. csdid 命令介绍

* 命令的安装
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 组样本。默认是不包括后者;
  • longlong2 对于冲击发生前的样本期的样本,使用 long gaps 而非 shortgaps;
  • 估计方法的选项包括:增进型双重稳健估计 (drimp)、基于 ipw 的双重稳健估计量 (dripw)、outcome regression 估计量 (reg)、标准化 ipw 估计量 (stdipw)、ipw 估计量 (ipw);
  • 标准误的选项默认是 robust and asymptotic standard errors,同时也可以选择 wbootcluster 等;
  • 不同组别政策效应加总 agg(aggtype)。加总方式 (aggtype) 包括:加总所有组所有时期 (simple)、每一组或者一群在所有时期加总 (group)、所有组在每一个时期分别加总 (calendar)、采用事件发生法在每一个时期分别加总 (event) 四大类。

3. csdid 基本用法

我们采用作者在帮助文档里提供的数据 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 分解方法来理解背后原因。

4. csdid 进阶用法

除了更加稳健,以及能有更大概率获得无偏估计结果,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

5. 进一步讨论

在帮助文件里,作者还给出了一些注意事项,主要有以下三点:

  • drdid 命令类似,csdid 假定协变量不随时间变化。当使用面板数据回归时,即使协变量是时变的,csdid 命令也只使用协变量基期值来估计;
  • 要使用 csdid 命令来获得 DID 无偏估计量,一般只能使用 never-treated 和 (或) not-yet-treated 的数据作为控制组,而不能使用 always-treated 组。否则由于异质性处理效应,平行趋势假设并不能满足,从而估计结果仍旧有偏;
  • 当使用面板数据时,与 drdid 命令不同,csdid 命令并没有要求数据为强平衡面板数据,即并没有要求没有缺失值。但是,在估计不同组的 ATT 时,存在缺失值的样本并不会进入回归。这和 TWFE 里对于缺失值的处理方式是一致的。

6. 相关推文

Note:产生如下推文列表的 Stata 命令为:
lianxh did, m
安装最新版 lianxh 命令:
ssc install lianxh, replace

  • 专题:倍分法DID
    • Stata+R:合成DID原理及实现-sdid
    • DID的陷阱和注意事项
    • Stata:事件研究法的稳健有效估计量-did_imputation
    • DID最新进展:异质性处理条件下的双向固定效应DID估计量 (TWFEDD)
    • 公开课:双重差分方法的新进展——交错型DID
    • 倍分法:交错DID与Stata操作
    • Stata倍分法新趋势:did2s-两阶段双重差分模型
    • Stata-DID:不同处理时点不同持久期的倍分法(flexpaneldid)
    • DIDM:多期多个体倍分法-did_multiplegt
    • Stata:双重差分的固定效应模型-(DID)
    • 多期DID之安慰剂检验、平行趋势检验
    • DID边际分析:让政策评价结果更加丰满
    • Big Bad Banks:多期 DID 经典论文介绍
    • Stata:多期倍分法 (DID) 详解及其图示
    • 倍分法DID详解 (二):多时点 DID (渐进DID)
    • 倍分法DID详解 (三):多时点 DID (渐进DID) 的进一步分析

课程推荐:因果推断实用计量方法
主讲老师:丘嘉平教授
🍓 课程主页https://gitee.com/lianxh/YGqjp

New! Stata 搜索神器:lianxhsongbl  GIF 动图介绍
搜: 推文、数据分享、期刊论文、重现代码 ……
👉 安装:
. ssc install lianxh
. ssc install songbl
👉  使用:
. lianxh DID 倍分法
. songbl all

🍏 关于我们

  • 连享会 ( www.lianxh.cn,推文列表) 由中山大学连玉君老师团队创办,定期分享实证分析经验。
  • 直通车: 👉【百度一下: 连享会】即可直达连享会主页。亦可进一步添加 「知乎」,「b 站」,「面板数据」,「公开课」 等关键词细化搜索。


您可能也对以下帖子感兴趣

文章有问题?点此查看未经处理的缓存