柚子快報激活碼778899分享：R語言統(tǒng)計分析——回歸診斷2

http://yzkb.51969.com/

參考資料：R語言實戰(zhàn)【第2版】

????????R語言的car包提供的大量函數(shù)，大大增強了擬合和評價回歸模型的能力。如下：

函數(shù)目的qqPlot()分位數(shù)比較圖durbinWatsonTest()對誤差自相關(guān)性做Durbin-Watson檢驗crPlots()成分與殘差圖ncvTest()對非恒定的誤差方差做得分檢驗spreadLevelPlot()分散水平檢驗outlierTest() Bonferroni離群點檢驗 avPlots()添加的變量圖形influencePlot()回歸影響圖scatterplot()增強的散點圖scatterplotMatrix()增強的散點圖矩陣vif()方差膨脹因子

? ? ? ? 另外，gvlma包提供了對所有線性模型假設(shè)進行檢驗的方法。

1、正態(tài)性

? ? ? ? 與基礎(chǔ)包中的plot()函數(shù)相比，qqPlot()函數(shù)提供了更為精確的正態(tài)假設(shè)檢驗方法，它畫出了在n-p-1個自由度的圖分布下的學(xué)生化殘差（studentized residual，也稱學(xué)生化刪除殘差或折疊化殘差）圖形，其中n是樣本大小，p是回歸參數(shù)的數(shù)目（包括截距項）。如下：

# 加載car包

library(car)

# 獲取數(shù)據(jù)

states<-as.data.frame(state.x77[,c("Murder","Population",

"Illiteracy","Income","Frost")])

# 擬合多元線性模型

fit<-lm(Murder~Population+Illiteracy+Income+Frost,data=states)

# 繪制Q-Q圖

qqPlot(fit,

simulate=TRUE,

main="Q-Q Plot")

? ? ? ? ?當simulate=TRUE時，將會使用參數(shù)自助法生成95%的置信區(qū)間。

? ? ? ? 從上圖中可以看出，除了Nevada，所有的點都離直線很近，并都落在置信區(qū)間內(nèi)，這表明正態(tài)性假設(shè)符合得很好。作為一個異常值，我們也應(yīng)該關(guān)注Nevada，它有一個很大的正殘差值（真實值-預(yù)測值），表明模型低估了該地的謀殺率。

# 查看Nevada數(shù)據(jù)

states["Nevada",]

# 查看Nevada的擬合數(shù)據(jù)（預(yù)測值）

fitted(fit)["Nevada"]

# 查看Nevada的殘差

residuals(fit)["Nevada"]

# 查看Nevada的學(xué)生化殘差

rstudent(fit)["Nevada"]

注：學(xué)生化殘差是殘差除以它的標準差后得到的數(shù)值，用以直觀地判斷誤差項服從正態(tài)分布這一假定是否成立，若假定成立，學(xué)生化殘差的分布也應(yīng)服從正態(tài)分布?！緦W(xué)生化殘差_百度百科】

? ? ? ? 可視化誤差還有其他方法，比如使用residplot()函數(shù)生成學(xué)生化殘差柱狀圖（即直方圖），并添加正態(tài)曲線、核密度曲線和軸須圖。

resiplot<-function(fit,nbreaks=10){

z<-rstudent(fit)

hist(z,breaks=nbreaks,freq=FALSE,

xlab="Studentized Residual",

main="Distribution of Error")

rug(jitter(z),col="brown")

curve(dnorm(x,mean=mean(z),sd=sd(z)),

add=TRUE,col="blue",lwd=2)

lines(density(z)$x,density(z)$y,

col="red",lwd=2,lty=2)

legend("topright",

legend=c("Normal Curve","Kernel Density Curve"),

lty=1:2,col=c("blue","red"),cex=0.7)

}

resiplot(fit)

? ? ? ? 如上圖所示，除了一個明顯的離群點，誤差基本上都很好的服從了正態(tài)分布。

2、誤差的獨立性

? ? ? ? 判斷因變量值（或殘差）是否相互獨立，最好的方法是依據(jù)收集數(shù)據(jù)方式的先驗知識。例如，時間序列數(shù)據(jù)通常呈現(xiàn)自相關(guān)——相隔時間越近的觀測相關(guān)性大于相隔越遠的觀測。car包提供了一個可做Durbin-Watson檢驗的函數(shù)，能夠檢測誤差的序列相關(guān)性。繼續(xù)使用上面的案例進行Durbin-Watson檢驗：

durbinWatsonTest(fit)

? ? ? ? 結(jié)果中，p值不顯著（p=0.242）說明無自相關(guān)性，誤差之間相互獨立。滯后項（lag=1）表明數(shù)據(jù)集中每個數(shù)據(jù)都是與其中后一個數(shù)據(jù)進行比較的。該檢驗適用于時間獨立的數(shù)據(jù)，對于非聚集型的數(shù)據(jù)并不使用。注意，durbinWatsonTest()函數(shù)使用的自助法來導(dǎo)出p值。如果添加了選項simulate=TRUE，則每次運行測試時獲得的結(jié)果都將略有不同。

3、線性

? ? ? ? 通過成分殘差圖（component plus residual plot）也稱偏殘差圖（partial residual plot），我們可以查看因變量與自變量之間是否呈非線性關(guān)系，也可以查看是否有不同于已設(shè)定線性模型的系統(tǒng)偏差，圖形可用car包中的crPlots()函數(shù)繪制。繼續(xù)使用上面的案例數(shù)據(jù)，如下：

# 導(dǎo)入car包

library(car)

# 繪制成分殘差圖

crPlots(fit)

? ? ? ? 如果成分殘差圖中的圖形存在非線性，則說明我們對自變量的函數(shù)形式建模不夠充分。那么就需要添加一些曲線成分，如多項式，或?qū)σ粋€或多個變量進行變換（如log(x)代替x），或用其他回歸變體形式而不是線性回歸。

? ? ? ? 從上圖中可以看出，4個自變量的成分殘差圖可以認為是線性的，尤其是前兩個自變量。

4、同方差性

? ? ? ? car包提供了兩個有用的函數(shù)，可以判斷誤差方差是否恒定。

????????ncvTest()函數(shù)生成一個計分檢驗，零假設(shè)為誤差方差不變，備擇假設(shè)為誤差方差隨著擬合值水平的變化而變化。若檢驗顯著，則說明存在異方差性（誤差方差不恒定）。

? ? ? ? spreadLevelPlot()函數(shù)創(chuàng)建一個添加了最佳擬合曲線的散點圖，展示標準化殘差決對值與擬合值得關(guān)系。

# 導(dǎo)入car包

library(car)

ncvTest(fit)

spreadLevelPlot(fit)

? ? ? ? 根據(jù)ncvTest()的計分檢驗結(jié)果顯示不顯著（p=0.186），說明滿足方差不變假設(shè)。?

? ? ? ? spreadLevelPlot()函數(shù)作圖顯示，數(shù)據(jù)點在水平的最佳擬合曲線（虛線）周圍呈水平隨機分布。如果違反了該假設(shè)，我們將看到一個非水平的曲線。而運行代碼結(jié)果顯示“Suggested power transformation: 1.209626”，建議進行1.2次冪的冪次轉(zhuǎn)換。因接近1，所以不需要進行變換。

柚子快報激活碼778899分享：R語言統(tǒng)計分析——回歸診斷2

http://yzkb.51969.com/

精彩文章