Stata双重差分（DID）实战：从数据清洗到安慰剂检验的完整流程（附代码）
Stata双重差分（DID）实战：从数据清洗到安慰剂检验的完整流程
当我们需要评估某项政策或事件的影响时，双重差分法（Difference-in-Differences, DID）是最常用的因果推断方法之一。这种方法通过比较处理组和对照组在政策实施前后的变化，来识别政策的净效应。本文将带你走过一个完整的DID分析流程，从最基础的数据准备到高级的稳健性检验，每个步骤都配有可立即执行的Stata代码。
1. 数据准备与清洗
任何严谨的实证分析都始于高质量的数据准备。在开始DID分析前，我们需要确保数据格式正确、变量定义清晰，并进行必要的清洗处理。
首先导入数据，假设我们有一个Excel格式的原始数据文件：
import excel "D:/research_data/policy_impact.xlsx", sheet("Sheet1") firstrow case(lower) clear导入后，检查数据结构和变量类型至关重要。使用describe命令可以快速浏览数据概况。常见问题包括数值变量被识别为字符串、日期格式混乱等。对于字符串型数值变量，可以使用destring进行转换：
destring firm_id revenue employees, replace force日期处理是另一个常见痛点。如果数据中有类似"2015-12-31"的日期字符串，我们需要先统一格式：
gen date = date(original_date_var, "YMD") format date %td异常值处理直接影响估计结果的稳健性。Winsorization是常用的方法：
ssc install winsor2 winsor2 revenue profit, replace cuts(1 99)对于面板数据，必须正确定义时间维度和个体维度：
xtset firm_id year2. DID核心变量构建
构建正确的DID模型变量是整个分析的核心。我们需要明确定义三个关键要素：处理组标识、政策时间点和交互项。
假设我们研究的是2018年实施的一项政策，影响部分城市的企业：
// 定义处理组（1=受政策影响的城市，0=不受影响的城市） gen treated = (city_code == "BJ" | city_code == "SH" | city_code == "GZ" | city_code == "SZ") // 定义政策后时期 gen post = (year >= 2018) // 生成DID核心交互项 gen did = treated * post在实际研究中，处理组的定义可能更复杂。例如，政策可能分批次实施，这时需要构建更精细的处理变量：
gen treat_year = . replace treat_year = 2018 if city_code == "BJ" replace treat_year = 2019 if city_code == "SH" replace treat_year = 2020 if city_code == "GZ" gen post = (year >= treat_year) if !missing(treat_year) replace post = 0 if missing(post) gen did = (treat_year <= year) if !missing(treat_year) replace did = 0 if missing(did)3. 基准回归模型设定
DID分析的核心是正确设定回归模型。最基本的DID方程包括处理组虚拟变量、时间虚拟变量和它们的交互项。
reghdfe revenue did treated post, absorb(year) vce(cluster city_code)更严谨的做法是控制个体固定效应和时间固定效应：
reghdfe revenue did, absorb(firm_id year) vce(cluster city_code)在实际应用中，我们通常还需要加入其他控制变量：
global controls size age leverage roa reghdfe revenue did $controls, absorb(firm_id year) vce(cluster city_code)模型结果输出可以使用esttab生成专业表格：
ssc install estout eststo clear eststo: reghdfe revenue did, absorb(firm_id year) vce(cluster city_code) eststo: reghdfe revenue did $controls, absorb(firm_id year) vce(cluster city_code) esttab using "results.rtf", replace b(3) t(3) star(* 0.1 ** 0.05 *** 0.01) nogaps4. 平行趋势检验
DID方法的关键假设是处理组和对照组在政策前具有平行趋势。验证这一假设通常有两种方法：事件研究法和图示法。
事件研究法通过构建政策前后的时间虚拟变量来检验趋势：
forvalues y = 2015/2020 { gen event`y' = (year == `y') * treated } drop event2018 // 以政策前一年为基准期 reghdfe revenue event2015 event2016 event2017 event2019 event2020, absorb(firm_id year) vce(cluster city_code)图示法更直观，可以清晰展示处理组和对照组在政策前后的趋势：
preserve collapse (mean) revenue, by(treated year) twoway (line revenue year if treated==1, lcolor(blue)) /// (line revenue year if treated==0, lcolor(red)), /// xline(2018, lpattern(dash)) legend(label(1 "处理组") label(2 "对照组")) restore5. 稳健性检验与安慰剂检验
为确保DID结果的可靠性，我们需要进行一系列稳健性检验，其中安慰剂检验是最有力的验证方法之一。
5.1 随机分配处理组的安慰剂检验
这种检验通过随机分配"虚假"处理组来验证我们的基准结果是否可能由偶然因素导致：
matrix b = J(500,1,.) matrix se = J(500,1,.) matrix p = J(500,1,.) forvalues i = 1/500 { preserve gen random_treated = runiform() > 0.5 gen random_post = (year >= 2018) gen random_did = random_treated * random_post reghdfe revenue random_did, absorb(firm_id year) vce(cluster city_code) matrix b[`i',1] = _b[random_did] matrix se[`i',1] = _se[random_did] matrix p[`i',1] = 2*ttail(e(df_r), abs(_b[random_did]/_se[random_did])) restore } svmat b, names(coef) svmat se, names(se) svmat p, names(pvalue) twoway (histogram coef1, fraction color(gs12)) /// (scatter pvalue1 coef1, msymbol(Oh) mcolor(red)), /// xline(0.042, lcolor(blue)) // 基准回归中的真实DID系数 legend(label(1 "随机系数分布") label(2 "P值"))5.2 改变政策时间的安慰剂检验
另一种方法是假设政策发生在实际时间之前，检验是否会出现"虚假"效应：
gen placebo_post = (year >= 2016) // 假设政策提前两年实施 gen placebo_did = treated * placebo_post reghdfe revenue placebo_did, absorb(firm_id year) vce(cluster city_code)5.3 排除其他政策干扰
如果样本期间有其他重大政策变化，我们需要验证这些政策是否影响了我们的结果：
gen other_policy = (year >= 2019) // 假设2019年有其他政策 gen other_policy_effect = treated * other_policy reghdfe revenue did other_policy_effect, absorb(firm_id year) vce(cluster city_code)6. 异质性分析与动态效应
了解政策对不同群体的差异化影响，可以丰富研究结论和政策含义。
6.1 分组回归分析
// 按企业规模分组 sum size, detail gen large = (size > r(p50)) if !missing(size) reghdfe revenue did if large==1, absorb(firm_id year) vce(cluster city_code) reghdfe revenue did if large==0, absorb(firm_id year) vce(cluster city_code) // 按行业分组 reghdfe revenue did if industry=="制造业", absorb(firm_id year) vce(cluster city_code) reghdfe revenue did if industry=="服务业", absorb(firm_id year) vce(cluster city_code)6.2 动态处理效应
考察政策效果随时间的变化：
forvalues y = 2018/2020 { gen post`y' = (year == `y') * treated } reghdfe revenue post2018 post2019 post2020, absorb(firm_id year) vce(cluster city_code)7. 结果可视化
专业的结果可视化能极大提升研究的传播效果。对于DID结果，常用的图形包括：
7.1 系数图
coefplot, baselevels keep(post2018 post2019 post2020) /// coeflabels(post2018="2018" post2019="2019" post2020="2020") /// vertical yline(0) xline(0) /// ytitle("政策效应大小") xtitle("年份") /// title("动态处理效应") /// addplot(line @b @at) ciopts(recast(rcap))7.2 平行趋势图
event_plot, stub(event*) default_look /// together graph_opt(xtitle("年份") ytitle("处理效应") /// xlabel(2015(1)2020) title("事件研究图") /// legend(label(1 "点估计") label(2 "95%置信区间")))7.3 安慰剂检验分布图
use placebo.dta, clear sum coef1, detail local true_effect = 0.042 // 替换为实际估计值 twoway (kdensity coef1) /// (function y = 0, range(`=r(min)' `=r(max)') lcolor(black)), /// xline(`true_effect', lcolor(red) lpattern(dash)) /// legend(label(1 "安慰剂系数分布") label(2 "")) /// xtitle("估计系数") ytitle("密度") /// title("安慰剂检验: 真实效应与随机分布比较")8. 实际应用中的注意事项
在完成上述分析后，还有一些关键问题需要考虑：
- 样本选择偏差：处理组和对照组在可观测特征上是否平衡？可以使用
psmatch2进行倾向得分匹配 - 溢出效应：处理组的变化是否影响了对照组？考虑使用空间计量方法检验
- 长期效应：政策效果是暂时的还是持久的？需要足够长的观测期
- 机制分析：政策通过什么渠道产生影响？加入中介变量进行分析
// 倾向得分匹配示例 ssc install psmatch2 psmatch2 treated size age leverage, logit neighbor(3) caliper(0.1) common reghdfe revenue did [pw=_weight], absorb(firm_id year) vce(cluster city_code) // 机制分析示例 reghdfe investment did, absorb(firm_id year) vce(cluster city_code) reghdfe revenue did investment, absorb(firm_id year) vce(cluster city_code)最后，记得保存完整的工作流程：
log using "DID_analysis.log", replace // 所有分析代码 log close