如何在R语言中实现多期DID

发布时间:2021-12-27 阅读 4068

Stata连享会   主页 || 视频 || 推文 || 知乎 || Bilibili 站

温馨提示: 定期 清理浏览器缓存,可以获得最佳浏览体验。

New! lianxh 命令发布了:
随时搜索推文、Stata 资源。安装:
. ssc install lianxh
详情参见帮助文件 (有惊喜):
. help lianxh
连享会新命令:cnssc, ihelp, rdbalance, gitee, installpkg

课程详情 https://gitee.com/lianxh/Course

课程主页 https://gitee.com/lianxh/Course

⛳ Stata 系列推文:

PDF下载 - 推文合集

作者:张梓瑶 (南开大学)
邮箱nk_zhangziyao@163.com


目录


1. 简介

双重差分模型 (Difference in Difference,DID) 是政策评估的常用方法。其中,经典的 DID 方法通常拥有两期 (处理前和处理后) 和两组 (处理组和对照组) 数据。在具体应用中,考虑到个体和时间异质性,我们常用 DID 固定效应模型对平均处理效应 (ATE) 进行估计。

但这种经典的 DID 方法在实际应用中往往面临挑战。例如,当存在多个政策时点冲击时,DID 固定效应模型对处理效应的估计时可能会产生严重偏误。为此,本文将介绍一个基于 R 语言的 did 包来解决上述问题,以及进行平行趋势检验。

2. 理论背景

2.1 两期 DID 固定效应模型的问题

基本的两期 DID 设定模型为:

其中, Gi 为分组虚拟变量 (处理组=1,控制组=0),Dt 为分期虚拟变量 (政策实施后=1,政策实施前=0),交互项 Gi×Dt 系数 β 为 DID 模型重点考察的平均处理效应 (Average Treatment Effect,ATE)。

上述基本的 DID 模型本质上是一种固定效应模型,常数项 α 是对照组的固定效应,γ 是处理组与对照组的固定效应差异,δ 是事件发生前后的固定效应差异,β 是事件固定效应 (处置效应)。

当然,如果我们有更多时间信息的面板数据,就可以通过细化组别和时期的固定效应来提高上述模型的精度。具体地,将处理组和对照组的固定效应 (α,γ) 细化为个体效应 μi,以及将事件前后的时期固定效应 δ 细化为每年的固定效应 λt,细化后的模型则可以较好地降低估计系数的方差。

两期 DID 的固定效应模型设定如下:

但两期 DID 固定效应模型估计可能会受到多期政策冲击的挑战。由于多期 DID 政策冲击时点不同因而存在多个不同的对照组,比如后接受政策冲击的处理组可能是先接受政策冲击的对照组,故存在一种动态的平均处理效应估计,导致估计偏误。

那么我们如何识别这种动态的平均处理效应呢?

2.2 多期 DID 平行趋势假设

DID 进行估计的前提是平行趋势假设,即处理组和对照组在政策冲击之前应具有相同的趋势。两期 DID 平行趋势假设条件为:

多期 DID 的平行趋势假设条件为:对 g=2,...,τt=2,...,τ,当 gt 时,

考虑控制变量后的多期 DID 平行趋势假设条件为:对 g=2,...,τt=2,...,τ,当 gt 时,

在上式中, Yit(0) 表示未接受处理的个体 i 在 t 时期的潜在结果;Yit(g) 表示个体 i 在 t 时期的潜在结果,该个体在 g 时期接受政策冲击;Gi 表示个体 i 接受处理的时间段;Ci 表示个体是否为对照组 (对照组随 t 变化,为简单起见,这里将从未接受过处理的个体作为对照组);Yit 表示个体 i 在 t 时期观察到的结果;Xi 表示参与回归的控制变量。

2.3 处理效应估计

