用R進行多元線性迴歸分析建模
阿新 • • 發佈:2019-02-09
概念:多元迴歸分析預測法,是指通過對兩個或兩個以上的自變數與一個因變數的相關分析,建立預測模型進行預測的方法。當自變數與因變數之間存在線性關係時,稱為多元線性迴歸分析。
下面我就舉幾個例子來說明一下
例一:謀殺率與哪些因素有關
變數選擇
我們可以明顯的看出謀殺率與人口,文盲率相關性較大states<-as.data.frame(state.x77[,c('Murder','Population','Illiteracy','Income','Frost')]) cor(states)#檢視變數相關係數 Murder Population Illiteracy Income Frost Murder 1.0000000 0.3436428 0.7029752 -0.2300776 -0.5388834 Population 0.3436428 1.0000000 0.1076224 0.2082276 -0.3321525 Illiteracy 0.7029752 0.1076224 1.0000000 -0.4370752 -0.6719470 Income -0.2300776 0.2082276 -0.4370752 1.0000000 0.2262822 Frost -0.5388834 -0.3321525 -0.6719470 0.2262822 1.0000000
將它們的關係視覺化
library(car)
scatterplotMatrix(states,spread=FALSE)
還可以這麼看
fit<-lm(Murder~Population+Illiteracy+Income+Frost,data = states) summary(fit) 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
還可以這麼看
#install.packages('leaps')
library(leaps)
leaps<-regsubsets(Murder~Population+Illiteracy+Income+Frost,data = states,nbest = 4)
plot(leaps,scale = 'adjr2')
最大值0.55是隻包含人口,文盲率這兩個變數和截距的。
還可以這樣,比較標準迴歸係數的大小
zstates<-as.data.frame(scale(states))#scale()標準化 zfit<-lm(Murder~Population+Illiteracy+Income+Frost,data = zstates) coef(zfit) (Intercept) Population Illiteracy Income Frost -2.054026e-16 2.705095e-01 6.840496e-01 1.072372e-02 8.185407e-03
通過這幾種方法,我們都可以明顯的看出謀殺率與人口,文盲率相關性較大,與其它因素相關性較小。
迴歸診斷
> 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
標記異常值qqPlot(fit,labels = row.names(states),id.method = 'identify',simulate = T)
圖如下,點一下異常值然後點finish就可以了檢視它的實際值11.5與擬合值3.878958,這條資料顯然是異常的,可以拋棄
> states['Nevada',]
Murder Population Illiteracy Income Frost
Nevada 11.5 590 0.5 5149 188
> fitted(fit)['Nevada']
Nevada
3.878958
> outlierTest(fit)#或直接這麼檢測離群點
rstudent unadjusted p-value Bonferonni p
Nevada 3.542929 0.00095088 0.047544
car包有多個函式,可以判斷誤差的獨立性,線性,同方差性library(car)
durbinWatsonTest(fit)
crPlots(fit)
ncvTest(fit)
spreadLevelPlot(fit)
綜合檢驗#install.packages('gvlma')
library(gvlma)
gvmodel<-gvlma(fit);summary(gvmodel)
檢驗多重共線性
根號下vif>2則表明有多重共線性
> sqrt(vif(fit))
Population Illiteracy Income Frost
1.115922 1.471682 1.160096 1.443103
都小於2所以不存在多重共線性
例二:女性身高與體重的關係
attach(women)
plot(height,weight)
通過圖我們可以發現,用曲線擬合要比直線效果更好
那就試試唄
fit<-lm(weight~height+I(height^2))#含平方項
summary(fit)
Call:
lm(formula = weight ~ height + I(height^2))
Residuals:
Min 1Q Median 3Q Max
-0.50941 -0.29611 -0.00941 0.28615 0.59706
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 261.87818 25.19677 10.393 2.36e-07 ***
height -7.34832 0.77769 -9.449 6.58e-07 ***
I(height^2) 0.08306 0.00598 13.891 9.32e-09 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.3841 on 12 degrees of freedom
Multiple R-squared: 0.9995, Adjusted R-squared: 0.9994
F-statistic: 1.139e+04 on 2 and 12 DF, p-value: < 2.2e-16
效果是很不錯的,可以得出模型為
把擬合曲線加上看看
lines(height,fitted(fit))
非常不錯吧
還可以用car包的scatterplot()函式
library(car)
scatterplot(weight~height,spread=FALSE,pch=19)#19實心圓,spread=FALSE刪除了殘差正負均方根在平滑曲線上
展開的非對稱資訊,聽著就不像人話,你可以改成TRUE看看到底是什麼,我反正不明白。
例三:含互動項
<strong>attach(mtcars)
fit<-lm(mpg~hp+wt+hp:wt)
summary(fit)
Call:
lm(formula = mpg ~ hp + wt + hp:wt)
Residuals:
Min 1Q Median 3Q Max
-3.0632 -1.6491 -0.7362 1.4211 4.5513
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 49.80842 3.60516 13.816 5.01e-14 ***
hp -0.12010 0.02470 -4.863 4.04e-05 ***
wt -8.21662 1.26971 -6.471 5.20e-07 ***
hp:wt 0.02785 0.00742 3.753 0.000811 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 2.153 on 28 degrees of freedom
Multiple R-squared: 0.8848, Adjusted R-squared: 0.8724
F-statistic: 71.66 on 3 and 28 DF, p-value: 2.981e-13</strong>
其中的hp:wt就是互動項,表示我們假設hp馬力與wt重量有相關關係,通過全部的三個星可以看出響應/因變數mpg(每加侖英里)與預測/自變數都相關,也就是說mpg(每加侖英里)與汽車馬力/重量都相關,且mpg與馬力的關係會根據車重的不同而不同。