UA MATH571A R語(yǔ)言回歸分析實(shí)踐 多元回歸2 醫(yī)療費(fèi)用的決定

• 系數(shù)的推斷與模型預(yù)測(cè)
• 模型診斷

這一講展示一下一元回歸中的模型診斷的手段怎么用在多元回歸中，同時(shí)介紹一下多元回歸做推斷和預(yù)測(cè)的方法。

系數(shù)的推斷與模型預(yù)測(cè)

先來(lái)看看系數(shù)的99%置信區(qū)間，顯著（顯著性水平1%）的系數(shù)的99%置信區(qū)間應(yīng)該不會(huì)包含0：

> confint(mreg01.lm,level=0.99)0.5 % 99.5 % (Intercept) -16116.04048 -10391.0546 Age 226.77096 288.1101 Sex -988.77965 728.8666 BMI 259.16501 401.8577 Chil 121.91230 832.2460 Smoke 22758.69080 24884.2267 Region -99.23445 956.4028

用confint函數(shù)可以計(jì)算置信區(qū)間，第一個(gè)輸入是模型對(duì)象，第二個(gè)輸入是置信水平。可以發(fā)現(xiàn)性別和地區(qū)的99%置信水平是包含0的，這與上一講發(fā)現(xiàn)這兩個(gè)解釋變量的系數(shù)在1%的顯著性水平下不顯著是一致的。需要注意的是這個(gè)置信區(qū)間是單個(gè)系數(shù)的置信區(qū)間，也就是說(shuō)當(dāng)我們關(guān)注Age的時(shí)候，是只考慮Age，沒有考慮其他系數(shù)的。因此更合理的做法是計(jì)算這些系數(shù)的聯(lián)合置信區(qū)間，可以采用Bonferroni調(diào)整或者WHS置信帶方法。比如計(jì)算Bonferroni調(diào)整系數(shù)

> n <- 7 # intercept and 6 coefficients > B <- qt(1-0.01/(2*n),mreg01.lm$df) > B [1] 3.195517

然后根據(jù)
$$[{\hat{\beta }}_i?Bse({\hat{\beta }}_i),{\hat{\beta }}_i+Bse({\hat{\beta }}_i)]$$
計(jì)算系數(shù)${\beta }_i$的Bonferroni調(diào)整的置信區(qū)間。或者計(jì)算WHS調(diào)整系數(shù)，

> W <- sqrt(n*qf(1-0.01,n,mreg01.lm$df)) > W [1] 4.309184

然后根據(jù)
$$[{\hat{\beta }}_i?Wse({\hat{\beta }}_i),{\hat{\beta }}_i+Wse({\hat{\beta }}_i)]$$
計(jì)算系數(shù)${\beta }_i$的WHS調(diào)整的置信區(qū)間，實(shí)際應(yīng)用的時(shí)候取$B$和$W$中的較小者。

接下來(lái)介紹一下多元線性回歸中怎么應(yīng)用predict函數(shù)。先討論擬合的用法，以第一個(gè)樣本為例：

predict(mreg01.lm,newdata = data.frame(Age = 19,Sex = 0,BMI = 27.9,Chil = 0,Smoke = 1,Region = 1),interval = "conf",level = 0.95) ## output ## fit lwr upr ## 1 25109.13 24004.59 26213.67

輸入第一項(xiàng)是模型對(duì)象，這一點(diǎn)與一元回歸相同，輸入第二項(xiàng)是一個(gè)數(shù)據(jù)框，數(shù)據(jù)框要包含所有的解釋變量，需要輸入要擬合的樣本每一個(gè)解釋變量對(duì)應(yīng)的值，然后data.frame會(huì)將其自動(dòng)合成一個(gè)數(shù)據(jù)框。如果要擬合多個(gè)樣本的值，以前三個(gè)樣本為例，可以用c()將每個(gè)解釋變量的值一起輸入到數(shù)據(jù)框中

predict(mreg01.lm,newdata = data.frame(Age = c(19,18,28),Sex = c(0,1,1),BMI = c(27.9,33.77,33),Chil = c(0,1,3),Smoke = c(1,0,0),Region = c(1,2,2)),interval = "conf",level = 0.95) ## output ## fit lwr upr ## 1 25109.132 24004.594 26213.670 ## 2 3746.041 3016.915 4475.167 ## 3 7020.111 6241.170 7799.053

如果不想手動(dòng)輸入數(shù)字可以用之前定義好的變量：

predict(mreg01.lm,newdata = data.frame(Age = Age[1:3],Sex = Sex[1:3],BMI = BMI[1:3],Chil = Chil[1:3],Smoke = Smoke[1:3],Region = Region[1:3]),interval = "conf",level = 0.95)

它會(huì)得到和手動(dòng)輸入數(shù)字一樣的結(jié)果，甚至我們還可以用讀入的原始數(shù)據(jù)來(lái)構(gòu)造數(shù)據(jù)框：

predict(mreg01.lm,newdata = data.frame(Age = Data[1:3,1],Sex = Data[1:3,2],BMI = Data[1:3,3],Chil = Data[1:3,4],Smoke = Data[1:3,5],Region = Data[1:3,6]),interval = "conf",level = 0.95)

做預(yù)測(cè)只需要把interval="conf"改成interval="pred"即可。

模型診斷

與一元的類似，只是要畫的圖有點(diǎn)多（對(duì)擬合值，每個(gè)解釋變量），這里舉三個(gè)例子：

plot(BMI,resid(mreg01.lm)) abline(h=0)

plot(fitted(mreg01.lm),resid(mreg01.lm)) abline(h=0)

qqnorm(resid(mreg01.lm))

從這三張圖看得出來(lái)，多元線性回歸顯然并不是個(gè)好模型，殘差呈現(xiàn)出來(lái)的模式是非常奇怪的，并且QQ圖明顯不是線性的。在做進(jìn)一步的診斷前，下一篇博客討論是不是每一個(gè)解釋變量都要添加到模型中，以及是因?yàn)橛靡淮芜€是二次的形式。

總結(jié)

以上是生活随笔為你收集整理的UA MATH571A R语言回归分析实践 多元回归2 医疗费用的决定的全部?jī)?nèi)容，希望文章能夠幫你解決所遇到的問(wèn)題。

如果覺得生活随笔網(wǎng)站內(nèi)容還不錯(cuò)，歡迎將生活随笔推薦給好友。