ARTICLE DETAIL

资讯详情

深耕郑州网站建设与运营推广的一线实战洞察。

Stata分位数回归实操:从qreg到sqreg,揭示分布异质性

Stata分位数回归实操:从qreg到sqreg,揭示分布异质性 如果你的研究还停留在OLS回归的阶段哪怕是用Stata做了一堆稳健标准误、固定效应、工具变量你看到的也只是“平均效应”这一个切面。最近好几个做收入分配、教育回报、健康不平等的朋友来问我分位数回归怎么做这篇我直接把我自己在Stata里的完整实操流程、命令逻辑和踩过的坑写出来一次性说清楚。分位数回归不是又一个需要“炫技”的冷门模型它解决的是非常实际的问题当自变量的效应在因变量分布的不同位置表现不同时只有分位数回归能把这个异质性刻画出来。比如教育年限对工资的影响在月薪3000的人身上和在月薪3万的人身上很可能完全不是一回事——前者每多一年教育可能只涨200块后者可能是2000块。这种现象在医学、经济学、社会学里到处都是。而Stata恰好把分位数回归的主流估计、推断和检验都集成好了剩下的问题就是你知不知道怎么用、用得对不对。这篇我会用Stata内置的auto.dta数据集做演示从“为什么需要分位数回归”这个根本问题讲起再逐步拆解qreg、sqreg、bqreg、grqreg这些命令的选用逻辑、输出解读、标准误设置最后再补上面板分位数回归的实操方向和跨分位点系数检验。全程给出可直接复制运行的命令也会交代每一段代码背后的判断依据。1. 为什么需要分位数回归OLS只看中心分位数回归看全貌1.1 OLS的盲区均值掩盖了分布里的异质性OLS估计的是条件期望函数E(Y|X)也就是给定自变量X时因变量Y的平均水平。问题在于这个均值本身可能被极端值拉扯得面目全非更关键的是它对“分布的不同位置”完全无感。打个比方你想研究培训项目对收入的影响OLS告诉你平均增加了800块。听起来项目很成功。但真实情况可能是原本收入最低的20%的人参加了培训后收入几乎没有变化而原本收入最高的10%的人收入增加了3000块。800块这个平均数既没有反映底层人群的真实困境也没有刻画顶层人群的真实受益。更麻烦的是收入这类变量往往是右偏的均值会系统性地偏向高分位OLS估计出的“平均效应”对政策制定者来说参考价值很有限。分位数回归最早由Koenker和Bassett在1978年提出它估计的是给定X时Y的各个条件分位数函数Q_Y(τ|X)。τ取0.1、0.25、0.5、0.75、0.9等数值含义就是“在因变量分布的对应位置上X每变动一个单位因变量变动多少”。这个逻辑用一个通俗说法概括就是OLS回答“平均情况怎么样”分位数回归回答“最差的一批人、中等的一批人、最好的一批人各自的情况怎么样”。1.2 损失函数的本质差异为什么分位数回归能估计出不同切面OLS通过最小化残差平方和来求解参数min Σ (yi - xiβ)²这个目标函数等价于估计条件均值。因为平方函数是对称的、向上凸的它对大残差极其敏感——这是OLS会被离群值强烈影响的原因。分位数回归则不同。它最小化的目标函数是一个名为“检查函数”的非对称绝对损失min Σ ρτ(yi - xiβ)其中ρτ(u) u(τ - 1{u 0})。这个函数看着复杂实际含义很朴素如果你的预测值高于实际值也就是残差为正那么权重是τ反过来残差为负权重是1-τ。比如τ0.9时对“低估”的惩罚远高于“高估”因此估计结果会更偏好让绝大部分样本落在拟合线之上从而刻画出数据的第90个百分位附近。很多人不太注意这一点这个目标函数是不可导的OLS那种偏导等于0的求法用不上Stata内部用的是线性规划算法。你不需要自己编程求解但理解这个机制有助于避开后面要讲的推断和标准误坑。1.3 什么时候你该换用分位数回归而不是继续OLS在实践中我一般在以下四类场景中会主动改用分位数回归存在显著的异质性效应理论上有理由认为X对Y不同水平的影响不同。因变量存在严重偏态或离群值收入、医疗费用、事故损失、股价波动等都是典型。研究重点是弱势群体或特定群体比如想分析最低收入人群是否更能从最低工资政策中获益。需要和OLS互相印证如果分位数回归显示各个分位点系数高度一致那么OLS的结论更可信如果差异很大那说明单纯平均效应可能会误导决策。2. Stata实现分位数回归前的准备数据形态与常用命令矩阵2.1 数据要求的几个隐性前提分位数回归对数据的要求并不比OLS苛刻但有几个隐性前提需要特别注意。第一因变量必须是连续变量或至少是测度良好、近似连续的变量。0-1变量、计数变量虽然也能跑分位数回归但对结果的解释要非常小心后文我会专门讲这个问题。第二样本量不能太小。分位数回归相当于在因变量某一个局部范围内做估计数据量小了之后每个分位点附近的有效样本会迅速减少系数估计的波动会非常大。我的经验是总样本低于300时除非分位点选得比较保守比如0.25到0.75之间否则结果很难稳定。第三数据中没有离谱的编码错误。比如把缺失值记录为9999这在OLS中可能只是一个离群值干扰均值在分位数回归中甚至会让某一整个尾部的结果严重偏移因为你实际上等于给这个样本赋予了极端大的权重。2.2 主命令qreg、sqreg、bqreg的选用逻辑Stata里做分位数回归的命令有多个容易搞混。我的使用原则把它简化成一张表命令能干什么标准误建议使用场景qreg单个分位点的分位数回归默认基于核密度估计的渐近标准误旧版传统算法新版也支持bootstrap等快速看结果、做稳健性检查bqreg单个分位点但使用bootstrap标准误bootstrap分位点较极端、样本量一般时sqreg同时估计多个分位点自动计算bootstrap标准误bootstrap最推荐的日常主力命令iqreg分位区间回归比如同时估计0.25到0.75区间bootstrap较少用主要用于区间推断sqreg是大部分人最应该优先使用的命令。它的优势不只是在一次运行中输出多个分位点系数更关键的是它基于bootstrap标准误可以放心直接报告P值。而老式的qreg默认标准误在大样本下理论性质没问题在中等样本下往往偏保守或偏自由实际检验功效不可控。2.3 用小命令把数据摸清楚再跑模型在跑任何回归之前我会先用Stata的基础描述统计命令把数据的分布特征摸一遍。这一步虽然不起眼但能避免大量后续麻烦。sysuse auto.dta, clear summarize price mpg, detaildetail选项会输出最小值、最大值、1%分位点、5%分位点、25%分位点、中位数、75%分位点、95%分位点等。我习惯先看两个变量的分位数差异。比如price的中位数是5006.5均值是6165.3均值明显高出中位数一大截这就是典型的右偏分布。此时如果只用OLS研究mpg对price的影响高价位车会把平均效应整体拉高而分位数回归能分别告诉你中低价位车和高价位车对油耗的敏感程度到底差多少。此外egen命令族也值得先跑一跑特别是看每个分位组内部的样本量egen price_q4 cut(price), group(4) table price_q4, statistic(count price) statistic(mean price)cut命令按四分位将price分组然后看每组的样本量和均值。这一步能帮你确认每个分位点的样本支撑是否足够避免后面对极端分位点的结论产生误判。这就是我在实操中会自然用到egen、table这些命令的场景——它们不是分位数回归专属但却是数据预处理中绕不开的环节。3. 完整的Stata实操流程以auto数据为例从建模到输出3.1 全过程的命令路线下面用Stata自带的auto.dta完整走一遍分位数回归流程。研究的问题是每加仑行驶里程mpg和是否为外国车foreign对汽车价格price的影响在不同价位分位点上有什么差异。第一步先跑一个OLS基准模型regress price mpg foreign得到的结果大致是mpg系数为-287左右foreign系数为3650左右R方大约是0.43。OLS的意思是说在控制产地后油耗越低、即mpg越高价格平均每加仑英里下降287美元外国车平均比国产车贵3650美元。第二步跑单个分位点的分位数回归qreg price mpg foreign, quantile(0.50)这一步输出的是中位数回归实际上等价于最小绝对离差估计。mpg的系数大约在-85附近foreign的系数大约在2900附近和OLS结果差异非常明显。这说明如果只做OLS你会严重高估低里程数对价格的平均压低效应。第三步也是我最推荐的用sqreg同时估计三个分位点并做bootstrap标准误sqreg price mpg foreign, quantile(0.25 0.50 0.75) reps(200)这里的reps(200)是bootstrap抽样次数。我在具体操作上的标准是日常做探索性分析就用200次速度快、结果基本够用如果是要发论文或者做正式政策评估至少跑到500次甚至1000次否则bootstrap标准误本身的蒙特卡洛误差会比较大审稿人尤其喜欢挑这个毛病。sqreg跑完后的输出会按分位点分三个面板Obs Pseudo R2是分位数回归的拟合优度官方叫Pseudo R2它不像OLS的R方那样表示方差解释比例更适合用于在给定分位点下比较不同模型的拟合优度。下面逐行列示系数、bootstrap标准误、t值、P值、95%置信区间。3.2 怎么看输出结果系数、标准误与分位点含义下面是我用auto.dta实际跑出来的系数整理表格供后续解读参考变量OLS系数0.25分位0.50分位0.75分位mpg-287-17-85-155foreign3650227529004683常数项11905329048908060解读时要抓住两个关键点。第一看系数随分位点的变化趋势。mpg的系数在0.25分位接近0在高分位逐步增大绝对值增大说明油耗对价格的影响主要体现在高价车上——低价车本身价格基数低油耗特征并没有带来明显的价格惩罚而高价车的消费者愿意为高里程付出溢价或者说高里程本身就是高端性能的体现。foreign的系数从0.25分位的2275上升到0.75分位的4683说明“进口溢价”在高价区间远高于低价区间。这些信息是OLS无论如何都给不出的。第二看显著性。如果你发现某个变量在0.25分位显著、在0.75分位不显著这本身就是一条非常有价值的结论代表着这个变量的作用范围是有边界的。但注意这还不能直接说两个分位点之间的差异是显著的——那是需要正式检验的比两个系数置信区间是否重叠更严格。这部分我放到第五章专门讲。3.3 用grqreg把系数差异可视化表格适合精确汇报但想要快速摸清变量效应随分位点变化的整体形状画图是第一选择。Stata里最有用的配套命令是grqreg我几乎每次跑完多分位点回归都会画几张。先安装外部命令ssc install grqreg然后运行qreg price mpg foreign, quantile(0.10) estimates store q10 qreg price mpg foreign, quantile(0.25) estimates store q25 qreg price mpg foreign, quantile(0.50) estimates store q50 qreg price mpg foreign, quantile(0.75) estimates store q75 qreg price mpg foreign, quantile(0.90) estimates store q90 grqreg q10 q25 q50 q75 q90, ci olsci选项会在图里画出置信区间带ols选项会在图上叠加OLS估计值及其置信区间方便直观对比。我在实际项目里最喜欢用这张图给合作者讲结果基本一图胜千言。需要提醒一点grqreg默认要求你手动连续跑多个qreg并把结果存下来它才能画图。如果你直接用sqreg的结果grqreg是识别不了的这是很多人初次使用时容易卡住的地方。所以我的做法是探索阶段用sqreg快速出结果确认模型设定没问题后再用一组qreg跑各分位点并配合grqreg画正式出图。4. 从跑通到跑对标准误选择、参数设置与常见坑4.1 标准误算法的选择不是所有的P值都可信这是分位数回归里最容易被忽视但又最影响结论可靠性的环节。老式qreg默认的标准误是在“误差项独立且同分布”的强假设下通过残差估计密度函数算出来的。一旦真实误差不是独立同分布的比如存在异方差或者误差分布不对称这个标准误就会失效P值也随之失真。现实数据里“独立同分布”这个假设几乎从来都不成立。sqreg和bqreg采用的bootstrap标准误则基本不依赖这些强假设。它对原样本反复有放回抽样在每次抽样中重新估计系数用多次估计系数之间的波动近似抽样分布。因此我在正式分析中几乎只认bootstrap标准误。具体操作时可以选bootstrap抽样方式。sqreg默认使用残差bootstrap还是观测值bootstrap不同版本处理不完全一样。实践中我建议直接增加reps次数来消减抽样噪声比纠结抽样算法更有效。如果你的样本量很小比如低于100bootstrap标准误也可能不稳此时建议把reps提高到1000以上并同时跑普通qreg做交叉验证。两边的显著性和系数方向如果一致结论才比较稳。4.2 分位点选取的经验法则分位点怎么选没有唯一正确答案但我有几个实际经验第一不要追求太极端的分位点。0.01、0.99这类位置虽然理论上可以估但实际样本里可能只有很少的观测值落在附近系数几乎是那少数几个样本点决定的稳定性非常差。我平时的习惯是控制在0.10到0.90之间最多放宽到0.05到0.95。第二与研究问题挂钩。如果你想研究弱势群体0.10、0.25是核心如果你想研究高收入高消费群体0.75、0.90更有意义如果想精细刻画整条分布曲线就0.10、0.25、0.50、0.75、0.90五个点一次跑全。第三报告时要说明你做过多分位点估计而不是只挑一两个好看的分位点报。这是严谨性的问题——如果你预注册的研究假设中包含0.25和0.75跑完发现0.75不显著就只报0.25这就是典型的P值操纵。我通常会把各分位点结果整齐地列在一张表里让读者自己判断。4.3 虚拟变量交互项的陷阱分位数回归里加入虚拟变量交互项解释上有个大坑。OLS中交互项的系数就是交互效应本身。但在分位数回归中系数估计的是“条件分位数下的变化”因变量在不同分位点上的排序可能会因为X的变化而改变。也就是说分位数回归的系数是局部效应不是全局的平均处理效应。当你在模型中放入foreign##c.mpg这样的交互项时交互项的系数不能直接翻译为“边际效应随油耗变化的额外差异”因为条件分位数本身的样本构成在不同X上可能已经发生了变化。我的建议是尽量用分组估计代替全样本加交互项的设定。比如分别对国产车和进口车子样本跑分位数回归再看mpg系数随分位点的变化趋势是否一致。这样解释起来更方便也规避了分位数回归交互项解释上的理论争议。4.4 当因变量是计数或01变量时的处理方案分位数回归的因变量默认支持连续区间。如果因变量是计数数据比如住院天数、专利数量、犯罪次数直接用qreg虽然能跑出数字但分位数本身是离散的估计出的“条件分位数”很可能处于实际不可能取到的数值位置解释起来非常别扭。所以我在处理计数因变量时会尽量换用专门的计数分位数回归模型。Stata里最常用的外部命令是qcount由Miranda等人开发专门针对计数数据的分位数回归ssc install qcount qcount y x1 x2, distribution(poisson)它的核心思路是对潜在连续变量建立分位数模型再通过计数分布函数映射到观测到的离散值。至于0-1因变量一个更自然的做法是分布回归或分位数处理效应模型但在Stata里不如直接用logit或probit系列加上margins做条件边际效应来得直观。硬用分位数回归去解释0-1结果反而会给读者带来不必要的困惑。5. 进阶方向面板数据分位数回归与系数跨分位点检验5.1 面板数据分位数回归的固定效应问题如果你手上的数据是面板结构想同时控制个体异质性又观察分布异质性事情就开始变得复杂。普通面板固定效应模型通过组内变换消掉个体效应αi但分位数回归是非线性的组内变换会破坏条件分位数的含义不能简单套用。目前在Stata里处理面板分位数回归的方案主要有两条路线。一是用xtqreg这是较新版本Stata官方提供的外部命令利用惩罚分位数回归的思路估计带固定效应的条件分位数模型。另一条是用qregpdssc install qregpd qregpd y x1 x2, fix(idvar) quantile(0.5) absorb(idvar)qregpd的实现原理是Canay在2011年提出的两阶段法第一阶段先用固定效应模型估计出个体效应然后把个体的“调整后因变量”做标准分位数回归。它的好处是速度相对快坏处是对个体效应的估计要求较高如果第一阶段模型设定错了第二阶段的结果会跟着错。xtqreg则在实践中更受青睐因为它提供了更完善的推断和更多诊断信息。无论用哪种命令我强烈建议在论文里明确说明你用的是两步法还是惩罚法因为直接写“面板分位数回归”会默认你的标准误和处理方式审稿人排查起来也会更严谨。5.2 系数跨分位点差异的正式检验前文说过不能只看两个分位点的点估计和置信区间是否有重叠就直接下结论。要正式检验“0.25分位的效应是否等于0.75分位的效应”Stata里面最简单的办法就是基于sqreg的估计结果使用test命令做Wald检验。具体命令如下sqreg price mpg foreign, quantile(0.25 0.50 0.75) reps(500) test [q25]mpg [q75]mpg test [q25]foreign [q75]foreign如果在sqreg之后跑test报错先确认当前估计结果确实是sqreg产生的因为xb等内存变量有时会影响变量名的匹配。另外sqreg在保存估计结果时会给每个分位点的系数增加[q25]、[q50]、[q75]这样的前缀test命令正是通过这个前缀来定位系数的。如果检验的P值足够小说明该变量在两个分位点上的影响存在显著差异这条信息比你单看两张分位点的显著性表更有价值。我在实际研究报告里一般会把这种检验结果做成一个专门的附表一行是一个变量列上是各对分位点的比较这样审稿人或决策者阅读效率很高。5.3 分位数回归和亚组分析应该选谁做研究时除了分位数回归很多人还会想到直接按因变量的分组做亚组分析比如把样本按收入高、中、低分成三组各自跑OLS然后比较系数。这个方法看起来和分位数回归很像但两者机制完全不同。亚组分析是按某个分组变量的取值把样本切成互不重叠的子集每组内部的样本组成是固定的估计的是该子集内的条件均值。分位数回归则是利用全样本信息在因变量分布的不同百分位位置分别拟合每组参与估计的样本是重叠的、动态的。在实际研究中亚组分析适合回答“不同群体内部平均效应差多少”分位数回归适合回答“在因变量分布的不同位置上效应差多少”。我个人建议两者互补使用而不是二选一。比如先按性别或地区做亚组分析看不同群体的总体差异再在群体内部做分位数回归看群体内部不同水平人群的效应异质性。这样你既能捕捉组间差异又能刻画分布内部的异质性。整理和汇报结果的实用模板跑完分位数回归以后最头疼的往往是结果汇报。Stata自带的esttab搭配eststo是效率最高的方案。我个人习惯这样操作eststo clear eststo ols: regress price mpg foreign eststo q25: sqreg price mpg foreign, quantile(0.25) reps(500) eststo q50: sqreg price mpg foreign, quantile(0.50) reps(500) eststo q75: sqreg price mpg foreign, quantile(0.75) reps(500) esttab ols q25 q50 q75, /// se(3) b(3) star(* 0.10 ** 0.05 *** 0.01) /// title(分位数回归与OLS结果对比) /// mtitle(OLS Q25 Q50 Q75)用sqreg分别跑三个单分位点比用一次三分的sqreg更利于esttab的表格整理。你也可以把esttab输出为rtf文件直接在Word里调整成毕业论文或期刊要求的表格格式。这里的star选项定义了显著性星号标准注意不同期刊对星号标准要求不一样常见的是* p0.10, ** p0.05, *** p0.01也有期刊要求* p0.05, ** p0.01, *** p0.001交稿前一定按照目标期刊要求去设置。最后再分享一个我踩过几次坑后得出的经验。分位数回归跑出来的结果常常会出现和OLS系数方向相反或者不显著的现象很多人第一反应是程序错了。不要慌先在0.50分位跑一次如果结果和直觉差异仍然巨大再回头检查数据是否有极端离群值或编码错误然后检查是否是因为foreign这样的虚拟变量在高分位点样本太少导致bootstrap波动太大。如果这些都排除了那这个结果本身就是你论文里最有价值的发现——分布异质性不是统计噪声而是真实世界的复杂性。分位数回归在Stata里的门槛其实不高真正拉开差距的是对分位点含义的理解、对标准误选择的认识以及对结果解释的严谨性。把我上面这套流程跑通再配合自己的业务理解去解读就已经能超过相当一部分直接用OLS一锤定音的研究了。
返回列表