查看原文
其他

双重差分法(DID)安慰剂检验的做法:随机抽取500次?

江河JH 功夫计量经济学 2022-12-31
“安慰剂”(placebo)一词来自医学上的随机实验,比如要检验某种新药的疗效。此时,可将参加实验的人群随机分为两组,其中一组为实验组,服用真药;而另一组为控制组,服用安慰剂(比如,无用的糖丸),并且不让参与者知道自己服用的究竟是真药还是安慰剂,以避免由于主观心理作用而影响实验效果。

双重差分法(DID)安慰剂检验的核心思想就是虚构处理组或者虚构政策时间进行估计,如果虚构情况下“伪政策虚拟变量”的系数依然显著,那么就说明原来的估计结果很有可能出现了偏误,我们的被解释变量y的变动很有可能是受到了其他政策或者随机性因素的影响。

说到虚构,那么自然是可以随机虚构,也可以不随机虚构(作者自己设定)。当然,我更推荐的还是随机虚构处理组或者是政策时间的方法。由于我们使用的数据基本都是“大N小T”型的短面板数据,所以随机虚构政策时间没什么意义,文献一般做法都是将政策年份统一提前2年或3年重新进行回归,看看政策虚拟变量系数是否依然显著。我们更多地还是随机虚构处理组,具体做法就是随机选取个体作为处理组,重复500次或者1000次,看看“伪政策虚拟变量”的系数是否显著

数据来源

石大千等(2018)发表在《中国工业经济》的论文《智慧城市建设能否降低环境污染》使用DID方法评估了智慧城市建设对城市环境污染的影响,《中国工业经济》期刊官网公布了这篇论文使用的数据和代码。接下来,我就使用这篇论文的数据,给大家分享一下双重差分法(DID)安慰剂检验中随机虚构处理组这种方法的Stata操作。

原文信息

石大千,丁海,卫平,刘建江.智慧城市建设能否降低环境污染[J].中国工业经济,2018(06):117-135.

随机虚构处理组的Stata操作

双重差分法(DID)安慰剂检验的一般做法就是随机选取个体作为处理组,重复500次或者1000次,看看“伪政策虚拟变量”的系数是否显著。在石大千等(2018)这篇论文中,处理组有32个城市,控制组有165个城市,所以我们需要从197个城市中随机选取32个城市作为“伪处理组”,假设这32个城市是智慧城市试点,其他城市为控制组,然后生成“伪政策虚拟变量”(交互项)进行回归。
重复进行500次,我们就会得到500次回归结果(包含“伪政策虚拟变量”估计系数、标准误和p值),最后我们可以绘制出500个“伪政策虚拟变量”估计系数的分布及相应的p值,直观展示安慰剂检验的结果。
说起来很简单,但操作起来还是蛮困难的,我们要解决如下两个关键问题:
(1)如何从197个城市中随机选取32个城市作为“伪处理组”,然后生成“伪政策虚拟变量”?因为是面板数据,所以我们首先只保留一期数据,然后使用sample命令从中随机抽取32个样本,保留所抽取样本的id编号,与原数据进行匹配,最终匹配上的就是我们的“伪处理组”样本,未匹配上的就是控制组样本。
(2)如何存储500次回归中的估计系数、标准误和p值?在Stata中,我们可以生成三个500行1列的矩阵(矩阵中的元素初始值为0),来分别存储500个“伪政策虚拟变量”的估计系数、标准误和p值。每回归一次,就将估计系数、标准误和p值分别存储到三个矩阵的对应位置。

最终代码

*-安慰剂检验-虚构处理组
mat b = J(500,1,0) //* 系数矩阵
mat se = J(500,1,0) //* 标准误矩阵 
mat p = J(500,1,0) //* P值矩阵

*循环500次
forvalues i=1/500{
 use smart_city2018.dta, clear
 xtset id year  //面板数据声明
 keep if year==2005   //保留一期数据
 sample 32, count   //随机抽取32个城市
 keep id  //得到所抽取样本的id编号
 save match_id.dta, replace   //另存id编号数据
 merge 1:m id using smart_city2018.dta //与原数据匹配
 gen treat = (_merge == 3) //将所抽取样本赋值为1,其余为0,得到政策分组虚拟变量
 gen period = (year >= 2012) //生成政策时间虚拟变量
 gen did = treat*period
 reghdfe lnrso did $xlist ,absorb(id year) vce(cluster id) 

 * 将回归结果赋值到对应矩阵的对应位置
 mat b[`i',1] = _b[did]
 mat se[`i'
,1] = _se[did]

 * 计算P值并赋值于矩阵
 mat p[`i',1] = 2*ttail(e(df_r), abs(_b[did]/_se[did]))
}

* 矩阵转化为向量
svmat b, names(coef)
svmat se, names(se)
svmat p, names(pvalue)
* 删除空值并添加标签
drop if pvalue1 == .
label var pvalue1 p值
label var coef1 估计系数
keep coef1 se1 pvalue1 

绘图代码

twoway (kdensity coef1) (scatter pvalue1 coef1, msymbol(smcircle_hollow) mcolor(blue)), ///
title("Placebo Test") ///
xlabel(-0.4(0.1)0.4) ylabel(,angle(0)) ///
xline(-0.171, lwidth(vthin) lp(shortdash)) xtitle("Coefficients") ///
yline(0.1,lwidth(vthin) lp(dash)) ytitle(p value) ///
legend(label(1 "kdensity of estimates") label( 2 "p value")) ///
plotregion(style(none)) ///无边框
graphregion(color(white)) //白底
上图展示了500个“伪政策虚拟变量”估计系数的分布及相应的p值,其中X轴表示“伪政策虚拟变量”估计系数的大小,Y轴表示密度值和p值大小,曲线是估计系数的核密度分布,蓝色圆点是估计系数对应的p值,垂直虚线是DID模型真实估计值-0.171,水平虚线是显著性水平0.1。
从图中可以看出,估计系数大都集中在零点附近,大多数估计值的p值都大于0.1(在10%的水平上不显著),这表明我们的估计结果不太可能是偶然得到的,因而不太可能受到了其他政策或者随机性因素的影响。不过,这个安慰剂检验的结果并不是那么理想,可能还需要在之前的一些环节做些微调,最好能让我们从图中看到DID模型真实估计值是明显的异常值。
这里,我只是提供了一种主流的安慰剂检验的做法,希望大家不要用“定式思维”来看待安慰剂检验。除此之外,还有其他各式各样的安慰剂检验做法,要根据政策情况、样本情况等方面进行具体分析。需要本篇推文使用的数据和代码的朋友,请在公众号后台对话框内回复关键词“placebo”。


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

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