正在加载
请稍等

菜单

Home 码农菜园 数据分析 R学习笔记-8 回归
Home 码农菜园 数据分析 R学习笔记-8 回归

R学习笔记-8 回归

数据分析 by   阅读量 13,354

回归分析是统计学的核心,通常指那些用一个或多个预测变量(自变量或解释变量)来预测响应变量(因变量、效标变量或结果变量)的方法。

8.1 回归的多面性

回归的种类很多(简单线性、多项式、多元线性、多变量、Logistic、泊松、Cox比例风险、时间序列、非线性、非参数、稳健),R中也提供了200多种用于拟合回归模型的函数,而我们重点研究对象是普通最小二乘(OLS)回归,包括简单线性回归、多项式回归和多元线型回归。

8.2 OLS回归

OLS回归拟合基于最小二乘法,关于最小二乘法的原理此处不赘述。

为了能够恰当地解释OLS模型的系数,数据必须满足以下统计假设:

  1. 正态性:对于固定的自变量值,因变量值呈正态分布;
  2. 独立性:各观测值之间相互独立;
  3. 线性:因变量和自变量之间为线性相关
  4. 同方差性:因变量的方差不随自变量的水平不同而变化。

拟合回归模型

在R中,拟合线性模型最基本的函数就是lm(),格式为:

formula为要拟合的模型形式,data为数据框,包含用于拟合模型的数据。拟合结果(myfit)存储至一个列表中,包含了所拟合模型的大量信息。formula形式如下,左边为响应变量,右边为各个预测变量。

除了lm()外,以下列出了其他一些对简单或多元回归分析有用的函数。拟合模型后,将这些函数应用于lm()返回的对象,可以得到更多额外模型信息。

  • summary():展示拟合模型的详细结果;
  • coefficients():列出拟合模型的模型参数(截距和斜率);
  • confint():提供模型参数的置信区间(默认95%);
  • fitted():列出拟合模型的预测值;
  • residuals():列出拟合模型的残差值;
  • anova():生成一个拟合模型的方差分析表,或者比较两个或更多拟合模型的方差分析表;
  • vcov():列出模型参数的协方差矩阵;
  • AIC():输出赤池信息统计量;
  • plot():生成评价拟合模型的诊断图;
  • predict():用拟合模型对新的数据集预测响应变量值。

简单线性回归

基础安装中的数据集women提供了15个年龄在30~39岁间女性的身高和体重,以下代码希望获得一个等式,使得可以通过身高来预测体重,并分辨出过重或过瘦的个体。

得到的预测等式为:Weight=-87.52+3.45×Height

回归系数(3.45)显著不为0,表示身高每增高1英寸,体重将预期增加3.45磅。R平方项(0.991)表明模型可以解释体重99.1%的方差,这也是实际和预测值之间的相关系数。残差标准误(1.53)是模型用身高预测体重的平均误差,F统计量检验所有的预测变量预测响应变量是否都在某个几率水平之上,由于简单回归只有一个预测变量,此处F检验等同于身高回归系数的t检验。

多项式回归

使用以下代码来添加二次项,以提高以上例子的回归预测精度。

I(height^2)表示向预测等式添加一个身高的平方项,以下代码给出了详细过程。

新的预测等式为:Weight=261.88-7.35×Height+0.083×Height²

可以看到,模型的方差解释率已经增加到了99.9%,二次项显著性(t=13.89,p<0.001)表明包含二次项提高了模型的拟合度。

一般来说,n次多项式生成一个n-1个弯曲的曲线,若要拟合三次多项式,可用以下代码,不过使用三次以上的项几乎没有必要。

使用car包中的scatterplot()函数可以很容易且方便地绘制二元关系图,提供了散点图、线性拟合曲线和平滑拟合曲线,还在相应边界展示了各个变量的箱线图。

多元线性回归

当预测变量不止一个时,简单线性回归就变成了多元线性回归。以state.x77数据集为例,我们来研究一个州的犯罪率和其他因素的关系,包括人口、文盲、平均收入和结霜天数。

考虑到lm()函数的输入是数据框(state.x77数据集为矩阵),需要首先进行如下转换:

在进行多元回归分析中,第一步最好先检查一下变量间的相关性。cor()函数提供了两变量之间的相关系数,car包中的scatterplotMatrix()则会生成散点图矩阵。

scatterplotMatrix()函数默认在非对角线区域绘制变量间的散点图,并添加平滑和线性拟合曲线,在对角线区域绘制每个变量的密度图和轴须图。可以看出,谋杀率是双峰的曲线,每个预测变量都一定程度上出现了偏斜。谋杀率随着人口和文盲率的增加而增加,随着收入水平和结霜天数的增加而下降,同时,越冷的州文盲率越低、收入水平越高。

Scatter Plot Matrix