ATE (Average Treatment Effects) = E[Yit(1)Yit(0)] 指的是从样本中随机抽取某个体的期望处理效应,即平均处理效应;ATT (Average Treatment Effects on the Treated) = E[Yit(1)Yit(0)|Di=1] 仅考虑处理组的平均处理效应,即处理组处理效应。

在随机分组的条件下 Di 独立于 [Yit(1),Yit(0)],此时 ATE = ATT。由于 ATE 既包括了处理组也包括了对照组个体,很多学者批评其定义过于宽泛。而对于政策制定者来说,ATT 可能是更加重要的考虑因素,因为他直接衡量了处理组中个体接受政策冲击前后的收益。

多期 DID 下的 ATT 的计算公式为:

其中,当 g=2,t=2 时,ATT(g=2,t=2) 表示两期 DID 处理组的平均处理效应。此外,考虑多期 DID 情况下,ATT(g=2,t=3) 表示个体在第 2 期接受处理在 t=3 时的处理效应;ATT(g=3,t=3) 表示个体在第 3 期接受处理在 t=3 时的处理效应。

下面我们表示出平均处理效应 (ATE) 的估计公式:

上式表示每组分别接受政策冲击后的 ATT;

上式表示所有处理组个体在参与处理后的总体处理效应;

上式中的 e 表示处理组已经暴露在政策冲击下的时间,θD(e) 可以看作一种动态效应估计指标。

3. R 语言实操

3.1 命令简介

did 是 Callaway 和 Sant’Anna 创建的用于多期 DID 处理效应估计的 R 命令包,其主要操作命令包括:

命令 功能
att_gt() Group-Time Average Treatment Effects
aggte() Aggregate Group-Time Average Treatment Effects
conditional_did_pretest() Pre-Test of Conditional Parallel Trends Assumption

关于上述命令的详细介绍,可点击对应命令查看。另外 did 包的其他命令介绍,请参考「Difference-in-Differences Methods」

did 包用于估计处理组处理效应 (ATT) 具有以下几个特点:

  • 考虑了多时点政策冲击的处理效应;
  • 考虑了在存在控制变量下的平行趋势假设条件;
  • 当存在多时点政策冲击时,ATT 估计避免了采用双向固定效应回归可能存在的估计偏误。

3.2 R 程序实操

下面我们用 R 语言来演示 did 包的使用。数据来自 Callaway 和 Sant’Anna (2020),研究问题是各州提高最低工资对县级青少年就业率影响。

样例数据集包含 2003 年至 2007 年 500 个县级青少年就业率的数据,其中一些州在 2004 年首次接受治疗,也有一些在 2006 年或 2007 年接受治疗。

主要变量 lemp 是县级青少年就业率的对数值,first.treat 是某个州首次提高最低工资的时间,可能是 2004 年、2006 年或 2007 年;year 是年份; countryreal 是每个县的 ID。

首先,在 R 中安装did 包。

writeLines('PATH="${RTOOLS40_HOME}\\usr\\bin;${PATH}"', con = "~/.Renviron")  #创建路径配置文件.Renviron
Sys.which("make")  #测试配置路径是否成功
options(repos=structure(c(CRAN="https://mirrors.tuna.tsinghua.edu.cn/CRAN/"))) #设置镜像
install.packages("did") #从CRAN上下载did包
library(did)  #调用did包
data(mpdta)   #调入样例数据

did 包使用 att_ gt 命令估计处理组平均处理效应 (ATT)。

out <- att_gt(yname = "lemp",   
              gname = "first.treat",  
              idname = "countyreal",
              tname = "year",
              xformla = ~1,
              data = mpdta,
              est_method = "reg"
              )
summary(out)

