Stata分位数回归I：理解边际效应和条件边际效应
👇 连享会 · 推文导航 | www.lianxh.cn
🍎 Stata：Stata基础 | Stata绘图 | Stata程序 | Stata新命令 📘 论文：数据处理 | 结果输出 | 论文写作 | 数据分享 💹 计量：回归分析 | 交乘项-调节 | IV-GMM | 时间序列 | 面板数据 | 空间计量 | Probit-Logit | 分位数回归 ⛳ 专题：SFA-DEA | 生存分析 | 爬虫 | 机器学习 | 文本分析 🔃 因果：DID | RDD | 因果推断 | 合成控制法 | PSM-Matching 🔨 工具：工具软件 | Markdown | Python-R-Stata 🎧 课程：公开课-直播 | 计量专题 | 关于连享会
连享会寒假班
作者：曹钰潇 (中山大学)
邮箱：caoyx25@mail2.sysu.edu.cn
编者按：本文主要摘译自下文，特此致谢！
Source：Quantile regressions I -Link-
目录
1. 简介
2. 从线性回归模型开始
3. 三种边际效应解释
3.1 个体效应——对 “我” 来说
3.2 条件效应——对 “像我这样的人” 来说
3.3 无条件效应——对 “所有的人” 来说
4. 总结
5. 相关推文
温馨提示： 文中链接在微信中无法生效。请点击底部「阅读原文」。或直接长按/扫描如下二维码，直达原文：
1. 简介
分位数回归 (Quantile regression，QR) 是回归分析的方法之一，最早由 Roger Koenker 和 Gilbert Bassett 于 1978 年提出。自其出现以来，已经发展出了大量细分模型和估计方法。
今天，我们将集中于三种分位数回归模型：条件分位数回归 (Koenker 和 Bassett，1978)，无条件分位数回归 (Firpo Fortin 和 Lemieux，2009)，以及分位数处理效应 (Firpo，2007)，它们之间的区别和解释常常让人混淆。
那么，它们之间的区别是什么？什么导致了这些区别呢？
正如分位数回归模型主要测量因变量的分布那样，以上三种模型的区别也主要与我们对因变量的测量类型有关。接下来，我们将在线性回归模型的基础上，区分个体边际效应、条件边际效应以及无条件边际效应，以帮助我们更好地理解 QR 模型的含义。
2. 从线性回归模型开始
假设你对分析自变量 怎样影响因变量 感兴趣，同时其他因素 也可能会影响 ，因此你会从这样的模型开始：
现在让我们做一个强的假设，即这些变量之间的关系可以写成一个线性函数。具体来说，假设数据生成过程 (DGP) 为：
同时假设误差项 完全独立于 ：
这个 “基本模型” 比 “入门级” 计量经济学书的单变量模型更先进：首先，模型系数是线性的，但与解释变量没有关系。其次，观察到的因子 和未观察到的因子 相互作用以确定 。由于 不可观测，我们可以认为这个模型实际上是异方差的。
接着从个体层面重写一下模型：
我们利用 Stata 生成数据：
. * 设定种子, 生成 250 个观察值
. clear
. set seed 102
. set obs 250
. * 假设 x 服从卡方分布, v 服从标准正态分布, 同时设定前 5 个观察值以便研究
. * 为作图方便, 利用缺失值将 x 限制在一个范围内
. gen x = rchi2(3)/2
. replace x=. if x>5
. replace x=.75 in 1
. replace x=2 in 2/3
. replace x=4 in 4
. replace x=4.5 in 5
. gen v = rnormal()*.5
. * 模拟数据
. gen y = 5 + 3*x -0.5*x^2+ v * (1+x)
(10 missing values generated)
. * 设置id以便作图
. gen id = _n
. * 生成数据图形
. set scheme s1color
. twoway (scatter y x, color(navy%20)) (scatter y x in 1/5, color(navy%80) mlabel(id) ///
> mlabcolor(black)), ytitle("Outcome y") xtitle("Independent var X") legend(off)
. graph export qr_fig1.png
上图展现了模拟样本的分布。按照 DGP， 和 之间的关系是非线性的，从图中可以看出其异方差性。同时为了方便起见，我们标记 5 个点作为样本中随机观察的例子。
那么，我们该如何解释这个模型呢？
3. 三种边际效应解释
3.1 个体效应——对 “我” 来说
解释的第一种方法是对数据中的每一个人或观测点进行解释。然而，这是一种很少会看到的解释方式，因为它要求我们知道 DGP 的每一个信息，这不仅包括特征 ，还包括不可观察因素 (或者模型是同方差的)。
通过将 看做一个连续变量，我们可以使用方程 (3) 来获得关于 的偏导，这样可以得到对个人 来说 对 的边际影响：
这里有几点需要注意：
只有当我们能观察到不可观察因素 ，或者模型是同方差，即 时，才能得到个人边际效应； 该效应不是恒定的，它取决于 和 ； 这是一个局部 (线性) 的 对 边际影响的近似值。如果 变化很大，这个近似值可能与模型预测的 的变化有很大差别。然而，在大多数情况下，这个近似值被认为是足够好的。
我们也可以根据 DGP 得到边际效应的准确值，如果 增加 单位，个人 的结果是：
接下来，我们假设 ，在 Stata 中绘制出展示边际效应的图。
* 数据清理
. gen dx=1
. gen xdx=x+dx
. gen xdxr=x+dx+0.1
. gen y_dy = 5 + 3*(xdx) -0.5*(xdx)^2+ v * (1+(xdx))
. * 画图
. twoway (scatter y x in 1/5, mlabel(id) mlabcolor(black) mlabpos(12) color(navy%80)) ///
> (scatter y_dy xdx in 1/5, mlabcolor(black) mlabpos(12) ms(d) color(maroon%80)) ///
> (pcarrow y x y_dy xdx in 1/5, color(black%40)) ///
> (rcap y y_dy xdxr in 1/5, color(black%40)) ///
> (scatter y x, color(navy%10)) (scatter y_dy xdx, color(maroon%10) ms(d)), ///
> ytitle("Outcome y") xtitle("Indep var X") xlabels(0/5) ///
> legend(order(1 "f(x,v)" 2 "f(x+dx,v)" 4 "{&Delta}y") col(3))
. graph export qr_fig2.png
如图所示，不同观测点的 的变化取决于它们的特征值 ，而且它们似乎紧跟总体趋势。然而，注意观察点 2 和 3，我们会发现两者有相同的 和不同的 ，这就导致了它们 的变化不同，即一种未观察到的异质性。换句话说，一个单位 的变化对每个人的影响是特定的。
我们也可以把边际效应的近似值与准确值进行比较：
. * 近似值
. gen mfx_aprox=3-x+v
. label var mfx_aprox "LL approximation"
. * 准确值
. gen mfx_exact=(y_dy-y)
. label var mfx_exact "Exact effect"
. list x v mfx_aprox mfx_exact in 1/5, sep(0)
+-----------------------------------------+
| x v mfx_aprox mfx_exact |
|-----------------------------------------|
1. | .75 .5645621 2.814562 2.314562 |
2. | 2 .1400093 1.140009 .6400089 |
3. | 2 -.1684819 .8315181 .3315182 |
4. | 4 .0549622 -.9450378 -1.445037 |
5. | 4.5 -.5240271 -2.024027 -2.524027 |
+-----------------------------------------+
对于观察点 1 ，增加 1 单位 将使其增加 2.315 个单位。对于观察点 2 和 3，同样的变化将使点 2 的结果增加 0.64 个单位，而点 3 的结果只增加 0.33 个单位，这是因为他们有不同的不可观测因素 。
3.2 条件效应——对 “像我这样的人” 来说
第一种解释表明，同样的 “政策” 适用于具有相同观察特征 的人，可能会因为未观察到的异质性 而产生不同结果。然而，这种方法是没有用的，因为 从来没有被观察过。因此，我们无法解释对任何特定个体的影响。
相比之下，第二种解释是更实际的，它不是试图量化 ，而是将不可观察因素通过平均化的方法从方程中剔除！
首先回顾一下 DGP：
如果我们把期望值的条件 看做是一个特定的值 ，可以得到：
使用零条件均值假设 ，可以得到我们更为熟悉的表达式：
此时，如果我们采取与第一种解释相同的方法，可以得到方程 (6) 相对于 (不是 ) 的偏导数，以求出一单位 的增加对 的 “平均/预期” 影响：
接下来，我们运用 Stata 得到以 为条件的 的预期均值，以及如果 变化 1 单位，条件均值将如何变化。
. * E(y|X)
. gen yy = 5 + 3*x -0.5*x^2
. * E(y'|(X+dx))
. gen yy_dy = 5 + 3*(x+dx) -0.5*(x+dx)^2
. * 画图
. twoway (scatter y x in 1/5, color(navy%25)) ///
> (scatter y_dy xdx in 1/5, ms(d) color(maroon%25)) ///
> (pcarrow y x y_dy xdx in 1/5, color(black%40)) ///
> (function 5 + 3*x -0.5*x^2, range(0 6) color(black%60)) ///
> (scatter y x, color(navy%10)) (scatter y_dy xdx, color(maroon%10) ms(d)) ///
> (scatter yy x if inlist(id,1,3,4,5), ///
> mlabel(id) mlabcol(black) mlabpos(12) color(navy%80)) ///
> (scatter yy_dy xdx if inlist(id,1,3,4,5), ///
> mlabel(id) mlabcol(black) mlabpos(6) color(maroon%80) ms(d) ), ///
> ytitle("Outcome y") xtitle("Indep var X") xlabels(0/5) ///
> legend(order(1 "f(x,v)" 2 "f(x+dx,v)" 4 "E(y|X)") col(3))
. graph export qr_fig3.png
如图所示，变化前后所有的预期均值都落在同一条 “线” 上——条件均值函数。不同于第一种解释适用于 “我”，解释回归系数的第二种方式适用于 “像我这样的人”。这意味着，这些影响并不适用于任何特定的观察值，而是发生在与 “我” 有相同特征的观察值中的平均情况，这就是条件效应的性质。
例如，观察值点 2 和点 3，它们都设定了特征 ，因此会出现同样的 0.5 个单位的 变化。
换句话说，第二种解释指的是在所有具有相同特征 的人中，预期效果或者条件均值上观察到的变化。如上，我们也可以解释为：如果一个 的人经历了一单位 的增加，它们的结果 预计会增加 0.5 个单位。同样也可以这样解释：比较两组具有完全相同特征的人，除了其中一组 ，另一组 ，第二组的平均结果 将高出 0.5 个单位。
因此再次强调：第二种解释需要从群体 (由其特征定义) 的角度思考，而不是从个人的角度思考。
3.3 无条件效应——对 “所有的人” 来说
第二种解释是我们最熟悉的一种，因为这也是大多数教科书所采用的解释方法。它也可能是我们最感兴趣的一种，因为它是我们在 “微观” 层面上最接近测量影响的方法。而第三种解释可以被认为是更多的 “宏观” 层面的解释。
从政策制定者的角度来看，第三种解释可能更实用，因为可以了解到由于一项改变了所有人 的政策，对人口的平均结果 会产生什么影响。比如说：如果人口的受教育年限平均增加 1 年，玻利维亚的贫困率将如何变化？
首先仍然考虑一下 DGP：
或者更好的是条件期望形式：
因为我们对均值的无条件影响感兴趣，所以应该用无条件均值的方式来写这个方程：
为了进一步简化，我们还要使用：
最后，无条件均值结果方程就可以写成：
从方程中可以发现：
人口的无条件均值 取决于 的无条件均值，； 人口的无条件均值 也取决于 的无条件方差，； 只有当我们知道 的一些分布特性时，我们才能预测 的无条件均值。
为简化方程，我们需要一个影响 ，而不是方差的所有 的改变。例如，每一个观察值都经历了一单位的增加。这就是 Firpo 等 (2000) 所说的 “位置转移” 效应：分布改变了“位置”，但没有改变形状。
因此，方程可以简化为只与 有关的无条件边际效应：
通过此方程我们可以直接计算出平均变化，不过需要注意：
解释一表明，由于观察到的 和未观察到的 的异质性，个体边际效应对每个人都是独一无二的； 解释二表明，条件边际效应 (均值) 的变化只是因为观察到 的异质性； 解释三表明，由于对应整个人口，无条件边际效应 (均值) 是恒定的。
接下来我们利用 Stata 来获得结果。
. * X 的均值
. /*
> 因为 2*X~Chi(3), E(X)=0.5*k =1.5, Var(X)=0.5*k=1.5
> 所以 E(y)=5 + 3*E(x) -0.5*E(x)^2-0.5*Var(X)
> E(Y) = 7.625 E(Y') = 8.625
> */
. twoway (scatter y x, color(navy%15)) (scatter y_dy xdx, ms(d) color(maroon%15)) ///
> (function 5 + 3*x -0.5*x^2, range(0 6) color(black%60)) ///
> (scatteri 7.625 1.5 "Before",mcolor(navy) msize(small) mlabcolor(black)) ///
> (scatteri 8.625 2.5 "After", mcolor(maroon) msize(small) mlabcolor(black)) ///
> (pcarrowi 7.625 1.5 8.625 2.5, color(black)), ///
> xline(1.5, lcolor(navy%50)) xline(2.5, lcolor(maroon%50)) ///
> yline(7.625, lcolor(navy%50)) yline(8.625, lcolor(maroon%50)) ///
> ylabel( 0(5)15 7.625 8.625, angle(0)) ///
> legend(order(1 "f(x,v)" 2 "f(x+dx,v)" 4 "E(y)" 5 "E(y')") col(4)) ///
> ytitle("Outcome y") xtitle("Indep var X")
. graph export qr_fig4.png
此图中，蓝点表示 变化前的数据，红点表示变化后的数据。请注意，点 没有落在 “条件均值” 线上，但是其变化是与它平行的。那么，该如何解释呢？
如果每个观察值 都增加 1 个单位，那么 就会增加 1 个单位，但方差保持不变，从而导致 的无条件均值 将增加 1 个单位，从 7.625 增加到 8.625。
4. 总结
在线性回归模型的框架内，本文介绍了一种不同的方式来思考个人、条件和无条件的边际效应，这些解释之间的差异取决于我们在进行解释时的对象是谁。
在数学上，边际效应可以被描述如下：
5. 相关推文
Note：产生如下推文列表的 Stata 命令为：
lianxh 边际, m
安装最新版lianxh命令：
ssc install lianxh, replace
专题：Stata教程 Stata-Python交互-5：边际效应三维立体图示 专题：数据处理 Stata：边际处理效应及其可视化-mtefe-T309 专题：回归分析 Stata：边际效应知多少？f_able命令（上） Stata：边际效应知多少？f_able命令（下） Stata: 手动计算和图示边际效应 Stata: 边际效应分析 专题：倍分法DID DID边际分析：让政策评价结果更加丰满 专题：交乘项-调节 Logit-Probit中的交乘项及边际效应图示 Stata：图示连续变量的连续边际效应 专题：Probit-Logit Logit-Probit：非线性模型中交互项的边际效应解读 专题：生存分析 mrsprep：边际相对生存分析的数据准备
New！ Stata 搜索神器：
lianxh和songblGIF 动图介绍
搜： 推文、数据分享、期刊论文、重现代码 ……
👉 安装：
. ssc install lianxh
. ssc install songbl
👉 使用：
. lianxh DID 倍分法
. songbl all
🍏 关于我们
连享会 ( www.lianxh.cn，推文列表) 由中山大学连玉君老师团队创办，定期分享实证分析经验。 直通车： 👉【百度一下：连享会】即可直达连享会主页。亦可进一步添加 「知乎」,「b 站」,「面板数据」,「公开课」 等关键词细化搜索。