devbox
菜鸟追梦旅行
这个屌丝很懒，什么也没留下！
关注作者
热门标签
热门文章
当前位置:   article > 正文

拓端tecdat|R语言广义线性模型GLM、多项式回归和广义可加模型GAM预测泰坦尼克号幸存者_gam glm

作者：菜鸟追梦旅行 | 2024-06-05 15:48:31

赞

踩

gam glm

最近我们被要求撰写关于GLM的研究报告，包括一些图形和统计输出。

 相关视频：非线性模型原理与R语言多项式回归、局部平滑样条、 广义相加模型GAM分析

非线性模型原理与R语言多项式回归、局部平滑样条、 广义相加模型GAM分析

，时长05:41

视频：R语言广义相加模型（GAM）在电力负荷预测中的应用

拓端tecdat：R语言广义相加模型（GAM）在电力负荷预测中的应用

本文通过R语言建立广义线性模型(GLM)、多项式回归和广义可加模型（GAM）来预测谁在1912年的泰坦尼克号沉没中幸存下来。

  1.  
  2. 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.  
  2. Coefficients:
  3. (Intercept)  
  4.     -0.4733  
  5.  
  6. Degrees of Freedom: 890 Total (i.e. Null);  890 Residual
  7. Null Deviance:	    1187 
  8. Residual Deviance: 1187 	AIC: 1189

给出β0的值，并且由于生存概率为

我们通过考虑

  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.  
  2. Coefficients:
  3. (Intercept)      Pclass2      Pclass3  
  4.      0.5306      -0.6394      -1.6704  
  5.  
  6. Degrees of Freedom: 890 Total (i.e. Null);  888 Residual
  7. Null Deviance:	    1187 
  8. Residual Deviance: 1083 	AIC: 1089

由于第1类是参考类，因此我们照旧考虑

  1. exp(0.5306)/(1+exp(0.5306))
  2. [1] 0.629623

predict预测：

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

我们可以尝试概率回归，我们得到的结果是一样的，

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

卡方独立性测试 ：在生存与否之间的检验统计量是多少？

卡方检验的命令如下

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

我们有一个列联表，如果变量是独立的，我们有，然后是统计量，我们可以看到，对测试的贡献

  1.  
  2.             1         2         3
  3.   0 -4.601993 -1.537771  3.993703
  4.   1  5.830678  1.948340 -5.059981

这给了我们很多信息：我们观察到两个正值，分别对应于“幸存”和“头等舱”与“无法幸存”和“三等舱”之间的强（正）关联，以及两个很强的负值，对应于“生存”和“第三等”之间的强烈负相关，以及“无法幸存”和“头等舱”。我们可以在下图上可视化这些值

  1.  
  2. ass(table( Survived, Pclass), shade = TRUE, las=3)

然后我们必须进行逻辑回归，并预测两名模拟乘客的生存概率

假设我们有两名乘客

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

让我们对所有变量进行简单回归，

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

两个变量并不重要。我们删除它们

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

我们有年龄这样的连续变量时，我们可以进行多项式回归

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

但是解释参数变得很复杂。我们注意到三阶项在这里很重要，因此我们将手动进行回归

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

可以看到，p值是相同的。简而言之，将年龄转换为年龄的非线性函数是有意义的。可以可视化此函数

  1.  
  2. plot(xage,yage,xlab="Age",ylab="",type="l")

实际上，我们可以使用样条曲线。广义可加模型（ gam ）是完美的可视化工具

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

我们可以看到年龄变量的变换，

并且我们发现变换接近于我们的3阶多项式。我们可以添加置信带，从而可以验证该函数不是真正的线性

我们现在有三个模型。最后给出了两个模拟乘客的预测，

  1. predict(reg,newdata=newbase,type="response")
  2.         1         2 
  3. 0.9605736 0.1368988 
  4. predict(reg3,newdata=newbase,type="response")
  5.         1         2 
  6. 0.9497834 0.1218426 
  7. predict(regam,newdata=newbase,type="response")
  8.         1         2 
  9. 0.9497834 0.1218426

可以看到莱昂纳多·迪卡普里奥（  Leonardo DiCaprio） 有大约12％的幸存机会（考虑到他的年龄，他有三等票，而且船上没有家人）。

最受欢迎的见解

1.用SPSS估计HLM层次线性模型模型

2.R语言线性判别分析（LDA），二次判别分析（QDA）和正则判别分析（RDA）

3.基于R语言的lmer混合线性回归模型

4.R语言Gibbs抽样的贝叶斯简单线性回归仿真分析

5.在r语言中使用GAM（广义相加模型）进行电力负荷时间序列分析

6.使用SAS，Stata，HLM，R，SPSS和Mplus的分层线性模型HLM

7.R语言中的岭回归、套索回归、主成分回归：线性模型选择和正则化

8.R语言用线性回归模型预测空气质量臭氧数据

9.R语言分层线性模型案例

声明：本文内容由网友自发贡献，不代表【wpsshop博客】立场，版权归原作者所有，本站不承担相应法律责任。如您发现有侵权的内容，请联系我们。转载请注明出处：https://www.wpsshop.cn/w/菜鸟追梦旅行/article/detail/677124
推荐阅读
相关标签

Copyright © 2003-2013 www.wpsshop.cn 版权所有，并保留所有权利。

  

闽ICP备14008679号