Stata模拟:控制变量!控制变量!Good-Controls-Bad-Controls

发布时间:2022-09-08 阅读 328

Stata连享会   主页 || 视频 || 推文 || 知乎 || Bilibili 站

温馨提示: 定期 清理浏览器缓存,可以获得最佳浏览体验。

New! lianxh 命令发布了:
随时搜索推文、Stata 资源。安装:
. ssc install lianxh
详情参见帮助文件 (有惊喜):
. help lianxh
连享会新命令:cnssc, ihelp, rdbalance, gitee, installpkg

课程详情 https://gitee.com/lianxh/Course

课程主页 https://gitee.com/lianxh/Course

⛳ Stata 系列推文:

PDF下载 - 推文合集

作者:李烨阳 (浙江大学)
邮箱li_yeyang95@163.com

编者按:本文主要整理自下文,特此致谢!
Source:Cinelli, C., A. Forney, J. Pearl, 2022, A crash course in good and bad controls, Sociological Methods & Research, Forthcoming. -PDF- -Link- -R- -Stata-


目录


1. 导言

在连享会推文「控制变量!控制变量!Good-Controls-Bad-Controls」中,我们已经介绍了如何使用因果图 (有向无环图,DAG) 筛选控制变量的理论知识,以及几种社会科学研究中较为典型的良好控制、糟糕控制和中性控制情形。在本文中,我们将使用 Stata 软件,对前文提到的 18 中模型进行模拟分析,以检验理论是否成立。

2. 良好的控制变量

为了方便读者阅读,本文模拟分析的顺序与「控制变量!控制变量!Good-Controls-Bad-Controls」一致,模型的编号与原文一致。

Stata 社区发布的 dag 程序包提供了对有向无环图进行模拟分析的外部命令。但由于这一程序包主要是面向流行病学研究设计的,它生成的模拟变量为 0-1 变量。而在社会科学研究中,我们常常需要使用到连续变量,同时为与原文的 R 语言模拟分析保持一致,本文并没有直接使用该程序包,感兴趣的读者可以自行尝试。

2.1 解释变量和被解释变量的混淆情形

模型 1-3 是存在同时影响解释变量和被解释变量的混淆的情形,我们已经知道,图中的 Z 都是良好的控制变量。

首先我们对模型 1 进行模拟分析。

第一步,先根据 DAG 图中变量间的关系,定义模型 1 的模拟程序:每次模拟我们随机生成 10000 个样本,z 是一个符合 N(0,1) 标准正态分布的随机变量,X 等于 Z 加上一个符合 N(0,1) 标准正态分布的随机变量,Y 等于 X 加 Z 再加上一个符合 N(0,1) 标准正态分布的随机变量。然后进行两次回归估计,第一次不控制变量 Z,第二次控制变量 Z,并分别记录每次回归 X 系数估计结果的点估计和 95% 置信区间的上界和下界。

*- step 1 定义模型的数值模拟程序

	* model 1 模拟程序
	capture program drop model_1 // 模型命名为model_1, 若有更换注意修改
	program model_1,rclass       // 若有更换注意修改
	clear
	set obs 10000

	* 此部分为针对不同模型结构需要修改的部分
	***************************
	* model 1: DAG结构
	gen z = rnormal()
	gen x = z + rnormal()
	gen y = x + z + rnormal()
	***************************

	* 不控制 z
	reg y x
	mat A=r(table)
	return scalar mean0=A[1,1]
	return scalar low0=A[5,1]
	return scalar up0=A[6,1]
    
	* 控制 z
	reg y x z 
	mat A=r(table)
	return scalar mean1=A[1,1]
	return scalar low1=A[5,1]
	return scalar up1=A[6,1]
	end

第二步,进行重复模拟分析。我们进行 100 次模拟,为了保证结果可复现,设定了种子值。

*- step 2 模拟模块
	clear
	set seed 11
	simulate mean0=r(mean0) low0=r(low0) up0=r(up0) ///
	mean1=r(mean1) low1=r(low1) up1=r(up1), reps(100): model_1 // 对 model_1 进行 100 次模拟, 若有更换注意修改
	gen id=_n