以下代码使用lm()函数拟合多元线性回归模型。

当预测变量不止一个时,回归系数的含义为:一个预测变量增加一个单位,其他预测变量保持不变时,因变量将要增加的数量。在本例中,文盲率的回归系数为4.14,表示控制人口、收入和温度不变时,文盲率上升1%,谋杀率将会上升4.14%,其系数在p<0.001的水平下显著不为0;相反,Frost的系数没有显著不为0(p=0.954),比表明当控制其他变量不变时,Frost和Murder不呈线性相关。总体来看,所有的预测变量解释了各州谋杀率57%的方差。

有交互项的多元线性回归

很多研究都会涉及交互项的预测变量,以mtcars数据框中的汽车数据为例,将汽车重量和马力作为预测变量,并包含交互项来拟合回归模型。

可以看到Pr(>|t|)栏中,马力和车重的交互项是显著的,说明响应变量与其中一个预测变量的关系依赖于另一个预测变量的水平,即每加仑汽油行驶英里数与汽车马力的关系依车重不同而不同。

8.3 回归诊断

虽然已经使用lm()函数拟合了OLS回归模型,通过summary()函数获取模型参数和相关统计量,但是为了确定模型是否合适、模型在多大程度上满足OLS模型统计假设,有必要对模型进行回归诊断。

使用confint()函数输出states多元回归问题的相关信息,结果表明,文盲率改变1%,谋杀率就在95%的置信区间[2.38, 5.90]中变化;Frost的置信区间包含0,说明当其他变量不变时,温度的改变和谋杀率无关。但是应当注意,这些结论都建立在数据满足统计假设的前提之上。

标准方法

最常见的检验方法就是对lm()函数返回的对象使用plot()函数,生成评价模型拟合情况的四幅图形。

women数据拟合分析

不妨来回顾一下OLS回归的统计假设:

  1. 正态性:当预测变量值固定时,因变量呈正态分布,其残差值应该是一个均值为0的正态分布。“正态Q-Q”(右上)是在正态分布对应的值下,标准化残差的概率图。若满足正态分布,则图上的点应落在呈45度角的直线上;
  2. 独立性:从以上四幅图中暂无法分辨出因变量值是否相互独立;
  3. 线性:若因变量和自变量线性相关,则残差值与预测值没有任何系统关系,即除了白噪声,模型应该包含数据中所有的系统方差。在“残差与拟合”(左上)中可以清楚地看到一个曲线关系,这意味着可能需要对回归模型加上一个二次项;
  4. 同方差性:若满足该性质,则在“位置尺度”(左下)中,水平线周围的点应随机分布(图中似乎满足此假设)。

最后一幅“残差与杠杆”(右下)指出了可能需要关注的单个观测点信息,包括离群点、高杠杆值点和强影响点。

  • 离群点:表明拟合回归模型对其预测效果不佳,产生了巨大的或正或负的残差;
  • 高杠杆值:表明其是一个异常的预测变量值组合,即在预测变量空间中,它是一个离群点,因变量值不参与计算一个观测点的杠杆值;
  • 强影响点:表明对模型参数的估计产生的影响过大,非常不成比例。

改进的方法

除了简单的使用plot()函数外,car包提供了大量拟合和评价回归模型的函数,gvlma包也提供了对所有线性模型假设进行检验的方法。

对于正态性,car包的qqPlot()函数提供了更为精确的正态假设检验方法。可以看出,除了Nevada,所有的点都离直线很近,并都落在置信区间内,表明正态性符合得很好。PS:点击某个点即可出现该点对应的标注。

car包的qqPlot()函数

以下自定义函数不需要加载car包,可以达到同样效果。

residplot函数分析结果

对于线性,可以使用成分残差图(也称偏残差图)来判断因变量和自变量之间是否呈线性关系,图形可用car包中的crPlots()函数绘制。

对于同方差性,可以使用car包的ncvTest()函数和spreadLevelPlot()函数进行检验。

另外,使用gvlma包的gvlma()函数,可以对线性模型假设进行综合验证,同时做偏斜度、峰度和异方差性的评价。

8.4 异常观测值

异常观测值包括之前提及的离群点、高杠杆值点和强影响点,这些数据点在一定程度上与其他观测点不同,可能对结果产生较大的负面影响。

离群点

在Q-Q图中,落在置信区间带外的点即可被认为是离群点。car包也提供了一种离群点的统计检验方法,outlierTest()函数,使用以下代码可以得出,Nevada被判定为离群点(p=0.048)。

高杠杆值点

高杠杆值观测点,即与其他预测变量有关的离群点,换言之,它们由许多异常的预测变量值组成,与响应变量没有关系。

