天天看點
傳回

拓端tecdat|R語言廣義線性模型GLM、多項式回歸和廣義可加模型GAM預測泰坦尼克号幸存者

R語言廣義線性模型GLM、多項式回歸和廣義可加模型GAM預測泰坦尼克号幸存者

本文通過R語言建立廣義線性模型(GLM)、多項式回歸和廣義可加模型（GAM）來預測誰在1912年的泰坦尼克号沉沒中幸存下來。

  1.  str(titanic)

資料變量為：

  • ​

    ​Survived​

     ​：乘客存活名額（如果存活則為1）
  • ​

    ​Pclass​

     ​：旅客艙位等級
  • ​

    ​Sex​

     ​：乘客性别
  • ​

    ​Age​

     ​：乘客年齡
  • ​

    ​SibSp​

     ​：兄弟姐妹/配偶人數
  • ​

    ​Parch​

     ​：父母/子女人數
  • ​

    ​Embarked​

     ​： 登船港口
  • ​

    ​Name​

     ​：旅客姓名

最後一個變量使用不多，是以我們将其删除，

titanic = titanic[,1:7]

現在，我們回答問題：

幸存的旅客比例是多少？

簡單的答案是

  1.  mean(titanic$Survived)
  2.  [1] 0.3838384

可以在下面的列聯表中找到

  1.  table(titanic$Survived)/nrow(titanic)
  2.  0 1
  3.  0.6161616 0.3838384

或此處幸存者的38.38 ％。也就是說，也可以通過不對任何解釋變量進行邏輯回歸來獲得（換句話說，僅對常數進行回歸）。回歸給出了：

  1.  Coefficients:
  2.  (Intercept)
  3.  -0.4733
  4.  Degrees of Freedom: 890 Total (i.e. Null); 890 Residual
  5.  Null Deviance: 1187
  6.  Residual Deviance: 1187 AIC: 1189

給出β0的值，并且由于生存機率為

拓端tecdat|R語言廣義線性模型GLM、多項式回歸和廣義可加模型GAM預測泰坦尼克号幸存者

我們通過考慮

拓端tecdat|R語言廣義線性模型GLM、多項式回歸和廣義可加模型GAM預測泰坦尼克号幸存者
  1.  exp(-0.4733)/(1+exp(-0.4733))
  2.  [1] 0.3838355

我們也可以使用​

​predict​