第三步,我们对模拟结果进行可视化。左图为不控制变量 Z 的情形,右图为控制变量 Z 的情形。横坐标是 100 次模拟的序列,纵坐标是 X 系数的估计值。蓝色附加线为 100 次模拟的 X 系数点估计值的均值,红色附加线为真实的数据模拟设定的 X 系数,本文中所有模型均设置为 1。

*- step 3 模拟结果可视化

	set scheme s1color //设定绘图模板
	sum mean0
	local b_mean `r(mean)'
	twoway (scatter mean0 id, ms(Oh) mc(gs5)) (rcap low0 up0 id, lc(gs5)), ///
	ytitle("E(β{sub:x}_hat)") legend(order(1 "E(β{sub:x}_hat)点估计" 2 "E(β{sub:x}_hat)置信区间")) ///
	yline(`b_mean', lcolor(blue)) 
	graph save model_1a.gph, replace //图名model_1a,若有更换注意修改

	sum mean1
	local b_mean `r(mean)'
	twoway (scatter mean1 id, ms(Oh) mc(gs5)) (rcap low1 up1 id, lc(gs5)), ///
	ytitle("E(β{sub:x}_hat|z)") legend(order(1 "E(β{sub:x}_hat|z)点估计" 2 "E(β{sub:x}_hat|z)置信区间")) ///
	yline(`b_mean', lcolor(blue))  yline(1, lcolor(red)) 
	graph save model_1b.gph, replace  //图名model_1b,若有更换注意修改
	
	graph combine model_1a.gph  model_1b.gph , /// //若有更换注意修改图名
	xsize(10) ysize(5) title("Model 1") //图形合并,图标题名Model 1,若有更换注意修改

结果如下图所示。可以看出,在不控制 Z 的情况下,100 次模拟的估计系数均值约为 1.5,且从图中可以看出,95% 置信区间均没有包括真实系数,估计系数是有偏的。而引入了控制变量 Z 之后,100 次模拟的估计系数均值与真实值非常接近,且绝大部分的单次模拟的 95% 置信区间均包含真实值。因此数值模拟结果支持了模型 1 中的 Z 是良好的控制变量这一观点。

为了节约篇幅,在后续的模型数值模拟介绍中,只给出步骤 1 定义模型的数值模拟程序步骤中 DAG 结构部分的代码,读者只需要将这一部分代码替换到上述例子的 step 1 中同样用一组星号分割开的位置中即可,并更改相应的模型命名、图名和图的标题名,其余部分毋须改动。当然,大家也可以通过点击「good_bad_controls_simulation.do」,下载完整的模拟代码。

模型 2 DAG 结构部分的代码为:

	***************************
	* model 2: DAG 结构
	gen u = rnormal()
	gen z = u + rnormal()
	gen x = z + rnormal()
	gen y = x + u + rnormal()
	***************************

模型 2 模拟结果为如下图所示,同样印证了 Z 是一个好的控制变量。

模型 3 DAG 结构部分的代码为:

	***************************
	* model 3: DAG结构
	gen u = rnormal()
	gen z = u + rnormal()
	gen x = u + rnormal()
	gen y = x + z + rnormal()
	***************************

其结果同样与理论一致。

2.2 解释变量和中介变量的混淆情形

模型 4、5、6 是存在同时影响解释变量和中介变量的混淆的情形,此三种情况同样是良好的控制变量。

下面分别给出模型 4、5、6 DAG 结构部分的代码,模拟分析结果均显示,不控制 Z 时存在混淆偏差,而控制了 Z 之后都可以得到 X 系数无偏的估计。

	***************************
	* model 4: DAG 结构 
	gen z = rnormal()
	gen x = z + rnormal()
	gen m = x + z + rnormal()
	gen y = m + rnormal()
	***************************

	***************************
	* model 5: DAG 结构 
	gen u = rnormal() 
	gen z = u + rnormal()
	gen x = z + rnormal()
	gen m = x + u + rnormal()
	gen y = m + rnormal()
	***************************

	***************************
	* model 6: DAG 结构 
	gen u = rnormal() 
	gen z = u + rnormal()
	gen x = u + rnormal()
	gen m = x + z + rnormal()
	gen y = m + rnormal()
	***************************

