在我们以前做线性回归的时候,总是公式上手直接套,并没有考虑到能够利用线性回归先要满足零均值、同方差、无自相关以及正态分布等假定,假设数据不符合同方差假定,那么残差的方差会被严重低估,从而使得伪错误的几率上升,整个模型的拟合精度会大幅下降。
那么我们如何去检测数据是否存在异方差呢?
在R语言里面,我们可以通过统计检验、图形法来找出答案。
例子数据是R内置的cars数据集,先利用lm函数来建立模型:
lmMod <- lm(dist ~ speed, data=cars) # 初始模型
#导入包
library(e1071)
library(caret) # 含box-cox 转换的包
library(lmtest) # 检验异方差的bp检验包
library(gvlma) # 线性模型的全局验证
bptest(lmMod) # Breusch-Pagan 检验
studentized Breusch-Pagan test
data: lmMod
BP = 3.2149, df = 1, p-value = 0.07297
这里结果中的p值大于0.05表示零假设(同方差检验)可以被拒绝,因此存在异方差,再通过线性模型假设(gvlma)的全局验证来确认这个结论。
gvlma(lmMod) # 验证线性回归的假设是否成立。
# Call:
gvlma(x = lmMod)
Value p-value Decision
Global Stat 15.801 0.003298 Assumptions NOT satisfied!
Skewness 6.528 0.010621 Assumptions NOT satisfied!
Kurtosis 1.661 0.197449 Assumptions acceptable.
Link Function 2.329 0.126998 Assumptions acceptable.
Heteroscedasticity 5.283 0.021530 Assumptions NOT satisfied!#验证上述结论
既然同方差不满足,来看看在图形方面数据是怎么表现的。
par(mfrow=c(2,2)) # 将图片切为四份
plot(lmMod)
图形如下:
左上角残差对拟合值和左下角标准化残差对拟合值的两幅图值得我们关注,如果同方差假设成立,我们会看到一个完全随机的分布,红线相对平坦,然而在这里明显是有一定“规律”的,错误(读取残差)不是完全随机的,并且可能包含一些信息,而这些信息又可以用于预测因变量。
既然知道假设不成立,我们可以通过以下两个办法来修正问题:
1、 通过添加创建的错误项(lmMod $ residuals),再重新构建模型
2、 通过转换因变量
Box-Cox变换
对于这个例子,我们将继续使用第二种方法来修改模型并预测'dist'。我们默认使用的转换是box-cox转换。
一旦找到最佳lambda值,我们就将使用它来变换为因变量,并使用新变换的变量构建模型。
distBCMod <- BoxCoxTrans(cars$dist)
print(distBCMod)
Box-Cox Transformation
50 data points used to estimate Lambda
Input data summary:
Min. 1st Qu. Median Mean 3rd Qu. Max.
2.00 26.00 36.00 42.98 56.00 120.00
Largest/Smallest: 60
Sample Skewness: 0.759
Estimated Lambda: 0.5
结果的最后一行显示,lambda的估计值为0.5,再进行变换var =(y ^ lambda - 1)/ lambda。 我们不必手动执行此操作,可以通过使用插入包中的“BoxCoxTrans”函数来完成。现在,我们将这个新的转换变量附加到汽车数据并重新构建模型。
cars <- cbind(cars, dist_new=predict(distBCMod, cars$dist)) # 增添转换值到数据列
head(cars) # 看下头6行数据
speed dist dist_new
1 4 2 0.8284271
2 4 10 4.3245553
3 7 4 2.0000000
4 7 22 7.3808315
5 8 16 6.0000000
6 9 10 4.3245553
重新构建模型,并检测同方差检验
lmMod_bc <- lm(dist_new ~ speed, data=cars)
bptest(lmMod_bc)
studentized Breusch-Pagan test
data: lmMod_bc
BP = 0.011192, df = 1, p-value = 0.9157
P值为0.91,我们不能拒绝零假设同方差假设。
gvlma(lmMod_bc) # 检查假设
Call:
lm(formula = dist_new ~ speed, data = cars)
Coefficients:
(Intercept) speed
0.5541 0.6448
ASSESSMENT OF THE LINEAR MODEL ASSUMPTIONS
USING THE GLOBAL TEST ON 4 DEGREES-OF-FREEDOM:
Level of Significance = 0.05
Call:
gvlma(x = lmMod_bc)
Value p-value Decision
Global Stat 3.03592 0.55183 Assumptions acceptable.
Skewness 2.81284 0.09351 Assumptions acceptable.
Kurtosis 0.04915 0.82455 Assumptions acceptable.
Link Function 0.09269 0.76078 Assumptions acceptable.
Heteroscedasticity 0.08124 0.77563 Assumptions acceptable.#这下好了
从上面的输出以及bptest(),我们都可以验证异方差性问题是否得到有效解决,再来看看图像。
plot(lmMod_bc)
图形如下:
在左上图中,残差变得更随机了,具有更直的平滑线,这意味着残差的方差不随着拟合值的增加而改变,这意味着去除了异方差性。顺利解决问题了。
鸡腿吃不腻
1楼 - 6 年,3 月 之前
这难道是伟大的原创(*´▽`)ノノ