​函數 

  1.  predict(glm(Survived~1, family=binomial,type="response")[1]
  2.  1
  3.  0.3838384

此外，在機率回歸中也适用，

  1.  reg=glm(Survived~1, family=binomial(link="probit"),data=titanic)
  2.  predict(reg,type="response")[1]
  3.  1
  4.  0.3838384

幸存的頭等艙乘客的比例是多少？

我們隻看頭等艙的人，

[1] 0.6296296

約63％存活。我們可以進行邏輯回歸

  1.  Coefficients:
  2.  (Intercept) Pclass2 Pclass3
  3.  0.5306 -0.6394 -1.6704
  4.  Degrees of Freedom: 890 Total (i.e. Null); 888 Residual
  5.  Null Deviance: 1187
  6.  Residual Deviance: 1083 AIC: 1089

由于第1類是參考類，是以我們照舊考慮

拓端tecdat|R語言廣義線性模型GLM、多項式回歸和廣義可加模型GAM預測泰坦尼克号幸存者
  1.  exp(0.5306)/(1+exp(0.5306))
  2.  [1] 0.629623

predict預測：

  1.  predict(reg,newdata=data.frame(Pclass="1"),type="response")
  2.  1
  3.  0.6296296

我們可以嘗試機率回歸，我們得到的結果是一樣的，

  1.  predict(reg,newdata=data.frame(Pclass="1"),type="response")
  2.  1
  3.  0.6296296

​卡方​獨立性測試 ：在生存與否之間的檢驗統計量是多少？

卡方檢驗的指令如下

  1.  chisq.test(table( Survived, Pclass))
  2.  Pearson's Chi-squared test
  3.  data: table( Survived, Pclass)
  4.  X-squared = 102.89, df = 2, p-value < 2.2e-16

我們有一個列聯表，如果變量是獨立的，我們有

拓端tecdat|R語言廣義線性模型GLM、多項式回歸和廣義可加模型GAM預測泰坦尼克号幸存者

，然後是統計量

拓端tecdat|R語言廣義線性模型GLM、多項式回歸和廣義可加模型GAM預測泰坦尼克号幸存者

，我們可以看到，對測試的貢獻

拓端tecdat|R語言廣義線性模型GLM、多項式回歸和廣義可加模型GAM預測泰坦尼克号幸存者
  1.  1 2 3
  2.  0 -4.601993 -1.537771 3.993703
  3.  1 5.830678 1.948340 -5.059981

這給了我們很多資訊：我們觀察到兩個正值，分别對應于“幸存”和“頭等艙”與“無法幸存”和“三等艙”之間的強（正）關聯，以及兩個很強的負值，對應于“生存”和“第三等”之間的強烈負相關，以及“無法幸存”和“頭等艙”。我們可以在下圖上可視化這些值

拓端tecdat|R語言廣義線性模型GLM、多項式回歸和廣義可加模型GAM預測泰坦尼克号幸存者
  1.  ass(table( Survived, Pclass), shade = TRUE, las=3)
拓端tecdat|R語言廣義線性模型GLM、多項式回歸和廣義可加模型GAM預測泰坦尼克号幸存者

然後我們必須進行邏輯回歸，并預測兩名模拟乘客的生存機率

假設我們有兩名乘客

newbase = data.frame(
 Pclass = as.factor(c(1,3)),
 Sex = as.factor(c("female","male")),
 Age = c(17,20),
 SibSp = c(1,0),
 Parch = c(2,0),

讓我們對所有變量進行簡單回歸，

Coefficients:
 Estimate Std. Error z value Pr(>|z|)
 (Intercept) 16.830381 607.655774 0.028 0.97790
 Pclass2 -1.268362 0.298428 -4.250 2.14e-05 ***
 Pclass3 -2.493756 0.296219 -8.419 < 2e-16 ***
 Sexmale -2.641145 0.222801 -11.854 < 2e-16 ***
 Age -0.043725 0.008294 -5.272 1.35e-07 ***
 SibSp -0.355755 0.128529 -2.768 0.00564 **
 Parch -0.044628 0.120705 -0.370 0.71159
 EmbarkedC -12.260112 607.655693 -0.020 0.98390
 EmbarkedQ -13.104581 607.655894 -0.022 0.98279
 EmbarkedS -12.687791 607.655674 -0.021 0.98334
 ---
 Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
  
 (Dispersion parameter for binomial family taken to be 1)
  
 Null deviance: 964.52 on 713 degrees of freedom
 Residual deviance: 632.67 on 704 degrees of freedom
 (177 observations deleted due to missingness)
 AIC: 652.67
  
 Number of Fisher Scoring iterations: 13

兩個變量并不重要。我們删除它們

Coefficients:
 Estimate Std. Error z value Pr(>|z|)
 (Intercept) 4.334201 0.450700 9.617 < 2e-16 ***
 Pclass2 -1.414360 0.284727 -4.967 6.78e-07 ***
 Pclass3 -2.652618 0.285832 -9.280 < 2e-16 ***
 Sexmale -2.627679 0.214771 -12.235 < 2e-16 ***
 Age -0.044760 0.008225 -5.442 5.27e-08 ***
 SibSp -0.380190 0.121516 -3.129 0.00176 **
 ---
 Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
  
 (Dispersion parameter for binomial family taken to be 1)
  
 Null deviance: 964.52 on 713 degrees of freedom
 Residual deviance: 636.56 on 708 degrees of freedom
 (177 observations deleted due to missingness)
 AIC: 648.56
  
 Number of Fisher Scoring iterations: 5

我們有年齡這樣的連續變量時，我們可以進行多項式回歸

Coefficients:
 Estimate Std. Error z value Pr(>|z|)
 (Intercept) 3.0213 0.2903 10.408 < 2e-16 ***
 Pclass2 -1.3603 0.2842 -4.786 1.70e-06 ***
 Pclass3 -2.5569 0.2853 -8.962 < 2e-16 ***
 Sexmale -2.6582 0.2176 -12.216 < 2e-16 ***
 poly(Age, 3)1 -17.7668 3.2583 -5.453 4.96e-08 ***
 poly(Age, 3)2 6.0044 3.0021 2.000 0.045491 *
 poly(Age, 3)3 -5.9181 3.0992 -1.910 0.056188 .
 SibSp -0.5041 0.1317 -3.828 0.000129 ***
 ---
 Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
  
 (Dispersion parameter for binomial family taken to be 1)
  
 Null deviance: 964.52 on 713 degrees of freedom
 Residual deviance: 627.55 on 706 degrees of freedom
 AIC: 643.55
  
 Number of Fisher Scoring iterations: 5

但是解釋參數變得很複雜。我們注意到三階項在這裡很重要，是以我們将手動進行回歸

Coefficients:
 Estimate Std. Error z value Pr(>|z|)
 (Intercept) 5.616e+00 6.565e-01 8.554 < 2e-16 ***
 Pclass2 -1.360e+00 2.842e-01 -4.786 1.7e-06 ***
 Pclass3 -2.557e+00 2.853e-01 -8.962 < 2e-16 ***
 Sexmale -2.658e+00 2.176e-01 -12.216 < 2e-16 ***
 Age -1.905e-01 5.528e-02 -3.446 0.000569 ***
 I(Age^2) 4.290e-03 1.854e-03 2.314 0.020669 *
 I(Age^3) -3.520e-05 1.843e-05 -1.910 0.056188 .
 SibSp -5.041e-01 1.317e-01 -3.828 0.000129 ***
 ---
 Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
  
 (Dispersion parameter for binomial family taken to be 1)
  
 Null deviance: 964.52 on 713 degrees of freedom
 Residual deviance: 627.55 on 706 degrees of freedom
 AIC: 643.55
  
 Number of Fisher Scoring iterations: 5

可以看到，p值是相同的。簡而言之，将年齡轉換為年齡的非線性函數是有意義的。可以可視化此函數

  1.  plot(xage,yage,xlab="Age",ylab="",type="l")
拓端tecdat|R語言廣義線性模型GLM、多項式回歸和廣義可加模型GAM預測泰坦尼克号幸存者

實際上，我們可以使用樣條曲線。廣義可加模型（ ​​gam​​ ）是完美的可視化工具

(Dispersion Parameter for binomial family taken to be 1)
  
 Null Deviance: 964.516 on 713 degrees of freedom
 Residual Deviance: 627.5525 on 706 degrees of freedom
 AIC: 643.5525
 177 observations deleted due to missingness
  
 Number of Local Scoring Iterations: 4
  
 Anova for Parametric Effects
 Df Sum Sq Mean Sq F value Pr(>F)
 Pclass 2 26.72 13.361 11.3500 1.407e-05 ***
 Sex 1 131.57 131.573 111.7678 < 2.2e-16 ***
 bs(Age) 3 22.76 7.588 6.4455 0.0002620 ***
 SibSp 1 14.66 14.659 12.4525 0.0004445 ***
 Residuals 706 831.10 1.177
 ---
 Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

我們可以看到年齡變量的變換，

拓端tecdat|R語言廣義線性模型GLM、多項式回歸和廣義可加模型GAM預測泰坦尼克号幸存者

并且我們發現變換接近于我們的3階多項式。我們可以添加置信帶，進而可以驗證該函數不是真正的線性

拓端tecdat|R語言廣義線性模型GLM、多項式回歸和廣義可加模型GAM預測泰坦尼克号幸存者
predict(reg,newdata=newbase,type="response")
 1 2
 0.9605736 0.1368988
 predict(reg3,newdata=newbase,type="response")
 1 2
 0.9497834 0.1218426
 predict(regam,newdata=newbase,type="response")
 1 2
 0.9497834 0.1218426
多項式回歸 r語言 GLM GAM 廣義線性模型
上一篇: punycode和中文互相轉換
下一篇: java xml/bean轉換

繼續閱讀