缘起
- 之前写了一个stata如何做placebo test的文章，本来只是写给自己看，记录自己的学习过程，没想到浏览上千，也是有点意外加激动。
- 有些朋友说根据我写的代码，改到自己的项目里，跑不出来结果，正好今天有学弟问我这个问题，于是我把代码重写了一次，这次用stata官方给的auto数据写的，应该没有bug了，有问题欢迎留言评论。
- 上次的文章链接：
https://www.jianshu.com/p/ef4f920494d
（该文粗略讲了安慰剂的原理，并提供一版自明性相对差点的stata代码，本文在这篇文章的基础上，修改了stata代码。）
实现
假设我的数据集是auto.dta
因变量为price
自变量为rep78
控制变量为headroom weight length
回归模型为OLS，即主回归代码为reg rep78 headroom weight length
下面是1000次placebo test的代码
forvalue i=1/1000{
sysuse auto, clear //调入数据
*- 思路：打乱rep78,即将rep78的全部取值拿出暂存，然后随机赋给每一个样本
*- 打乱rep78,即将rep78的全部取值拿出暂存
g obs_id= _n //初始样本序号
gen random_digit= runiform() //生成随机数
sort random_digit //按新生成的随机数排序
g random_id= _n //产生随机序号
preserve
keep random_id rep78 //保留虚拟的rep78
rename rep78 random_rep78
rename random_id id //重命名为id，以备与其他变量合并（merge）
label var id 原数据与虚拟处理变量的唯一匹配码
save random_rep78, replace
restore
drop random_digit random_id rep78 //删除原来的rep78
rename obs_id id //重命名为id，以备与random_rep78合并（merge）
label var id 原数据与虚拟处理变量的唯一匹配码
save rawdata, replace
*- 合并，回归，提取系数
use rawdata, clear
merge 1:1 id using random_rep78,nogen
reg price random_rep78 headroom weight length
g _b_random_rep78= _b[random_rep78] //提取x的回归系数
g _se_random_rep78= _se[random_rep78] //提取x的标准误
keep _b_random_rep78 _se_random_rep78
duplicates drop _b_random_rep78, force
save placebo`i', replace //把第i次placebo检验的系数和标准误存起来
}
*- 纵向合并1000次的系数和标准误
use placebo1, clear
forvalue i=2/1000{
append using placebo`i' //纵向合并1000次回归的系数及标准误
}
gen tvalue= _b_random_rep78/ _se_random_rep78
kdensity tvalue, xtitle("t值") ytitle("分布") saving(placebo_test)
*-删除临时文件
forvalue i=1/1000{
erase placebo`i'.dta
}
*-mark：后续有时间再将它封装成程序【完】
文献范例
- 周茂、陆 毅、杜 艳、 姚 星，2018：《开发区设立于地区制造业升级》，《中国工业经济》第3期。
image.png
