> states<- as.data.frame(state.x77[,c("Murder", "Population","Illiteracy","Income","Frost")])
> fit <- lm(Murder ~ Population + Illiteracy + Income + Frost, data=states)
> confint(fit)
? ? ? ? ? ? ? ? ? ? 2.5 %? ? ? 97.5 %
(Intercept) -6.552191e+00 9.0213182149
Population? 4.136397e-05 0.0004059867
Illiteracy? 2.381799e+00 5.9038743192
Income? ? ? -1.312611e-03 0.0014414600
Frost? ? ? -1.966781e-02 0.0208304170
> #結(jié)果表明，文盲率改變1%時，謀殺率就在95%的置信區(qū)間[2.38，5.90]中變化。另外，因為Frost的置信區(qū)間包含0，所以可以得出結(jié)論：當其他變量不變時，溫度的改變與謀殺率無關(guān)。
標準方法
>fit <- lm(weight~ height, data=women)
> par(mfrow=c(2,2))
> plot(fit) #結(jié)果如下圖
正態(tài)性 當預(yù)測變量值固定時，因變量成正態(tài)分布，則殘差值也應(yīng)該是一個均值為0的正態(tài)分布。正態(tài)“Q-Q，右上”是正態(tài)分布對應(yīng)的值下，標準化殘差的概率圖。拖滿足正態(tài)分布，那么圖上的點應(yīng)該落在呈45°的直線上；否則，就違反了正態(tài)性的假設(shè)；
線性? 若因變量和自變量線性相關(guān)，模型應(yīng)該包含數(shù)據(jù)中所有的系統(tǒng)方差。在“殘差圖與擬合圖，左上”中可以清楚地看到一個曲線關(guān)系，暗示著可能需要對回歸模型加上一個二次項。
同方差性? 若滿足不變方差假設(shè)，那么“位置尺度圖，左下”中，水平線周圍的點應(yīng)該隨機分布。
殘差與杠桿圖（右下），提供了可能關(guān)注的單個觀測點的信息，從圖形中可以鑒別出離群點，高桿杠點和強影響點。
第二次修正
> fit2 <- lm(weight ~ height + I(height^2), data=women)
Warning message:
顯示串列沒有完全被刷新
> par(mfrow=c(2,2))
> plot(fit2)
> #這兒的顯示基本還算理想，但是觀測點13和觀測點15對整個回歸的影響還是比較大，刪除之后再看結(jié)果
> newfit<- lm(weight ~ height +I(height^2),data=women[-c(13,15),])
> par(mfrow=c(2,2))
Warning message:
顯示串列沒有完全被刷新
> plot(newfit)
繼續(xù)來看運用函數(shù)進行修正
car包提供了大量函數(shù)，大大增強了擬合和評價回歸模型的能力
qqPlot()? ? ? ? ? ? ? ? ? ? ?分位數(shù)比較圖
durbinWatsonTest()? 對誤差自相關(guān)性做durbin-Watson檢驗
ncvTest()? ? ? ? ? ? ? ? ? ? ?對非恒定的誤差方法做的分檢驗
speardLevelPlot()? ? ? ?分散水平檢驗
outlierTest()? ? ? ? ? ? ? ? ? bonferronni離群點檢驗
avluencePlot()? ? ? ? ? ? ? 回歸影響圖
scatterplot()? ? ? ? ? ? ? ? ? 增強的散點圖
scatterplotMatrix()? ? ? ? 增強的散點圖矩陣
vif()? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?方差膨脹因子
舉例來說
> states <- as.data.frame(state.x77[,c("Murder","Population","Illiteracy","Income","Frost")])
> fit <- lm(Murder ~ Population + Illiteracy + Income + Frost, data=states)
> par(mfrow=c(2,2))
> plot(fit)
很明顯，結(jié)果中“Nevada”是對結(jié)果的影響就比較大，屬于噪點
可視化誤差：
> states["Nevada",]
? ? ? Murder Population Illiteracy Income Frost
Nevada? 11.5? ? ? ? 590? ? ? ? 0.5? 5149? 188
> fitted(fit)["Nevada"] #顯示擬合值
? Nevada
3.878958
> residuals(fit)["Nevada"] #殘差擬合值
? Nevada
7.621042
> rstudent(fit)["Nevada"]
? Nevada
3.542929
> #可見，Nevada有一個很大的正殘差值（真實值-預(yù)測值），表明該模型低估了該州的謀殺率。Nevada的謀殺率是11.5%，而模型預(yù)測的謀殺率是3.9%。
繪制學生化殘差圖函數(shù)
> residplot<- function(fit, nbreaks=10){
+ z<- rstudent(fit)
+ hist(z, breaks=nbreaks,fre=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)
+ lines(density(z)$x, density(z)$y,col="red", lwd=2, Ity=2)
+ legend("topright", legend=c("Normal Curve","Kernel Density Curve"), Ity=1:2, col=c("blue","red"), cex=.7)
+ }
> residplot(fit)
結(jié)果如下：
線性模型假設(shè)的綜合驗證
gvlma( )函數(shù)由pena和slate編寫，能對線性模型假設(shè)進行綜合驗證，同時還能做偏斜度、峰度和差異方差性的評價。
> 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 value Pr(>|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 standard error: 2.535 on 45 degrees of freedom
Multiple R-squared:? 0.567,? ? Adjusted R-squared:? 0.5285
F-statistic: 14.73 on 4 and 45 DF,? p-value: 9.133e-08
ASSESSMENT OF THE LINEAR MODEL ASSUMPTIONS
USING THE GLOBAL TEST ON 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.
Heteroscedasticity 0.4824? 0.4873 Assumptions acceptable.
> #該結(jié)果中的P=0.597，高于顯著性水平，因此接受原假設(shè)。
R語言回歸診斷的基本知識到這就結(jié)束了，咱們下期再見！O(∩_∩)O哈哈~
