分位数回归

第16章 分位数回归 方匡南朱建平姜叶飞传统的线性回归模型描述了因变量的条件均值分布与自变量X的关系，为了和分位数回归相区别，因此把传统的回归又称为均值回归（Mean Regression）。其中，OLS是估计回归系数的昀基本方法。如果模型的随机误差项来自均值为零、方差相同的分布，那么回归系数的OLS估计为昀佳线性无偏估计（BLUE）；如果随机误差项是正态分布，那么回归系数的OLS估计与MLE估计一致，均为昀小方差无偏估计（MVUE）。此时它具有无偏性、有效性等优良性质。  但实际中，假设不能够满足时，为了弥补普通昀小二乘法（OLS）在回归分析中的缺陷，1818年Laplace提出了中位数回归（Median Regression），利用昀小绝对偏差估计（Least absolute deviance, LAD）。在此基础上，1978年Koenker和Bassett把中位数回归推广到了一般的分位数回归（Quantile Regression）上。分位数回归是估计一组回归变量X与被解释变量Y的分位数之间关系的建模方法。本章主要介绍基本的分位数回归及其应用。 第16章  分位数回归 16.1问题的提出例16-1：恩格尔定律 德国统计学家恩格尔 ( Engel)使用收集的235个比利时家庭的收入与食物支出数据得 出其著名的恩格尔定律其著名的恩格尔定律:收入越高的家庭将其收入用于食物支出的比例越低。 No.收入消费No.收入消费No.收入消费1420.1577255.839411502.8390338.001421692.6397466.95832541.4117310.958712616.7168412.361322997.8770543.39693901.1575485.680013790.9225520.000623506.9995317.71984639.0802402.997414555.8786452.401524654.1587424.32095750.8756495.560815713.4412512.720125933.9193518.96176945.7989633.797816838.7561658.839526433.6813338.00147829.3979630.756617535.0766392.599527587.5962419.64128979.1648700.440918596.4408443.558628896.4746476.320091309.8789830.958619924.5619640.116429454.4782386.3602101492.3987815.360220487.7583333.839430584.9989423.2783表 16-1  恩格尔定律的部分数据 第16章  分位数回归 16.1问题的提出install.packages(quantreg)library(quantreg)data(engel)attach(engel)hist(foodexp)curve(density(foodexp),add=T)plot(income,foodexp,xlab=HouseholdIncome,ylab=FoodExpenditure,type=n,cex=.5)points(income,foodexp,cex=.5,col=blue) 图16-1是家庭消费支出的核密度函数图，从该图可以看出消费支出不符合正态性假设，是一个右偏的分布。另外，从图16-2收入消费散点图来看，消费支出和收入之间存在着异方差，即随着收入的增加，消费之间的差异在扩大。对于这样的数据，如果直接利用传统的均值回归方法会有问题，此处利用分位数回归是比较合适的 图16-1 消费核密度函数图 图16-2 收入消费散点图 16.2总体分位数和总体中位数第16章  分位数回归 定义16-1：对于一个连续随机变量y，其总体第τ分位数()yτ的定义是：y小于等于()yτ的概率是τ，即 ()()()()PyyFyτττ=≤=其中()P表示概率，()()Fyτ表示y的累积分布函数(cdf)。 比如(0.25)3y=，则意味着3y≤的累积概率是0.25，即(3)0.25Py≤=。并且，1()()()yFyττ−=，即()()Fyτ的反函数是()yτ。 当0.5τ=时，()yτ是y的中位数。0.75τ=时，()yτ是y的第3/4分位数，0.25τ=时，()yτ是y的第1/4分位数。若y服从标准正态分布，则(0.5)0y=，(0.95)1.645y=，(0.975)1.96y=。 16.2总体分位数和总体中位数第16章  分位数回归 另外，如果随机变量    分布是对称的，那么其均值与中位数是相同的。当其中位数小于均值时，分布是右偏的。反之，分布是左偏的。一般来讲，工资的分布是右偏的（如图16-3），所以如果单纯以平均工资来反映工资的话，这是很不恰当的，因此美国等一些国家除了公布平均工资外，还会同时公布工资的中位数和1/4、3/4分位数等。 y图16-3  美国2012家庭年收入分布 注：该图来自美国劳工与统计调查局 16.3样本（经验）分位数估计第16章  分位数回归 对一个离散的随机变量y，取其容量为T的样本序列(1,,TyyK)，计算第τ分位数的方法如下： 首先将数据从小到大排序，标号为i，1,2,,iT=K。然后利用表16-2所列的方法计算随机变量y的第τ分位数的排列序号的i；如果i为整数，则随机变量y的第τ分位数即为yi，如果i不是整数，则随机变量ty的第τ分位数为：()[][]1[]([])()iiiyyiiyyτ+=+−− 其中[i]表示不大于i的最大整数。给定一个具体的随机变量y，对于一个容量为T的样本，则y的第τ分位数的序号i的计算方法如下。在大样本情况下，各方法收敛到同一值。 连续样本的经验分位数利用表16-3的方法计算。 16.3样本（经验）分位数估计第16章  分位数回归 quantile(x, probs = seq(0, 1, 0.25), na.rm = FALSE, names = TRUE, type = 7, ...) 其中，type是1到9的取值，代表不同的经验分位数的算法，默认是第7种算法，具体的各种算法如下表。 （1）离散样本分位数 对于types 1, 2 and 3, Q[i](p)是关于p的离散函数,当i = 1和2时， m = 0，当i = 3时， m = -1/2. Type1经验分布函数的反函数.假如g=0,则γ=0；g取其他值时，γ取1。Type2与Type1类似，但是在非连续处取均值，当g=0,则γ=0.5；g取其他值时，γ取1。Type3SAS定义方法：昀近的偶数顺序统计量。假如g=0,则γ=0并且j是偶数，γ=0；其他情况，γ取1。表 16-2     离散样本经验分位数 第16章  分位数回归 （2）连续样本分位数 对于types 4-9, Q[i](p)是关于p的连续函数, 以及对应的gamma = g和 m 详见下表。样本（经验）分位数可以通过点(p[k],x[k])之间进行线性插值得到，其中x[k]是第k个顺序统计量。关于p[k]的具体表达式详见下表。 Type4m=0，p[k]=k/n.也就说经验分布函数的线性插值。Type5m=1/2，p[k]=(k-0.5)/n.这是一个分段线性回归函数。水文研究比较常用该方法。Type6m=p，p[k]=k/(n+1).因此，p[k]=E[F(x[k])].MinitabandSPSS用这种方法Type7m=1-p，p[k]=(k-1)/(n-1).此时,p[k]=mode[F(x[k])].S语言使用此方法.Type8m=(p+1)/3，p[k]=(k-1/3)/(n+1/3).则p[k]=~median[F(x[k])].这个分位数估计方法近似与中位数无偏，而不管x的分布。Type9m=p/4+3/8，p[k]=(k-3/8)/(n+1/4). 当x是正态分布是，该估计结果对期望顺序统计量是近似无偏的。表 16-3    连续样本经验分位数 quantile(x,probs=c(0.1,0.25,0.5,0.75),type=2)10%25%50%75%-1.24081576-0.682396630.062314510.71881867我们利用quantile()函数求经验分位数 16.3样本（经验）分位数估计16.4分位数回归原理第16章  分位数回归 离差绝对值              在中位数时取到昀小值。因此，中位数回归估计量可以通过昀小绝对离差法（least absolute deviation, LAD）估计 ∑−αytttyµ′=+Xβ∑′−(0.5)βXtˆtyβ(0.5)ˆβty(0.5)βXXtˆ)ˆ()5.0(=ty对于线性回归模型 。通过求 昀小，得到 的中位数回归系数估计量 ，从而得到 的中位数回归拟合值 。 16.4分位数回归原理第16章  分位数回归 现在我们把中位数回归推广到分位数回归。对于回归模型，被解释变量y对以X为条件的第τ分位数用函数()|yXτ表示，其含义是：以X为条件的y小于等于()|yXτ的概率是τ，即()()(|)(|)pyyXFyXτττ≤==，或者可以写成1()()|(|)ttyXFyXττ−= 其中()(|)FyXτ和1()(|)tFyXτ−分别是y在给定X条件下的累积概率分布函数(cdf)和其反函数。则()|yXτ称作被解释变量ty对X的条件分位数函数。而()()(|)(|)ttFyXfyXττ′=则称作分位数概率密度函数。 16.4分位数回归原理第16章  分位数回归  Koenker和Bassett(1978)证明，若用()ˆyτ表示y的τ分位数回归估计量，则对于以检查函数（check function）wτ为权数，y对任意值α的加权离差绝对值和∑−ατtyw只有在()ˆ=yτα时取得最小值。其中         ::(1)()()iiTTtttiytywyyyτααατατα≥−=−−−+−∑∑∑     (16-1) 其中(0,1)τ∈。因此，分位数回归可以通过加权的最小绝对离差和法（weighted least absolute deviation, WLAD）进行估计。 16.4分位数回归原理第16章  分位数回归 根据式(16-1)，对于线性回归模型tttyµ′=+Xβ，求第τ分位数回归方程系数的估计量)(ˆτβ的方法就相当于是求使得下式（目标函数）达到最小时的解， ()()()()()()ˆˆ00()()ˆˆ::ˆˆ(1)ˆˆ(1)()()ttttTTttuuTTtttyXtyXQuuyyττττττττββττττ≥′′≥=−−+′′=−−−+−∑∑∑∑ttXβXβ(16-2) 其中，tu)(ˆτ表示第τ分位数回归方程对应的残差。 第τ分位数的回归方程表达式是         ()()ˆˆtyττ′=tXβ 16.4分位数回归原理第16章  分位数回归 当0.5τ=时，式（16-2）变为         ∑∑∑=′≥′′−=′−+′−−=TttTXyttTXyttyyyQtt1)0.5(ˆ:)0.5(ˆ:)0.5(ˆ0.5)ˆ(0.5)ˆ(0.5)0.5()0.5(βXβXβXββ  (16-3) (0.5)(0.5)ˆˆty′=tXβ称作中位数回归方程，)0.5(ˆβ是中位数回归系数估计量。