模型 4-6 的模拟结果为:

2.3 处理后变量一定是坏的控制变量吗

虽然大部分情况下处理后变量是坏的控制变量,但也有例外。例如模型 15 中,如果以 W 为条件是给定的前提,此时 Z 是一个好的控制变量。

定义模型 15 的数值模拟程序为:

*- step 1 定义模型的数值模拟程序

	* model 15 模拟程序
	capture program drop model_15 
	program model_15,rclass 
	clear
	set obs 10000
	
	***************************
	* model 15: DAG 结构 
	gen x = rnormal() 
	gen u = rnormal() 
	gen z = x + rnormal()
	gen w = z + u + rnormal()
	gen y = x -2*u + rnormal()
	***************************

	* 控制 w, 不控制 z
	reg y x w
	mat A=r(table)
	return scalar mean0=A[1,1]
	return scalar low0=A[5,1]
	return scalar up0=A[6,1]

	* 控制 w 和 z
	reg y x w z 
	mat A=r(table)
	return scalar mean1=A[1,1]
	return scalar low1=A[5,1]
	return scalar up1=A[6,1]
	end

模拟结果如下图所示,左图为只控制 W,此时会造成选择性偏误,而当控制 W 的同时也控制 Z,由于控制对撞变量而开放的路径又重新截断了,可以获得无偏的估计。

3. 糟糕的控制变量

3.1 处理前变量一定是好的控制变量吗

模型 7 是一个 M 偏误的典型例子,这里的 Z 虽然是处理前变量,但却是糟糕的变量。

模型 7 (a) DAG 结构部分的代码为:

	***************************
	* model 7: DAG结构 
	gen u1 = rnormal() 
	gen u2 = rnormal() 
	gen z = u1 + u2 + rnormal()
	gen x = u1 + rnormal()
	gen y = x -4*u2 + rnormal()
	***************************

模拟结果如下图所示,由于 Z 是一个对撞变量,在不控制 Z 时可以获得无偏估计,控制 Z 反而会引入非因果的相关性。

模型 10 同样是一个控制变量虽然是处理前变量,但仍然是糟糕的控制变量的例子。

模型 10 DAG 结构部分的代码为:

	***************************
	* model 10: DAG结构 
	gen z = rnormal() 
	gen u = rnormal() 
	gen x = -2*z + u + rnormal()
	gen y = x + 2*u + rnormal()
	***************************

模拟结果如下图所示,在不增加控制变量 Z 时,由于有混淆会造成估计偏差,100 次模拟的估计系数均值约为 1.3,但是这里引入控制变量 Z 后,不仅不能截断后门路径,还会放大偏差,100 次模拟的估计系数均值约为 2.0。

3.2 警惕阻断因果路径

模型 11、12 是由于控制中介变量或其子变量造成过度控制偏差的典型范例。

模型 11、12 DAG 结构部分的代码分别为:

	***************************
	* model 11: DAG 结构 
	gen x = rnormal() 
	gen z = x + rnormal() 
	gen y = z + rnormal()
	***************************

    ***************************
	* model 12: DAG 结构 
	gen x = rnormal() 
	gen m = x + rnormal() 
	gen z = m + rnormal() 
	gen y = m + rnormal()
	***************************

模拟结果分别如下图所示,在不控制中介变量或其子变量时,估计是无偏的,但引入他们作为控制变量则会导致间接因果效应路径被截断,从而低估结果,因此中介变量或其子变量是糟糕的控制变量。

3.3 警惕打开对撞变量所在路径

模型 16 和模型 17 中的 Z 均为对撞变量,如果控制对撞变量将引入选择性偏差,因此是糟糕的控制变量。

