断点回归设计RDD是当前最热门的因果推断计量方法,最主要的原因在于它的透明性和强因果识别性,里面的每一步都可以成功运行出来,若需要do文件和数据dta的请进入计量经济圈社群直接提取(文末)。 gen y = outcome // 结果变量 gen d = running>0 // 处理变量(0/1种类) gen v = running // 分配变量或参考变量 gen vd = v*d // 交互项 local i=1 forvalues i=2/4 { gen v`i'=v^`i' gen v`i'd=v`i'*d } // 产生分配变量的三次方、四次方和他们与处理变量的交互项 qui tab year, gen(dyear) // 如果在面板数据中,想要控制年份可以产生虚拟变量 gen pop2 = pop^2 // 将来用在回归中作为协变量,pop的平方项 ———————————————— ************************************* *图形识别,提供三种方式 ************************************* **1.结果变量是不是在断点处跳跃---------
global sizebin 0.2 //根据你的那个running variable选择箱体,这个你自己设定参数 gen bin=floor(v/$sizebin) gen midbin=bin*$sizebin+0.5*$sizebin bys bin: egen mean=mean(y) reg y d v v2 vd v2d, robust predict fit predict fitsd, stdp gen upfit=fit+1.645*fitsd // 产生置信区间的上边界 gen downfit=fit-1.645*fitsd // 产生置信区间的下边界 preserve // 第一种方式绘制断点回归图 twoway (rarea upfit downfit v, sort fcolor(gs12) lcolor(gs12)) /// (line fit v if v<0, sort="" lcolor(green)="" lwidth(thick))=""> (line fit v if v>0, sort lcolor(red) lwidth(thick)) /// (scatter mean midbin, msize(large) mcolor(black) msymbol(circle_hollow)), /// ytitle('') xtitle('treatment, X (cutoff: X=0)') xline(0, lcolor(black)) /// legend(off) xlabel(-1(0.2)1) title('policy implementation') graph copy all, replace restore
cmogram y v,cut(0) scatter lineat(0) qfitci // 第二种方式绘制断点回归图形 rdplot y v, cut(0) nbins(10) // 第三种方式绘制断点回归图 /**通过图形识别,我们发现在断点处结果变量y发生了跳跃**/ ——————————————— ******************************** *估计结果,使用三种方式 *********************************** **1. 非参数估计-------------- rdrobust y v,c(0) kernel(uni) bwselect(mserd) all // 使用rdrobust进行的非参数估计 rdrobust y v, c(0) kernel(tri) bwselect(mserd) all // 这里使用的是triangular密度估计 rdrobust y v, c(0) kernel(epa) bwselect(mserd) all // 这里使用的是epanechnikov密度估计
**2. 非参数估计---------------------- rd y v, mbw(50 100 200) gr z0(0) kernel(tri) // 这个根据最优带宽计算了三个相应带宽,感觉比较方便 rd y v, mbw(50 100 200) gr z0(0) kernel(rec) // 这里使用的是rectangle密度估计
**3. 参数估计:局部线性回归------ rdbwselect y v, c(0) kernel(uni) bwselect(mserd) // 选择最优带宽 preserve keep if v>= -0.216 & v<= 0.216 =""> eststo x1: reg y d, robust // 面板的话选择xtreg,如果是2sls选择xtivregre eststo x2:reg y d##c.v, robust eststo x3:reg y d##c.(v v2), robust // 局部线性回归法,选择2阶多项式 eststo x4:reg y d##c.(v v2 v3), robust // 局部线性回归法,选择3阶多项式 eststo x5:reg y d##c.(v v2 v3 v4), robust // 局部线性回归法,选择4阶多项式 esttab x1 x2 x3 x4 x5 using y.rtf, star(* .1 ** .05 * .01) nogap nonumber replace /// se(%5.4f) ar2 aic(%10.4f) bic(%10.4f) //输出结果到rtf格式 restore
******************************** *稳健性检验 ********************************* **1. 加入协变量后看看回归结果是不是依然显著----- *1.1 非参估计加入协变量 rd y v, cov(pop pop2) mbw(50 100 200) z0(0) kernel(tri) // 加入协变量pop和pop2
*1.2 参数估计加入协变量 preserve eststo x11: reg y d pop pop2, robust // 加入协变量pop和它的平方项 eststo x21:reg y d##c.v pop pop2, robust eststo x31:reg y d##c.(v v2) pop pop2, robust eststo x41:reg y d##c.(v v2 v3)pop pop2, robust eststo x51:reg y d##c.(v v2 v3 v4) pop pop2, robust esttab x11 x21 x31 x41 x51 using y1.rtf, star(* .1 ** .05 * .01) nogap nonumber replace /// se(%5.4f) ar2 aic(%10.4f) bic(%10.4f) //输出加入协变量后的结果到rtf格式 restore ———————————————— **2.检验其中的协变量是不是在断点处连续-------
**2.1 绘制图形检验一下协变量pop是不是连续的 cmogram pop v,cut(0) scatter lineat(0) qfitci // 第二种方式绘制断点回归图形 rdplot pop v, cut(0) nbins(10) // 第三种方式绘制断点回归图
**2.2 使用估计方法估计出来具体系数看显著不 ** 非参数估计----------------- rdrobust pop v,c(0) kernel(uni) bwselect(mserd) all // 使用rdrobust进行的非参数估计
** 参数估计:局部线性回归------ rdbwselect pop v, c(0) kernel(uni) bwselect(mserd) // 最优带宽的选择 preserve keep if v>= -0.175 & v<= 0.175 =""> eststo xa:reg pop d, robust eststo xb:reg pop d##c.v, robust // 用协变量作为伪结果变量,进行断点回归,选择1阶多项式 eststo xb:reg pop d##c.(v v2), robust // 用协变量作为伪结果变量,进行断点回归,选择2阶多项式 eststo xc:reg pop d##c.(v v2 v3), robust // 用协变量作为伪结果变量,进行断点回归,选择3阶多项式 eststo xd:reg pop d##c.(v v2 v3 v4), robust // 用协变量作为伪结果变量,进行断点回归,选择4阶多项式 restore esttab x11 x21 x31 x41 x51 using m.rtf, star(* .1 ** .05 * .01) nogap nonumber replace /// se(%5.4f) ar2 aic(%10.4f) bic(%10.4f) //输出加入协变量后的结果到rtf格式 /**结果显示pop回归方程不是显著的,所以rdd是适用于此的**/ ———————————————— **3.Mccracy检验:操纵running variable检验--- net install DCdensity, from('http://www./DCdensity') // 安装McCrary检验命令
*注意:以下这个关于分配变量在断点处跳跃的操纵检验会随着下面的binsize和bandwidth设置而不同的
preserve DCdensity v, breakpoint(0) generate(Xj Yj r0 fhat se_fhat) b(0.2) h(0.216) // McCracy test gen upfhat=fhat+1.645*se_fhat gen lowfhat=fhat-1.645*se_fhat twoway (rarea upfhat lowfhat r0 if r0<0, sort="" fcolor(gs12)="" lcolor(gs12))=""> (rarea upfhat lowfhat r0 if r0>0, sort fcolor(gs12) lcolor(gs12)) /// (line fhat r0 if r0<0, lcolor(red))="" (line="" fhat="" r0="" if="" r0="">0, lcolor(blue)) /// (scatter Yj Xj if Yj>0, mcolor(gs4) msymbol(circle_hollow)), /// ytitle('Density') xtitle('') xline(0) legend(off) restore gen t= .079111002/.143889525 // 产生t值,这个需要你根据系数提取出来
display 2*ttail(2651, t) // 得到p值,2651是自由度
/**可以看出在5%显著性水平下实际上Mccrary检验是通不过的,证明没有操纵**/
** 把邻近断点处的那些密度分布放大一些看,这样可以更能清楚地看见是不是有操纵—- preserve DCdensity v, breakpoint(0) generate(Xj Yj r0 fhat se_fhat) b(0.2) h(0.216) // McCracy test local breakpoint 0 local cellmpname Xj local cellvalname Yj local evalname r0 local cellsmname fhat local cellsmsename se_fhat drop if `cellmpname' < -1="" |="" `cellmpname'=""> 0.5 // 把小于-1和大于0.5的部分都去掉 drop if `evalname' < -1="" |="" `evalname'=""> 0.5 tempvar hi quietly gen `hi' = `cellsmname' + 1.96*`cellsmsename' tempvar lo quietly gen `lo' = `cellsmname' - 1.96*`cellsmsename' gr twoway (scatter `cellvalname' `cellmpname', msymbol(circle_hollow) mcolor(gray)) /// (line `cellsmname' `evalname' if `evalname' < `breakpoint',="" lcolor(black)="" lwidth(medthick)) =""> (line `cellsmname' `evalname' if `evalname' > `breakpoint', lcolor(black) lwidth(medthick)) /// (line `hi' `evalname' if `evalname' < `breakpoint',="" lcolor(black)="" lwidth(vthin)) =""> (line `lo' `evalname' if `evalname' < `breakpoint',="" lcolor(black)="" lwidth(vthin)) =""> (line `hi' `evalname' if `evalname' > `breakpoint', lcolor(black) lwidth(vthin)) /// (line `lo' `evalname' if `evalname' > `breakpoint', lcolor(black) lwidth(vthin)), /// xline(`breakpoint', lcolor(black)) legend(off) restore
—————————————— ** 4.安慰剂检验-----------------------
**4.1 改变断点的位置----------------- rdplot y v if v<0, c(-0.25) ="">
rdplot y v if v>0, c(0.25) // 将原来的断点0改变为新的断点0.25
rdrobust y v,c(-0.25) kernel(uni) bwselect(mserd) all // 新断点处使用rdrobust进行的非参数估计 rdrobust y v,c(0.25) kernel(uni) bwselect(mserd) all // 新断点处使用rdrobust进行的非参数估计 /** 通过以上发现改变断点后不显著了,所以我们的断点选择是有道理的**/
**4.2 改变带宽----------------- rdrobust y v,c(0) kernel(uni) h(0.1) all // 改变带宽为0.1 rdrobust y v,c(0) kernel(uni) h(0.4) all // 改变带宽为0.4 /** 通过以上发现改变带宽并没有影响其显著性,因此我们识别的因果效应很稳健**/
|