高杠杆值的观测点可通过帽子统计量(hat statistic)判断,在以下代码所绘制的图形中,虚线以上(帽子均值2倍或3倍的位置)的点即为高杠杆值点,PS:点击某个点即可出现该点对应的标注。

强影响点

强影响点即对模型参数估计值影响有些比例失衡的点,当移除某个强影响点时,模型可能会随之发生巨大的改变。

使用Cook距离可以判断强影响点,以下代码判断出Alaska、Hawaii和Nevada是强影响点。

使用变量添加图(added variable plot)可以得到更直观的效果,car包中的avPlots()函数可提供变量添加图。图中的直线表示响应预测变量的实际回归系数,可以通过每个点对直线斜率的影响程度判断强影响点。

当然,还可以使用car包中的influencePlot()函数,将离群点、杠杆值和强影响点的信息整合到一幅图中。

结果显示,Nevada和Rhode Island是离群点,New York、California、Hawaii和Washington有高杠杆值,Nevada、Alaska和Hawaii为强影响点。

8.5 改进措施

对于回归诊断过程中发现的问题,有四种方法可以处理这些违背回归假设的问题:删除观测点、变量变换、添加或删除变量、使用其他回归方法。以下简要介绍删除观测点和变量变换。

删除观测点

删除离群点通常可以提高数据集对于正态假设的拟合度,而强影响点会干扰结果,通常也会被删除。删除最大的离群点或者强影响点后,模型需要重新拟合,若离群点或强影响点仍然存在,重复以上过程直至获得比较满意的拟合。

需要注意的是,删除观测点虽然可以剔除错误数据和无效数据,但是也有可能忽略掉异常数据,而这些异常点往往会包含一些非常有意义的信息。

变量变换

当模型不符合正态性、线性或同方差性假设时,一个或多个变量的变换通常可以改善或调整模型效果,即使用Y^λ替代Y。

car包中的powerTransform()函数通过λ的最大似然估计来正态化变量Y^λ,以下代码应用于数据states。

结果表明,可以用Murder^0.6来正态化Murder,考虑到0.6很接近0.5,故也可以尝试用平方根变换来提高模型正态性的符合程度。

当违反了线性假设时,对预测变量进行变换常常会比较有用,car包中的boxTidwell()函数即通过获得预测变量幂数的最大似然估计来改善线性关系。

8.6 选择“最佳”的回归模型

尝试获取一个回归方程时,便面对着从众多可能的模型中做选择的问题。一般来说,模型选择的目标有两点:准确和精简。

模型比较

可以使用基础安装中的anova()函数比较两个嵌套模型的拟合优度,即一个模型的项完全包含于另一个模型中。由于检验不显著(p=0.994),因此可以认为不需要将Income和Frost添加到线性模型中。

AIC()考虑了模型的统计拟合度以及用来拟合的参数数目,AIC值越小的模型越要优先选择,说明模型用较少的参数获得了足够的拟合度。不同于anova(),AIC()无需嵌套模型。

变量选择

以上方法仅适用于比较两个已知模型,如果涉及到的自变量较多,需要从中筛选出合适的预测变量,常用的方法有两种:逐步回归法和全子集回归。

逐步回归法包括以下三类:

  1. 向前逐步回归:从不包含任何预测变量开始,每次添加一个预测变量,直到添加变量不会使模型有所改进;
  2. 向后逐步回归:从包含所有预测变量开始,每次删除一个变量,直到会降低模型质量为止;
  3. 逐步回归:每次进入一个变量,重新评价后删除没有共享的变量,预测变量可能会被添加或删除好几次,直到获得最优。

MASS包中的stepAIC()函数给予AIC准则实现了逐步回归模型(向前、向后和双向)。

需要注意的是,逐步回归法虽然能找到一个局部的好模型,但是不能保证其就是最优的。

全子集回归克服了逐步回归的缺点,即所有可能的模型都会被检验,评判准则可以是R平方、调整R平方或者Mallows Cp统计量。

8.7 深层次分析

最后,介绍评价模型泛化能力和变量相对重要性的方法。

交叉验证

交叉验证是指将一定比例的数据挑选出来作为训练样本,另外的样本作为保留样本,先在训练样本上获取回归方程,然后在保留样本上做预测。

所谓k重交叉验证,即样本被分为k个子样本,然后轮流将k-1个子样本组合作为训练集,另外一个子样本作为保留集,得到k个预测方程并记录k个保留样本的预测表现结果,最后求其平均值。当k等于观测总数目时,k重交叉验证又称作刀切法。

相对重要性

为了获取各个预测变量重要性的相对排序,最简单的方法是比较其回归系数。注意在进行回归分析前,需要用scale()函数将数据标准化为均值为0、标准差为1的数据集,接下来用R回归即可获得标准化的回归系数。

22 2015-04

1条评论

  1. 匿名说道:

    真棒!!!

发表评论