#> Call:
#> att_gt(yname = "lemp", tname = "year", idname = "countyreal", 
#>     gname = "first.treat", xformla = ~1, data = mpdta, est_method = "reg")
#> 
#> Group-Time Average Treatment Effects:
#>  Group Time ATT(g,t) Std. Error [95% Simult.  Conf. Band]  
#>   2004 2004  -0.0105     0.0259       -0.0796      0.0586  
#>   2004 2005  -0.0704     0.0312       -0.1535      0.0127  
#>   2004 2006  -0.1373     0.0375       -0.2373     -0.0372 *
#>   2004 2007  -0.1008     0.0352       -0.1947     -0.0069 *
#>   2006 2004   0.0065     0.0229       -0.0547      0.0677  
#>   2006 2005  -0.0028     0.0203       -0.0568      0.0513  
#>   2006 2006  -0.0046     0.0179       -0.0522      0.0430  
#>   2006 2007  -0.0412     0.0210       -0.0972      0.0148  
#>   2007 2004   0.0305     0.0159       -0.0118      0.0728  
#>   2007 2005  -0.0027     0.0164       -0.0465      0.0410  
#>   2007 2006  -0.0311     0.0183       -0.0797      0.0176  
#>   2007 2007  -0.0261     0.0179       -0.0737      0.0215  
#> ---
#> Signif. codes: `*' confidence band does not cover 0
#> 
#> P-value for pre-test of parallel trends assumption:  0.16812
#> Control Group:  Never Treated,  Anticipation Periods:  0
#> Estimation Method:  Outcome Regression

从上面 summary 的输出结果可以看出,该命令不仅估计出 gt 的平均处理效应,也估计出了 gt 的平均处理效应用以检验平行趋势假设。这里我们使用 Wlad 检验的 p 值统计量来进行平行趋势检验。

多期 DID 的平均处理效应的图示命令如下:

ggdid(out, ylim = c(-.25,.1))

上图中的红色圆点是政策冲击前的平均治疗效果 (采用了 95% 的置信区间),可以看出红色圆点在 95% 的置信水平下趋近于 0,符合平行趋势假设。蓝色圆点表示接受政策冲击后的平均治疗效果,在本案例中表示提高最低工资水平对县级青少年就业率的影响。

es <- aggte(out, type = "dynamic")
summary(es)

#> Call:
#> aggte(MP = out, type = "dynamic")
#> 
#> Overall ATT:  
#>      ATT Std. Error     [95%  Conf. Int.]  
#>  -0.0772     0.0205   -0.1175      -0.037 *
#>  
#> Dynamic Effects:
#>  Event time Estimate Std. Error [95% Simult.  Conf. Band]  
#>          -3   0.0305     0.0153       -0.0102      0.0712  
#>          -2  -0.0006     0.0132       -0.0356      0.0345  
#>          -1  -0.0245     0.0145       -0.0630      0.0141  
#>           0  -0.0199     0.0120       -0.0518      0.0120  
#>           1  -0.0510     0.0175       -0.0977     -0.0042 *
#>           2  -0.1373     0.0396       -0.2426     -0.0319 *
#>           3  -0.1008     0.0346       -0.1929     -0.0087 *
#> ---
#> Signif. codes: `*' confidence band does not cover 0
#> 
#> Control Group:  Never Treated,  Anticipation Periods:  0
#> Estimation Method:  Outcome Regression

did 包允许使用 aggte 估计动态处理效应,命令如下:

es <- aggte(out, type = "dynamic")
summary(es)

#> Call:
#> aggte(MP = out, type = "dynamic")
#> 
#> Overall ATT:  
#>      ATT Std. Error     [95%  Conf. Int.]  
#>  -0.0772     0.0205   -0.1175      -0.037 *
#> 
#> Dynamic Effects:
#>  Event time Estimate Std. Error [95% Simult.  Conf. Band]  
#>          -3   0.0305     0.0153       -0.0102      0.0712  
#>          -2  -0.0006     0.0132       -0.0356      0.0345  
#>          -1  -0.0245     0.0145       -0.0630      0.0141  
#>           0  -0.0199     0.0120       -0.0518      0.0120  
#>           1  -0.0510     0.0175       -0.0977     -0.0042 *
#>           2  -0.1373     0.0396       -0.2426     -0.0319 *
#>           3  -0.1008     0.0346       -0.1929     -0.0087 *
#> ---
#> Signif. codes: `*' confidence band does not cover 0
#> 
#> Control Group:  Never Treated,  Anticipation Periods:  0
#> Estimation Method:  Outcome Regression

