8.3回归诊断
> fit<-lm(weight~height,data=women)
> par(mfrow=c(2,2))
> plot(fit)
为理解这些图形,我们来回顾一下oLs回归的统计假设。
口正态性当预测变量值固定时,因变量成正态分布,则残差值也应该是一个均值为0的正态分布。正态Q-Q图(Normal Q-Q,右上)是在正态分布对应的值下,标准化残差的概率图。若满足正态假设,那么图上的点应该落在呈45度角的直线上;若不是如此,那么就违反了正态性的假设。
口独立性你无法从这些图中分辨出因变量值是否相互独立,只能从收集的数据中来验证。上面的例子中,没有任何先验的理由去相信一位女性的体重会影响另外一位女性的体重。假若你发现数据是从一个家庭抽样得来的,那么可能必须要调整模型独立性的假设。
口线性若因变量与自变量线性相关,那么残差值与预测(拟合)值就没有任何系统关联。换句话说,除了自噪声,模型应该包含数据中所有的系统方差。在“残差图与拟合图”( Residuals vs Fitted,左上)中可以清楚的看到一个曲线关系,这暗示着你可能需要对回归模型加上一个二次项。
口同方差性若满足不变方差假设,那么在位置尺度图(Scale-Location Graph,左下)中,水平线周围的点应该随机分布。该图似乎满足此假设。最后一幅“残差与杠杆图”(Residuals vs Leverage,右下)提供了你可能关注的单个观测点的信息。从图形可以鉴别出离群点、高杠杆值点和强影响点。
8.3.2改进的方法
qqPlot() 分位数比较图
durbinWatsonTest()对误差自相关性做Durbin-Watson检验
crPlots()成分与残差图
ncvTest()对非恒定的误差方差做得分检验
spreadLevelPlot()分散水平检验
outlierTest()Bonferroni离群点检验
avPlots()添加的变量图形
inluencePlot()回归影响图
scatterplot()增强的散点图
scatterplotMatrix()增强的散点图矩阵
vif()方差膨胀因子
1.正态性
与基础包中的plot ( )函数相比,qqPlot()函数提供了更为精确的正态假设检验方法,它画出了在n-p-1个*度的t分布下的学生化残差(( studentized residual,也称学生化删除残差或折叠化残差)图形,其中n是样本大小,p是回归参数的数目(包括截距项)。
Eg:
> library(car)
> states=data.frame(state.region,state.x77)
> fit<-lm(Murder~Population+Illiteracy+Income+Frost,data=states)
> qqPlot(fit,labels=row.names(states),id.method="identify",simulate=TRUE,main="Q-QPlot")
绘制学生化残差图的函数
residplot<-function(fit,nbreaks=10){
z<-rstudent(fit)
hist(z,breaks=nbreaks,freq=FALSE,
xlab="Studentized Residual",
main="Distribution of Errors")
rug(jitter(z),col="brown")
curve(dnorm(x,mean=mean(z),sd=sd(z)),
add=TRUE,col="blue",lwd=2)
legend("topright",legend=c("NormalCurve","Kernel Density Curve"),
lty=1:2,col=c("blue","red"),cex=.7)
}
residplot(fit)
2.误差的独立性
car包提供了一个可做Durbin-Watson检验的函数,能够检测误差的序列相关性。
> durbinWatsonTest(fit)
lagAutocorrelation D-W Statistic p-value
1 -0.2006929 2.317691 0.284
Alternative hypothesis: rho != 0
3. 线性
通过成分残差图(component plus residual plot)也称偏残差图(partialresidual plot),你可以看看因变量与自变量之间是否呈非线性关系,也可以看看是否有不同于已设定线性模型的系统偏差,图形可用car包中的crPlots()函数绘制。
4. 同方差性
ncvTest()函数生成一个计分检验,零假设为误差方差不变,备择假设为误差方差随着拟合值水平的变化而变化。若检验显著,则说明存在异方差性(误差方差不恒定)。spreadLevelPlot()函数创建一个添加了最佳拟合曲线的散点图,展示标准化残差绝对值与拟合值的关系
检验同方差性
> library(car)
> ncvTest(fit)
Non-constant Variance Score Test
Variance formula: ~ fitted.values
Chisquare = 1.746514 Df = 1 p = 0.1863156
> spreadLevelPlot(fit)
Suggested power transformation: 1.209626
8.3.3 线性模型假设的综合验证
gvlma包中的gvlma()函数
> library(gvlma)
> gvmodel<-gvlma(fit)
> summary(gvmodel)
Call:
lm(formula = Murder ~Population + Illiteracy + Income + Frost,
data = states)
Residuals:
Min 1Q Median 3Q Max
-4.7960 -1.6495-0.0811 1.4815 7.6210
Coefficients:
Estimate Std. Error t valuePr(>|t|)
(Intercept)1.235e+00 3.866e+00 0.319 0.7510
Population 2.237e-04 9.052e-05 2.471 0.0173 *
Illiteracy 4.143e+00 8.744e-01 4.738 2.19e-05 ***
Income 6.442e-05 6.837e-04 0.094 0.9253
Frost 5.813e-04 1.005e-02 0.058 0.9541
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’1
Residual standarderror: 2.535 on 45 degrees of freedom
MultipleR-squared: 0.567, Adjusted R-squared: 0.5285
F-statistic: 14.73 on 4and 45 DF, p-value: 9.133e-08
ASSESSMENT OF THELINEAR MODEL ASSUMPTIONS
USING THE GLOBAL TESTON 4 DEGREES-OF-FREEDOM:
Level of Significance= 0.05
Call:
gvlma(x = fit)
Value p-value Decision
Global Stat 2.7728 0.5965 Assumptions acceptable.
Skewness 1.5374 0.2150 Assumptions acceptable.
Kurtosis 0.6376 0.4246 Assumptions acceptable.
Link Function 0.1154 0.7341 Assumptions acceptable.
Heteroscedasticity0.4824 0.4873 Assumptions acceptable.
8.3.4 多重共线性
多重共线性可用统计量VIF(Variance Inflation Factor,方差膨胀因子)进行检测。VIF的平
方根表示变量回归参数的置信区间能膨胀为与模型无关的预测变量的程度(因此而得名)。car
包中的vif()函数提供VIF值。一般原则下,vif >2就表明存在多重共线性问题
检测多重共线性:
> library(car)
> vif(fit)
PopulationIlliteracy Income Frost
1.245282 2.165848 1.345822 2.082547
> sqrt(vif(fit))>2
PopulationIlliteracy Income Frost
FALSE FALSE FALSE FALSE