分类目录归档:回归分析

回归分析、线性模型、非线性模型、广义线性模型

LARS算法简介

最近临时抱佛脚,为了讨论班报告Group Regression方面的文章,研究了Efron等人于2004年发表在Annals of Statistics里一篇被讨论的文章LEAST ANGLE REGRESSION。这篇文章很长,有45页。加上后面一些模型方面大牛的讨论的文章,一共有93页。对于这种超长论文,我向来敬畏。后来因为要报告的文章里很多东西都看不懂,才回过头来研读这篇基石性的文章。

所谓大牛,就是他能提出一种别人从来没有提出过的想法。大牛们看待问题的角度和常人不同。比如在回归中常用的逐步回归法。我们小辈们只知道向前回归,向后回归还有二者结合的一些最基本的想法。比如向前回归,就是先选择和响应最相关的变量,进行最小二乘回归。然后在这个模型的基础上,再选择和此时残差相关度最高的(也就是相关度次高)的变量,加入模型重新最小二乘回归。之后再如法继续,直到在某些度量模型的最优性准则之下达到最优,从而选取一个最优的变量子集进行回归分析,得到的模型是相比原模型更加简便,更易于解释的。这种方法,牺牲了模型准确性(预测有偏),但是提高了模型的精确度(方差变小)。大多数本科生对逐步回归的理解也就如此了。Efron看待这个问题时,比起常人更高了一个层次。他首先指出,逐步向前回归,有可能在第二步挑选变量的时候去掉和X1相关的,但是也很重要的解释变量。这是因为它每次找到变量,前进的步伐都太大了,侵略性太强。

因此在这个基础上,Efron提出了Forward stagewise。也就是先找出和响应最相关的一个变量,找到第一个变量后不急于做最小二乘回归,而是在变量的solution path上一点一点的前进(所谓solution path是指一个方向,逐步回归是在这个方向上进行),每前进一点,都要计算一下当前的残差和原有的所有变量的相关系数,找出绝对值最大的相关系数对应的变量。我们可以想像,刚开始,前进的步伐很小,相关系数绝对值最大的对应的变量一定还是第一步选入的变量。但是随着前进的进程不断向前,这个相关系数的绝对值是在慢慢减小的,直到找到另外一个变量X2,它和当前前残差的相关系数和第一个入选变量X1的相关系数绝对值相同,并列第一。此时把X2也加入回归模型中,此时回归模型在X1上的系数已经确定了,如果在X1的solution path上继续前进,则得到的与当前残差相关系数最大的变量一定是X2,所以不再前进,而是改为在X2的solution path上前进,直到找到第三个变量X3,使得X3的与当前残差的相关系数绝对值最大。这样一步一步进行下去。每一步都是很多小步组成。直到某个模型判定准则生效,停止这个步骤。在每一个solution path上的计算都是线性的。总体的solution path是分段线性的。这种算法是一种自动进行模型构建的方法。它和传统的Forward selection在本质上是一样的,都是选择一个变量,然后选择一个继续进行的solution path,在该方向上前进。这两种方法的solution path的选择方法是一样的,唯一的区别就是前进的步伐不一样,Forward selection的前进步伐很大,一次到头,而stagewise则是一小步一小步前进。这样比Forward selection要谨慎一些,会免于漏掉一些重要的变量。