上面 summary 的输出结果中,event time 表示个体第一次接受政策冲击时的处理效应,比如 event time=0 表示个体刚好处于政策冲击期的处理效应,event time=-1 表示个体在接受政策冲击前一年的处理效应。

上图与多期 DID 平均处理效应的图相似,但这里将政策冲击看作一个事件,横坐标轴表示事件发生时间。通过对每段时期的处理效应进行估计,上图表示个体参与治疗的总体处理效应。

总体平均处理效应 (ATE) 的计算程序如下:

group_effects <- aggte(out, type = "group")
summary(group_effects)

#> Call:
#> aggte(MP = out, type = "group")
#> 
#> Overall ATT:  
#>     ATT Std. Error     [95%  Conf. Int.]  
#>  -0.031     0.0125   -0.0554     -0.0066 *
#> 
#> Group Effects:
#>  Group Estimate Std. Error [95% Simult.  Conf. Band]  
#>   2004  -0.0797     0.0274       -0.1388     -0.0207 *
#>   2006  -0.0229     0.0165       -0.0585      0.0127  
#>   2007  -0.0261     0.0176       -0.0641      0.0120  
#> ---
#> Signif. codes: `*' confidence band does not cover 0
#> 
#> Control Group:  Never Treated,  Anticipation Periods:  0
#> Estimation Method:  Outcome Regression

从总体平均处理效应 (ATE) 的结果来看,提高最低工资会使青少年就业率下降 3.1%,在 10% 的置信水平下显著。

4. 参考资料

  • Callaway B, Sant’Anna P H C. Difference-in-differences with multiple time periods[J]. Journal of Econometrics, 2021, 225(2): 200-230. -PDF-
  • Difference-in-Differences -Link-
  • Difference-in-Differences Methodology -Link-
  • 邱嘉平. 因果推断实用计量方法[M]. 上海:上海财经大学出版社,2020.
  • 陈强. 高级计量经济学及 Stata 应用 (第二版)[M]. 教育出版社,2014.

5. 相关推文

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

相关课程

免费公开课

最新课程-直播课

专题 嘉宾 直播/回看视频
最新专题 文本分析、机器学习、效率专题、生存分析等
研究设计 连玉君 我的特斯拉-实证研究设计-幻灯片-
面板模型 连玉君 动态面板模型-幻灯片-
面板模型 连玉君 直击面板数据模型 [免费公开课,2小时]
  • Note: 部分课程的资料,PPT 等可以前往 连享会-直播课 主页查看,下载。

课程主页

课程主页

关于我们

  • Stata连享会 由中山大学连玉君老师团队创办,定期分享实证分析经验。
  • 连享会-主页知乎专栏,700+ 推文,实证分析不再抓狂。直播间 有很多视频课程,可以随时观看。
  • 公众号关键词搜索/回复 功能已经上线。大家可以在公众号左下角点击键盘图标,输入简要关键词,以便快速呈现历史推文,获取工具软件和数据下载。常见关键词:课程, 直播, 视频, 客服, 模型设定, 研究设计, stata, plus, 绘图, 编程, 面板, 论文重现, 可视化, RDD, DID, PSM, 合成控制法

连享会小程序:扫一扫,看推文,看视频……

扫码加入连享会微信群,提问交流更方便

✏ 连享会-常见问题解答:
https://gitee.com/lianxh/Course/wikis

New! lianxhsongbl 命令发布了:
随时搜索连享会推文、Stata 资源,安装命令如下:
. ssc install lianxh
使用详情参见帮助文件 (有惊喜):
. help lianxh