模型 16、17 DAG 结构部分的代码分别为:

	***************************
	* model 16: DAG 结构 
	gen x = rnormal() 
	gen u = rnormal() 
	gen z = x + u + rnormal() 
	gen y = x + 2*u + rnormal()
	***************************

	***************************
	* model 17: DAG 结构 
	gen x = rnormal() 
	gen y = x + rnormal()
	gen z = x + y + rnormal()
	***************************

模拟结果分别为:

控制模型 18 中的 Z 会导致案例控制偏差。

模型 18 DAG 结构部分的代码为:

	***************************
	* model 18: DAG 结构 
	gen x = rnormal() 
	gen y = x + rnormal()
	gen z = y + rnormal()
	***************************

模拟结果为:

4. 中性的控制变量

4.1 可能提高精度的情形

模型 8 中的 Z 是一个中性的控制变量,但控制它可以提高估计结果的精确度。

模型 8 DAG 结构部分的代码为:

	***************************
	* model 8: DAG 结构 
	gen z = rnormal() 
	gen x = rnormal()
	gen y = x + 2*z + rnormal()
	***************************

从模拟结果可以看出,虽然无论是否控制 Z,都可以获得无偏估计,但是控制 Z 可以提高估计精度,即系数估计值的标准误降低,系数置信区间范围缩小。

模型 13 中的 Z 同样是一个提高估计精度的中性的控制变量。

模型 13 DAG 结构部分的代码为:

	***************************
	* model 13: DAG 结构 
	gen z = rnormal() 
	gen x = rnormal()
	gen m = 2*z + rnormal()
	gen y = x + 2*m + rnormal()
	***************************

模拟结果为:

4.2 可能降低精度的情形

模型 9 和 14 中的 Z 均为中性的控制变量。但控制它们会减少 X 的变动性,从而降低估计结果的精确度。

模型 9 和 14 DAG 结构部分的代码分别为:

	***************************
	* model 9: DAG 结构 
	gen z = rnormal() 
	gen x = 2*z + rnormal()
	gen y = x + 2*rnormal()
	***************************

	***************************
	* model 14: DAG 结构 
	gen x = rnormal() 
	gen z = 2*x + rnormal()
	gen y = x + 2*rnormal()
	***************************   

从模拟结果可以看出,虽然无论是否控制 Z,都可以获得无偏估计,但是控制 Z 会估计精度,即系数估计值的标准误增加,系数置信区间范围增大。

5. 相关推文

Note:产生如下推文列表的 Stata 命令为:
lianxh 控制变量, m
安装最新版 lianxh 命令:
ssc install lianxh, replace

相关课程

免费公开课

最新课程-直播课

专题 嘉宾 直播/回看视频
最新专题 文本分析、机器学习、效率专题、生存分析等
研究设计 连玉君 我的特斯拉-实证研究设计-幻灯片-
面板模型 连玉君 动态面板模型-幻灯片-
面板模型 连玉君 直击面板数据模型 [免费公开课,2小时]
  • Note: 部分课程的资料,PPT 等可以前往 连享会-直播课 主页查看,下载。

课程主页

课程主页

关于我们

  • Stata连享会 由中山大学连玉君老师团队创办,定期分享实证分析经验。
  • 连享会-主页知乎专栏,700+ 推文,实证分析不再抓狂。直播间 有很多视频课程,可以随时观看。
  • 公众号关键词搜索/回复 功能已经上线。大家可以在公众号左下角点击键盘图标,输入简要关键词,以便快速呈现历史推文,获取工具软件和数据下载。常见关键词:课程, 直播, 视频, 客服, 模型设定, 研究设计, stata, plus, 绘图, 编程, 面板, 论文重现, 可视化, RDD, DID, PSM, 合成控制法

连享会小程序:扫一扫,看推文,看视频……

扫码加入连享会微信群,提问交流更方便

✏ 连享会-常见问题解答:
https://gitee.com/lianxh/Course/wikis

New! lianxhsongbl 命令发布了:
随时搜索连享会推文、Stata 资源,安装命令如下:
. ssc install lianxh
使用详情参见帮助文件 (有惊喜):
. help lianxh