从这个视角来看,我们可以选择另外一种solution path。Efron等人在这篇文章中,就提出了一种新的solution path。在已经入选的变量中,寻找一个新的路径,使得在这个路径上前进时,当前残差与已入选变量的相关系数都是相同的。直到找出新的与当前残差相关系数最大的变量。从几何上来看,当前残差在那些已选入回归集的变量们所构成的空间中的投影,是这些变量的角平分线。下面我简单的描述一下这个算法:

  • 第一步,我们初始的估计模型为0,那么当前的残差就是Y,我们找出X’Y中绝对值最大的那个对应的变量,记为X1,把它加入回归模型。这一步中X’Y是当前残差和所有变量的相关系数向量。(注意这里Y都已经中心化,X中心标准化过了)。
  • 第二步,在已选的变量的solution path上前进,solution path就是s1*X1,s1是X1与当前残差的相关系数的符号。在这个path上前进,直到另外一个变量出现,使得X1与当前残差的相关系数与它和当前残差的相关系数相同。记这个变量为X2,把它加入回归模型中。
  • 第三步,找到新的solution path。Efron在文章中提出了一种找出满足LARS条件的solution path的解法。solution path需要使得已选入模型变量和当前残差的相关系数均相等。因此这样的路径选择它的方向很显然就是$X_k(X_k’X_k)^{-1}1$的指向(因为$X_k'(X_k(X_k’X_k)^{-1})1$的元素都相同,保证了LARS的要求,当然这里或许会有一些其他的解,也能满足LARS的要求,有没有达人能想到或许证明这个解是唯一的)。只要再标准化这个向量,我们便就找到了solution path的方向。在这个方向上前进,直到下一个满足与当前残差相关系数绝对值最大的变量出现。如此继续下去。

LARS算法,保证了所有入选回归模型的变量在solution path上前进的时候,与当前残差的相关系数都是一样的。这一点,比起Forward stagewise要捷径一些,走得更快一些。

LARS算法已经在SAS和R中实现了。作为回归模型选择的一种重要的算法,LARS相比起传统的Forward selection和Forward stagewise,既不那么富于侵略性,又比较走捷径。LARS算法在lasso 估计的求解中也有非常好的应用。在Efron等人的同篇论文中有详细的讨论。关于lasso和它的LARS算法,笔者将在今后的文章中介绍。

 

从线性模型到广义线性模型(2)——参数估计、假设检验

1.GLM参数估计——极大似然法

为了理论上简化,这里把GLM的分布限定在指数分布族。事实上,实际应用中使用最多的分布就是指数分布族,所以这样的简化可以节省很多理论上的冗长论述,也不会限制实际应用。
如前文如述,指数分布族的概率密度函数可以统一地写为:

$f_Y(y;\theta,\Psi)=exp[(y\theta – b(\theta))/{\Psi} + c(y;\Psi)]$

这里为了在模型中体现散布参数(dispersion parameter)$\phi$,把上述密度函数中的$\Psi$记做
$\Psi=a_i(\phi)={\phi}/w_i$
从而响应变量的单个观测值的(加权)对数似然函数可以表示为:

$logL({\theta}_i,\phi;y_i)=w_i[(y_i{\theta}_i-b({\theta}_i))/{\phi}]+c(y_i,\phi)$

再结合观测值之间的独立性,全体观测值的对数似然函数可记做:$\sum_i logL({\theta}_i,\phi;y_i)$
一般情况下最大化上述的对数似然函数很难找到解析解(正态分布是特例之一),因而必须使用数值方法求解。McCullagh和Nelder(1989)证明了使用Newton-Raphson方法,结合Fisher scoring算法,上述对数似然函数的最大化等价于连续迭代的加权最小二乘法(iteratively weighted least squares, or IRWLS)。

广义线性模型的IRWLS算法如下:
1.设置线性估计量和响应变量的均值的初始估计值: $\hat {\eta}_0$和$\hat {\mu}_0$
这里$\hat {\mu}_0$是根据经验或是专家意见等信息对$\mu=E(Y)$的一个估计值,而$\hat {\eta}_0$可以利用模型建立时选用的联接函数来获得,即$\hat {\eta}_0=g(\hat {\mu}_0)$。这一函数关系也用于计算步骤2和3中$\eta$对$\mu$一阶导数。
2.构造调整的因变量(adjusted dependent variable):$z_0=\hat {\eta}_0+(y-{\hat \mu}_0){d\eta \over d\mu}|_{\hat {\eta}_0}$
3.构造权重:$w^{-1}_0={({d\eta \over d\mu})}^2|{\hat {\eta}_0V(\hat {\mu}_0)}$
这里$V(\hat {\mu}_0)$是利用方差函数(variance function)和$\hat {\mu}_0$构造的$Var(Y)$的估计值。
4.利用步骤2和3构造的调整的因变量和权重,拟合普通线性模型(ordinary linear model),预测/拟合(predict)新的线性估计量和均值: $\hat {\eta}_1$和$\hat {\mu}_1$
5.重复步骤2-4直到收敛(满足一定的迭代步数或是精度要求)。
此时得到的模型就是极大似然估计方法下的广义线性模型。IRWLS的算法思路也从另一个方面说明了广义线性模型是普通线性模型的推广。在广义线性模型的实际应用中,IRWLS算法是最常用的极大似然估计求解方法。对于特殊的案例,也有其他的特殊的参数估计方法。比如对于在精算学科中最常用的列联表(contigency table)数据或案例就有Bailey-Simon法、边际总和法(marginal totals)、最小二乘法(least squares)、直接法(direct method)等。

2.假设检验

2.1 空模型和全模型

一个极端的情况,所有自变量$x_i$对于响应变量$Y$都没有影响,也即是为所有的响应变量$Y$拟合一个共同的均值,即只有一个参数。这样的模型称为空模型(null model)。对于普通线性模型(正态分布下的GLM)而言,空模型的具体形式就是$y=\mu + \epsilon$。对于特殊的数据或案例类型,可能存在着其他的限制条件(constraints)从而空模型的参数个数大于1。比如非寿险精算中经常用到的列联表(contigency table)数据,其空模型就可能包含了行号、列号、对角线序号等限制。

相反的一个极端情况就是,所有自变量$x_i$的每一个观测值或称为数据的样本点(data points)对于响应变量$Y$都有影响,这样的模型称为全模型(full or saturated model)。一般可以通过构造阶数足够高的多项式或者把所有的量化观测值(quantitative)视为质化观测值(qualitive),并且引入适当数量的交叉项(interactions)来构造全模型。

统计建模的目的之一就是把样本数据划分为随机成分和系统成分两大部分。在这一点上,空模型认为响应变量的变动(variation)完全由随机性(random variation)造成,而全模型则认为响应变量的变动完全来自于系统成分(systematic)。一个直观地理解就是全模型是在现有的数据或样本的条件下,针对某一种分布所能拟合的最优模型,因而可以做为检验目标模型拟合优度的一个标准(measure)。

2.2 偏差(Deviance)

如果把全模型的对数似然函数记为$l(y,\phi|y)$,把目标模型的对数似然函数记为$l({\hat {\mu}},\phi|y)$,那么目标模型与全模型在拟合优度上的偏离的定义可写成$2(l(y,\phi|y)-l({\hat {\mu}},\phi|y))$。再结合观测值的独立性假设和指数散布族的假设,那么上述偏离的定义可以简化为:

$\sum_i 2w_i(y_i({\hat {\theta}_i} – {\tilde {\theta}_i}) – b({\tilde {\theta}_i}) + b({\hat {\theta}_i})) /{\phi}$

其中$a_i(\phi)={\phi}/w_i$,$\tilde {\theta}$是全模型下的参数估计值,$\hat {\theta}$是目标模型下的参数估计值。如果把上式写成$D(y,\hat {\mu})/{\phi}$,那么$D(y,\hat {\mu})$称为偏差(Deviance),$D(y,\hat {\mu})/{\phi}$则称为标准化偏差(scaled deviace)。
此外,皮尔逊卡方统计量(Pearson’s chi-square statistics):

$X^2={\sum_i (y_i – {{\hat \mu}_i})^2 \over Var({\hat {\mu}}_i)}$

也是衡量模型偏离程度(discrepancy)的统计量之一,在一些场合可以做为偏差的替代选择。

2.3 拟合优度检验

广义线性模型的假设检验可以分为两种:一是检验目标模型相对于数据或预测值的拟合有效性的检验(goodness of fit test);另外一种则是对“大”模型以及对“大”模型的参数施加一定的线性约束(linear restrictions)之后得到的“小”模型之间的拟合优度比较检验。直观上的理解就是,“大”模型具有更多的参数,即从参数的线性约束总可把一个或多个参数用其他参数的线性组合来表示,然后代入“大”模型,从而参数的个数减少,派生出所谓的“小”模型,也就是说“大”和“小”并非任意的,而是具有一种派生关系(nested models)。如果把全模型认为是“大”模型,而目标模型是“小”模型,那么上述两种检验的本质是相同的。因而假设检验的零假设(null hypothsis)可以统一且直观地设定为:“小”模型(目标模型)是正确的模型。

如果把大模型记做$\Omega$,把小模型记做$\omega$,其标准化偏差之差记做$D_{\omega} – D_{\Omega}$,其自由度之差记做$df_{\omega}-df_{\Omega}$,则构造如下的统计量:${(D_{\omega} – D_{\Omega})/(df_{\omega}-df_{\Omega})} \over {\phi}$。

当$\phi$是已知常数时,比如泊松和二项分布的情况下$\phi=1$,上述统计量在零假设下渐近地(asymptotically)服从卡方分布(正态分布时正好是卡方分布)。当$\phi$未知时,通常需要用估计值代替。最常用的估计值是$\hat {\phi}=X^2/(n-p)$这里n是数据中观测值的数量,p是目标模型的参数个数。此时上述的统计量在零假设下近似地(approximately)服从F分布(正态分布时严格服从F分布)。注意上述两种情况下,渐近和近似的区别。

对于某一个参数,可以使用其估计值的标准误(standard error)来构造一个z统计量来检验其显著性,即$z=\hat {\beta}/se(\hat {\beta})$。在零假设下,z统计量在普通线性模型,也就是正态分布下的广义线性模型中就是我们熟知的t统计量,严格服从t分布。在其他分布下的广义线性模型中,渐近地服从正态分布。z检验也称为Wald检验,在广义线性模型中效果不如上述的偏差检验,因而较少使用。

从线性模型到广义线性模型(1)——模型假设篇

在统计学里,对特定变量之间的关系进行建模、分析最常用的手段之一就是回归分析。回归分析的输出变量通常记做$ Y$,也称为因变量(dependent)、响应变量(response)、被解释变量(explained)、被预测变量(predicted)、从属变量(regressand);输入变量通常记做$ x_1$,…,$x_p$,也称为自变量(independent)、控制变量(control&controlled)、解释变量(explanatory)、预测变量(predictor)、回归量(regressor)。本文根据作者自己的一些学习心得和理解,简单且不严格地介绍在模型假设方面普通线性模型和广义线性模型的区别和联系/推广(generalization)。广义线性模型的拟合检验、推断、诊断等方面的方法和手段依赖于模型所采用的分布类型,难以一概而论,将在作者后续的学习心得文章里具体介绍。

1.普通线性模型的简单回顾

普通线性模型(ordinary linear model)可以用下式表示:

$ Y={\beta}_0+{\beta}_1x_1+{\beta}_2x_2+…+{\beta}_{p-1}x_{p-1}+\epsilon$                          (1.1)

这里$ {\beta}_i$,$ i=1$,…,$p-1$称为未知参数,$ {\beta}_0$称为截矩项。

普通线性模型的假设主要有以下几点:

1.响应变量$ Y$和误差项$\epsilon$正态性:响应变量$ Y$和误差项$\epsilon$服从正态分布,且$\epsilon$是一个白噪声过程,因而具有零均值,同方差的特性。

2.预测量$x_i$和未知参数${\beta}_i$的非随机性:预测量$x_i$具有非随机性、可测且不存在测量误差;未知参数${\beta}_i$认为是未知但不具随机性的常数,值得注意的是运用最小二乘法或极大似然法解出的未知参数的估计值$\hat{\beta}_i$则具有正态性。

3.研究对象:如前所述普通线性模型的输出项是随机变量$ Y$。在随机变量众多的特点或属性里,比如分布、各种矩、分位数等等,普通线性模型主要研究响应变量的均值$E[ Y]$。

4.联接方式:在上面三点假设下,对(1.1)式两边取数学期望,可得

$ E[Y]={\beta}_0+{\beta}_1x_1+{\beta}_2x_2+…+{\beta}_{p-1}x_{p-1}$                                        (1.2)

从 (1.2)式可见,在普通线性模型里,响应变量的均值$E[ Y]$与预测量的线性组合${\beta}_0+{\beta}_1x_1+{\beta}_2x_2+…+{\beta}_{p-1}x_{p-1}$通过恒等式(identity)联接,当然也可认为通过形为$f(x)=x$的函数(link function)联接二者,即

$ E[Y]=f({\beta}_0+{\beta}_1x_1+{\beta}_2x_2+…+{\beta}_{p-1}x_{p-1})={\beta}_0+{\beta}_1x_1+{\beta}_2x_2+…+{\beta}_{p-1}x_{p-1}$

2.广义线性模型的简单介绍

广义线性模型(generalized linear model)正是在普通线性模型的基础上,将上述四点模型假设进行推广而得出的应用范围更广,更具实用性的回归模型。

1.响应变量的分布推广至指数分散族(exponential dispersion family):比如正态分布、泊松分布、二项分布、负二项分布、伽玛分布、逆高斯分布。exponential dispersion family的详细定义限于篇幅这里不做详细介绍。

2.预测量$x_i$和未知参数${\beta}_i$的非随机性:仍然假设预测量$x_i$具有非随机性、可测且不存在测量误差;未知参数${\beta}_i$认为是未知且不具有随机性的常数。

3.研究对象:广义线性模型的主要研究对象仍然是响应变量的均值$E[ Y]$。

4.联接方式:广义线性模型里采用的联连函数(link function)理论上可以是任意的,而不再局限于$f(x)=x$。当然了联接函数的选取必然地必须适应于具体的研究案例。同时存在着与假设2.1里提及的分布一一对应的联接函数称为标准联接函数(canonical link or standard link),如正态分布对应于恒等式,泊松分布对应于自然对数函数等。标准联接函数的推导及其应用上的优点涉及到指数分散族的标准化定义,这里不做详述。

3.简单的例子

考虑这样一个简单的退保案例:一个保险产品一共卖出12份保单(当然了这在现实中不可能,这里仅为示例),在保单期限内一共有6人退保。那么采用这12个投保人的特征数据如收入、职业、年龄等做为预测变量对退保/退保率进行回归分析时,普通线性模型不再适用,因为这里退保这一事件不再服从正态分布,而是二项分布(当然了如果观测值的个数足够大,比如大于30,正态分布是一个很好的近似)。此时就可采用广义线性模型(目标分布采用二项分布)进行建模。

4.补充:指数分布族的简介

指数分布族(exponential dispersion family)实质上是对一类具有以下形式的概率密度函数或具有此类密度函数的分布的总括:

$f_Y(y;\theta,\Psi)=exp[(y\theta – b(\theta))/{\Psi} + c(y;\Psi)]$

这里$\Psi$和$\theta$是实参数,$b(.)$和$c(.;.)$是实函数,该密度函数的支集(support)$D_{\Psi}$是$R$的子集,且不依赖于$\theta$。满足$\theta=\eta=g(\mu)$的联接函数$g(\mu)$称为标准联接函数(standard or canonical link)。

一般情况下参数$\Psi$的值是未知常数(fixed and unknown),因此在许多GLM文献里指数分布族又被称为单参数指数族(one-parameter exponential family)。对于比较常用的分布,$\Psi$和$\theta$的取值具有特殊的形式:

正态分布$N(\mu,{\sigma}^2)$:$\Psi={\sigma}^2$和$\theta=\mu$

泊松分布$Poisson(\lambda)$:$\Psi=1$和$\theta=log\lambda$

二项分布$Binomial(m,p)$:$\Psi=1$和$\theta=log(p/(p-1))$

负二项分布$Negative Binomial(r,p)$:$\Psi=1$和$\theta=log(1-p)$

伽玛分布$Gamma(\alpha,\beta)$:$\Psi=1/{\alpha}$和$\theta=-{\beta}/{\alpha}$

漫谈相关与回归

老师不断提醒我要对统计学的基本概念、定义及背景反复思考,这样才不会本末倒置,迷失方向。但是这个做起来很难,因为那些概念定义等看起来实在”太简单”、”没什么东西”,可能还是不能够平心静气吧!

最近静下来看了David Freedman等著的《统计学》的”相关与回归”部分,以及一篇关于直方图的文章,不免有些感慨!其实统计学中的很多概念、工具、方法等的实际意义或作用可能要比我们认为的要大很多,同时,当我们从一些概念定义等中发现出一些新东西时我们总会欣喜若狂。世界上的很多事物又何尝不是如此,人们对事物的了解总易受到传统或他人的影响仅仅停留在表面,很少达到全面而深刻,而一旦我们获得了那种深刻的洞察力,才发现真实世界是何等的精彩!一直以为直方图很简单,无非是一些代表频数的柱状图的组合而已,感觉没什么作用,但是看了一篇关于直方图制作方面的论文时,才认识到直方图的威力。直方图其实是非参数统计中估计总体分布特征的一项重要工具,选择好适当的组距和边界点(组距和最小边界点是关键),随着样本量的增大,它可以非常接近地反映数据的真实分布情况。其实,在统计中使用一种工具方法的目的也应该是使现有的数据尽可能多地反映出真实的信息,而这项工作往往是一个无底洞(这时又要考虑到效率问题了)。

散点图亦是如此。散点图给出了所有数据点的信息,但是如何从这些数据中获得结论或拟合模型,甚至用来预测?面对一张散点图,相关和回归应该是最容易想到的吧!这里主要谈谈两个变量间的相关和回归。

在研究两个变量的关系时,一般会先看看它们的散点图,在图中两变量的关系还是比较直观的,大致可以判断是否线性相关及相关性大小如何,是否是非线性相关等。而到底什么是相关呢?相关其实就是知道一件事对了解另一件事的帮助的大小。实际中,如果对某一事物不太了解,但是对与其有一定联系的另一事物有所了解,如果这种联系很强,那我们对于那件不了解的事物就有了更多的信息,或者说对这个不了解的事物有了更大的自信去预测。其实这也是研究中的一种常用的方法。

关于两个变量间的相关系数的计算。我们都知道两个变量X与Y的相关系数的计算公式为$Cov(x,y)/(SD(x)*SD(y))$,然而这已经是一个结果性东西了,我更推崇David Freedman等著的《统计学》中计算方法:先分别对两个变量做标准化,比如对变量X做标准化$(x_i-\bar{x})/SD(x)$,然后对应的标准量相乘,最后加总再求平均. 这种求法反映到散点图中,相当于对散点图的坐标刻度标准化,从而使两个坐标轴具有了相同的刻度,同时在直观两个变量之间的相关性大小时不会受到各自的标准差大小的影响。这个新的坐标系把所有的点(数据对)分到了不同的象限,通过观察各个象限的点的个数和大致分布情况便可以对相关性的大小与正负有直观的了解,比如更多的点都分布在一、三象限且群集于一条直线周围,那么这两个变量的一般具有较强的正线性相关。

我们都知道相关系数是-1到1之间的一个实数,那么相关系数为0.8是不是表示百分之八十的点群集在一条直线的周围吗?当然不是,相关系数是基于全体数据的一个综合信息,它反映的是所有点与某一条直线的群集程度,而不是一部分的点。由此也不免想到,我们在用到一些概念或定义时,也必须清楚这个概念或定义是基于怎样的对象,或有哪些局限条件或假定,比如概率论中的”事件”, “事件”是基于特定条件的,在具体使用过程中大家对这些特定条件都太”熟悉”以至于很少关注到它们,所以一旦条件改变了,大部分的人认为还是指同一个东西。

相关与因果也是这样,我们都知道相关不能同因果划等号,但实际过程中人们总”自然而然”地得出一些结论。一般来说,体重和升高相关,那体重大是不是就因为身高高呢?除非有一个固定的身材标准,而所有人都是这个标准。(随便提个问题:如果两个变量严格线性相关,即相关系数为1,那是否可以说这两个变量中一个为因一个为果呢?)。其实两个变量的相关更经常的情况是它们同时受到另外的一个或多个因素的影响,在这里可以通过对照试验或观察研究来进一步研究。另外,相关是可逆的,而因果则不可以。所以我们分析相关时总是如此谨慎地说,某某变化,与此相关的某某”相应地”如何变化。研究发现,个人收入与教育水平相关,高教育水平是不是高收入的原因呢?实际情况是它们相互影响:教育水平高的人收入一般较高,收入高的一般也更有能力获得继续教育的机会。虽然相关不是因果,然而有时我们并不需要弄清所有的因果关系,盯住输入和输出,只要存在相关,即使不是因果关系也不妨碍人们利用这种关系来进行推断。比如利用公鸡打鸣来预报太阳升起,虽然公鸡打鸣绝对不是日出的原因(虽然打鸣发生在先)。

在对两变量的相关关系有一定了解后,接下来的自然想法便是拟合回归模型。”回归”这一词来自于高尔顿的父子两代身高的研究,身高较高的父亲其儿子的平均身高要比父亲矮些,身高较矮的父亲其儿子的平均身高要比父亲高些,用高尔顿的话说就是”回归到平常”。虽然现在统计学上的”回归”这一概念已经远远超出的当时的定义,但是回归的原始思想依然有着非常重要的作用。”回归”,个人认为其实就是向中心的回归。在知道某地区18-24岁男子的身高的大致情况时,如果没有其他信息,让我们估计该地区中某一特定区域18-24岁男子的平均身高时(当然不是侏儒或篮球运动员之类的人),自然是用平均数(包括中位数)去估计了,这便是回归,没有其它的辅助信息时我们总倾向于平均值,这当然是符合统计思想的。两个变量的相关系数绝对值为1时,那么知道一变量的值就立即知道了另一变量的值;相关系数为0时,那么知道一变量的值对预测另一变量没有任何意义,那么我们就估计其值为平均值;相关系数绝对值介于0与1之间时,相关程度越大,我们越不倾向于取平均值。其实回归模型也是基于平均意义的,让我们来看看回归的本质(暂以两个变量x和y为例),回归是对每一个x值的y的平均值的估计,所以用回归模型来预测或估计总是平均意义的(这也是回归的思想),而针对某个特别的个体的预测则就需要非常的慎重了。

有这样一个例子,某学前班在儿童入学和结业时均要做智商测验,结果发现前后两次测验的分数平均都接近于100分,标准差为15分。但是仔细观察发现入学分数低于平均值的儿童结业时分数平均提高了5分,相反入学分数高于平均值的儿童结业时分数平均降低了5分,难道学前班会使儿童的智商平均化?其实没那么夸张,这只是回归效应的一个表现,只要两次测验分数的散点图中所有点不在同一条直线(这条直线的斜率为1)上,那么就会存在回归效应。观察得到的数据并不是真实值,都有或大或小的、或正或负误差,在大多数对称的概率分布中,观察值大于平均值的往往是其真实值加上了一个正的机会误差,观察值小于平均值的往往是其真实值加上了一个负的机会误差。所以在那个学前班中,入学分数较平均分低的儿童其真实分数一般是大于观察值的,因此在结业时的分数一般是要比入学时高,因为在向观察值的平均值,即真实值回归。

相关与回归是一定范围内的相关与回归,超出范围没有任何意义(经常实践的人应该会很少犯此类毛病的吧)。回归其实并不能增加信息量,它是一种结论(结论的准确性还有待评价),或对数据以某一种方式的总结,超出范围的估计预测是没有任何意义的。收入与教育水平有关,无休止的教育显然不会带来收入的持续的增加,所以人为地改变一个变量,希望通过回归模型的”魔力”来改变另一个变量是很荒谬的。另外,变量也是有范围或区域限制的,因此在使用回归模型做预测时是要非常谨慎的。

现在研究的回归往往都是多元回归,往往比较复杂,其实这是符合实际情况的,因此往往要用多个变量作为因子来拟合,但是这些变量是不是考察某一方面的较好指标呢,比如收入与教育水平有关,还可能与父母的社会地位有关,那这个”父母的社会地位”这一因子又该如何度量呢?这又是一个问题,尽管多元回归是一种非常有用的技术,但是永远代替不了对数据间内在关系的了解。由此可见实践经验的重要性!

使用回归分析,样本过少时不妨好先作图看看

回归分析往往是学统计、学计量课程时接触的第一个统计模型了,甚至不少人可能认为回归分析理所当然成为计量的绝大部分内容——毕竟很多教材中提到统计模型的时候,往往就一个OLS为主的讲法。回归分析的内容当然很广泛,也在学科中占据相对基础的位置。

学会OLS,有人还明白了ML等方法的含义;现在学统计分析的时候,或多或少会安排统计软件的实践课程,于是大家学会了使用Excel,乃至SAS中如何来做经典的回归分析。看过不少的文献,很多都忽略了回归分析模型诊断这个环节——可能很多标准教科书没有强调,甚至是没有讲;这不能不说是一个遗憾。

回归分析使用最广泛,误用的情况也多了些。下面使用一个经典的例子,来“恶心”一下那些“过分钟爱”经典回归分析的人——我在很多课堂上都举过这个例子(Anscombe),作为从基础课程向中级乃至高级课程的开场白。

x1 x2 x3 x4    y1   y2    y3    y4
1  10 10 10  8  8.04 9.14  7.46  6.58
2   8  8  8  8  6.95 8.14  6.77  5.76
3  13 13 13  8  7.58 8.74 12.74  7.71
4   9  9  9  8  8.81 8.77  7.11  8.84
5  11 11 11  8  8.33 9.26  7.81  8.47
6  14 14 14  8  9.96 8.10  8.84  7.04
7   6  6  6  8  7.24 6.13  6.08  5.25
8   4  4  4 19  4.26 3.10  5.39 12.50
9  12 12 12  8 10.84 9.13  8.15  5.56
10  7  7  7  8  4.82 7.26  6.42  7.91
11  5  5  5  8  5.68 4.74  5.73  6.89

上面有四对x,y,分别做经典回归分析的话,结果如下:

[[1]]
Estimate Std. Error  t value    Pr(>|t|)
(Intercept) 3.0000909  1.1247468 2.667348 0.025734051
x1          0.5000909  0.1179055 4.241455 0.002169629
[[2]]
Estimate Std. Error  t value    Pr(>|t|)
(Intercept) 3.000909  1.1253024 2.666758 0.025758941
x2          0.500000  0.1179637 4.238590 0.002178816
[[3]]
Estimate Std. Error  t value    Pr(>|t|)
(Intercept) 3.0024545  1.1244812 2.670080 0.025619109
x3          0.4997273  0.1178777 4.239372 0.002176305
[[4]]
Estimate Std. Error  t value    Pr(>|t|)
(Intercept) 3.0017273  1.1239211 2.670763 0.025590425
x4          0.4999091  0.1178189 4.243028 0.002164602

这时候你会发现p值、s.e.值都好的很,截距项和斜率项也非常好,甚至连R square都简直一模一样。

往往讲到这个时候,机灵点的学生就开始皱眉头了:数据差别挺大的,怎么模型都是一样呢!?

学过回归诊断的人立刻就嚷嚷来做个回归诊断来看看,四种情况下是不是有的属于“模型前提假设”不满足的情况。(当然,似乎本科阶段,教师没有过分强调这点的,很少有学生能自己自主的想到这点,所以往往让人感叹教学过程中的不完善。)

紧接着我就给他们做出一幅图:

anscombe

这时候他们有点明白,能回答说左上是我们“常见”的回归;右上完全就是二次曲线嘛;左下存在一个“异常点”;而右下,完全就属于“诡异”情形。

所以可以归纳一下:做统计模型的时候,一定记得要模型诊断下——重点就分析下残差是不是符合假设——这个是标准套路。如果你想“偷懒”,在只有一个自变量的时候,不妨在样本量不多(现在的统计软件已经非常强大了,数万对点很容易搞定;如果数据太多,你也可以抽样来试试嘛。)这样,“一眼”就能看出个大概了——这恐怕就是人比机器厉害的地方,人比模型厉害的地方。