双重差分法（DID）安慰剂检验的做法：随机抽取500次？
“安慰剂”（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')) ///
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模型真实估计值是明显的异常值。
豫公网安备41160202000603
站长QQ:729038198
关于我们
投诉建议