R语言学习笔记(十一):广义线性模型
#Logistic 回归 install.packages("AER") data(Affairs,package="AER") summary(Affairs)
affairs gender age yearsmarried children
Min. : 0.000 female:315 Min. :17.50 Min. : 0.125 no :171
1st Qu.: 0.000 male :286 1st Qu.:27.00 1st Qu.: 4.000 yes:430
Median : 0.000 Median :32.00 Median : 7.000
Mean : 1.456 Mean :32.49 Mean : 8.178
3rd Qu.: 0.000 3rd Qu.:37.00 3rd Qu.:15.000
Max. :12.000 Max. :57.00 Max. :15.000
religiousness education occupation rating
Min. :1.000 Min. : 9.00 Min. :1.000 Min. :1.000
1st Qu.:2.000 1st Qu.:14.00 1st Qu.:3.000 1st Qu.:3.000
Median :3.000 Median :16.00 Median :5.000 Median :4.000
Mean :3.116 Mean :16.17 Mean :4.195 Mean :3.932
3rd Qu.:4.000 3rd Qu.:18.00 3rd Qu.:6.000 3rd Qu.:5.000
Max. :5.000 Max. :20.00 Max. :7.000 Max. :5.000
table(Affairs$affairs)
0 1 2 3 7 12
451 34 17 19 42 38
#构建因子
Affairs$ynaffair[Affairs$affairs>0]<-1
Affairs$ynaffair[Affairs$affairs==0]<-0
Affairs$ynaffair<-factor(Affairs$ynaffair,levels=c(0,1),labels=c("No","Yes"))
table(Affairs$ynaffair)
No Yes
451 150
#拟合Logistic回归
fit.full<-glm(ynaffair~gender+age+yearsmarried+children+religiousness+education+occupation+rating,data=Affairs,family=binomial())
summary(fit.full)
Call:
glm(formula = ynaffair ~ gender + age + yearsmarried + children +
religiousness + education + occupation + rating, family = binomial(),
data = Affairs)
Deviance Residuals:
Min 1Q Median 3Q Max
-1.5713 -0.7499 -0.5690 -0.2539 2.5191
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 1.37726 0.88776 1.551 0.120807
gendermale 0.28029 0.23909 1.172 0.241083
age -0.04426 0.01825 -2.425 0.015301 *
yearsmarried 0.09477 0.03221 2.942 0.003262 **
childrenyes 0.39767 0.29151 1.364 0.172508
religiousness -0.32472 0.08975 -3.618 0.000297 ***
education 0.02105 0.05051 0.417 0.676851
occupation 0.03092 0.07178 0.431 0.666630
rating -0.46845 0.09091 -5.153 2.56e-07 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 675.38 on 600 degrees of freedom
Residual deviance: 609.51 on 592 degrees of freedom
AIC: 627.51
Number of Fisher Scoring iterations: 4
去除不显著的变量,重新拟合
fit.reduced<-glm(ynaffair~age+yearsmarried+religiousness+rating,data=Affairs,family=binomial())
summary(fit.reduced)
Call:
glm(formula = ynaffair ~ age + yearsmarried + religiousness +
rating, family = binomial(), data = Affairs)
Deviance Residuals:
Min 1Q Median 3Q Max
-1.6278 -0.7550 -0.5701 -0.2624 2.3998
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 1.93083 0.61032 3.164 0.001558 **
age -0.03527 0.01736 -2.032 0.042127 *
yearsmarried 0.10062 0.02921 3.445 0.000571 ***
religiousness -0.32902 0.08945 -3.678 0.000235 ***
rating -0.46136 0.08884 -5.193 2.06e-07 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 675.38 on 600 degrees of freedom
Residual deviance: 615.36 on 596 degrees of freedom
AIC: 625.36
Number of Fisher Scoring iterations: 4
#检验fit.reduced模型,检验新full和reduce模型的效果是否一致
anova(fit.reduced,fit.full,test="Chisq")
Analysis of Deviance Table
Model 1: ynaffair ~ age + yearsmarried + religiousness + rating
Model 2: ynaffair ~ gender + age + yearsmarried + children + religiousness +
education + occupation + rating
Resid. Df Resid. Dev Df Deviance Pr(>Chi)
1 596 615.36
2 592 609.51 4 5.8474 0.2108 #结果不显著,说明这两个模型不独立,效果类似
#模型参数
coef(fit.reduced) #返回结果基于对数
(Intercept) age yearsmarried religiousness rating
1.93083017 -0.03527112 0.10062274 -0.32902386 -0.46136144
exp(coef(fit.reduced))
(Intercept) age yearsmarried religiousness rating
6.8952321 0.9653437 1.1058594 0.7196258 0.6304248 #大于0是正向影响,小于0是负向影响
#评价预测变量对概率的影响概率
testdata<-data.frame(rating=c(1,2,3,4,5),age=mean(Affairs\$age),yearsmarried=mean(Affairs\$yearsmarried),religiousness=mean(Affairs$religiousness))
testdata
rating age yearsmarried religiousness #针对rating进行预测
1 1 32.48752 8.177696 3.116473
2 2 32.48752 8.177696 3.116473
3 3 32.48752 8.177696 3.116473
4 4 32.48752 8.177696 3.116473
5 5 32.48752 8.177696 3.116473
testdata$prob<-predict(fit.reduced,newdata=testdata,type="response")
testdata
rating age yearsmarried religiousness prob #rating从1变为5时,婚外情发生的可能性则从0.5逐渐减小到0.15
1 1 32.48752 8.177696 3.116473 0.5302296
2 2 32.48752 8.177696 3.116473 0.4157377
3 3 32.48752 8.177696 3.116473 0.3096712
4 4 32.48752 8.177696 3.116473 0.2204547
5 5 32.48752 8.177696 3.116473 0.1513079
#过度离势(是指所欲估计参数的实际变异数大于原先期望的变异数,这使我们进行统计推论时出现谬误。)
deviance(fit.reduced)/df.residual(fit.reduced) #检测是否有离群点
[1] 1.03248 #结果为1左右,没有过度离势
#过度离势检测
fit<-glm(ynaffair~age+yearsmarried+religiousness+rating,family=binomial(),data=Affairs)
fit.od<-glm(ynaffair~age+yearsmarried+religiousness+rating,family=quasibinomial(),data=Affairs)
pchisq(summary(fit.od)$dispersion*
fit$df.residual
,fit$df.residual,lower=F)
[1] 0.340122 (P>0.05) #结果为:不显著,说明模型没有过度离势的问题
#泊松回归
install.packages("robust")
data(breslow.dat,package="robust")
names(breslow.dat)
summary(breslow.dat[c(6,7,8,10)])
Base Age Trt sumY
Min. : 6.00 Min. :18.00 placebo :28 Min. : 0.00
1st Qu.: 12.00 1st Qu.:23.00 progabide:31 1st Qu.: 11.50
Median : 22.00 Median :28.00 Median : 16.00
Mean : 31.22 Mean :28.34 Mean : 33.05
3rd Qu.: 41.00 3rd Qu.:32.00 3rd Qu.: 36.00
Max. :151.00 Max. :42.00 Max. :302.00
opar<-par(no.readonly=TRUE)
par(mfrow=c(1,2))
#attach(breslow.dat)
hist(sumY,breaks=20,xlab="Seizure Count",main="Distribution of Seizures")
boxplot(sumY~Trt,xlab="Treatment",main="Group Comparisons")
par(opar)
#拟合泊松回归 fit<-glm(sumY~Base+Age+Trt,data=breslow.dat,family = poisson()) summary(fit)
Call:
glm(formula = sumY ~ Base + Age + Trt, family = poisson(), data = breslow.dat)
Deviance Residuals:
Min 1Q Median 3Q Max
-6.0569 -2.0433 -0.9397 0.7929 11.0061
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 1.9488259 0.1356191 14.370 < 2e-16 ***
Base 0.0226517 0.0005093 44.476 < 2e-16 ***
Age 0.0227401 0.0040240 5.651 1.59e-08 ***
Trtprogabide -0.1527009 0.0478051 -3.194 0.0014 **
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Dispersion parameter for poisson family taken to be 1)
Null deviance: 2122.73 on 58 degrees of freedom
Residual deviance: 559.44 on 55 degrees of freedom
AIC: 850.71
Number of Fisher Scoring iterations: 5
#解释模型参数
coef(fit)
(Intercept) Base Age Trtprogabide
1.94882593 0.02265174 0.02274013 -0.15270095
exp(coef(fit))
(Intercept) Base Age Trtprogabide
7.0204403 1.0229102 1.0230007 0.8583864
#过度离势
deviance(fit)/df.residual(fit)
[1] 10.1717 #大于1,存在过度离势现象
Overdispersion test Obs.Var/Theor.Var Statistic p-value
poisson data 62.87013 3646.468 0 #结果小于0.05,显著,说明存在过度离势现象
#存在过度离势,重新拟合泊松回归
fit.od<-glm(sumY~Base+Age+Trt,data=breslow.dat,family=quasipoisson())
summary(fit.od)
Call:
glm(formula = sumY ~ Base + Age + Trt, family = quasipoisson(),
data = breslow.dat)
Deviance Residuals:
Min 1Q Median 3Q Max
-6.0569 -2.0433 -0.9397 0.7929 11.0061
Coefficients:
Estimate Std. Error t value Pr(>|t|) #重新拟合后,估计参数没有变化,但是标准误差变大了不少
(Intercept) 1.948826 0.465091 4.190 0.000102 ***
Base 0.022652 0.001747 12.969 < 2e-16 ***
Age 0.022740 0.013800 1.648 0.105085
Trtprogabide -0.152701 0.163943 -0.931 0.355702
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Dispersion parameter for quasipoisson family taken to be 11.76075)
Null deviance: 2122.73 on 58 degrees of freedom
Residual deviance: 559.44 on 55 degrees of freedom
AIC: NA
Number of Fisher Scoring iterations: 5
转载于:https://www.cnblogs.com/GhostBear/p/7760466.html
R语言学习笔记(十一):广义线性模型相关推荐
- R语言学习笔记——高级篇:第十四章-主成分分析和因子分析
R语言 R语言学习笔记--高级篇:第十四章-主成分分析和因子分析 文章目录 R语言 前言 一.R中的主成分和因子分析 二.主成分分析 2.1.判断主成分的个数 2.2.提取主成分 2.3.主成分旋转 ...
- R语言学习笔记(1~3)
R语言学习笔记(1~3) 一.R语言介绍 x <- rnorm(5) 创建了一个名为x的向量对象,它包含5个来自标准正态分布的随机偏差. 1.1 注释 由符号#开头. #函数c()以向量的形式输 ...
- r语言c函数怎么用,R语言学习笔记——C#中如何使用R语言setwd()函数
在R语言编译器中,设置当前工作文件夹可以用setwd()函数. > setwd("e://桌面//") > setwd("e:\桌面\") > ...
- R语言学习笔记 07 Probit、Logistic回归
R语言学习笔记 文章目录 R语言学习笔记 probit回归 factor()和as.factor() relevel() 案例11.4复刻 glm函数 整理变量 回归:Logistic和Probit- ...
- R语言学习笔记——入门篇:第一章-R语言介绍
R语言 R语言学习笔记--入门篇:第一章-R语言介绍 文章目录 R语言 一.R语言简介 1.1.R语言的应用方向 1.2.R语言的特点 二.R软件的安装 2.1.Windows/Mac 2.2.Lin ...
- R语言学习笔记——入门篇:第三章-图形初阶
R语言 R语言学习笔记--入门篇:第三章-图形初阶 文章目录 R语言 一.使用图形 1.1.基础绘图函数:plot( ) 1.2.图形控制函数:dev( ) 补充--直方图函数:hist( ) 补充- ...
- R语言学习笔记 06 岭回归、lasso回归
R语言学习笔记 文章目录 R语言学习笔记 比较lm.ridge和glmnet函数 画岭迹图 图6-4 <统计学习导论 基于R语言的应用>P182 图6-6<统计学习导论 基于R语言的 ...
- R语言学习笔记(八)--读写文件与网络爬虫
R语言学习笔记(八) 1 工作路径 2 保存R对象 3 Scan函数 3-1 从控制台读取数据 3-2 从txt文件读取数据 3-3 从url读取数据 4 按行读写文本文件 5 读取文本文件(txt. ...
- R语言学习笔记(三)多元数据的数据特征、相关分析与图形表示
文章目录 写在前面 独立性检验 χ2\chi^2χ2独立性检验 Fisher独立性检验 Cochran-Mantel-Haenszel χ2\chi^2χ2独立性检验 相关性分析 相关性检验 相关性检 ...
- R语言学习笔记(六)回归分析
文章目录 写在前面 普通最小二乘(OLS)回归法 正态假设 简单线性回归 多项式回归 多元线性回归 有交互项的多元线性回归 小结 回归诊断 标准方法 综合验证方法 多重共线性 广义线性回归--Logi ...
最新文章
- 2022-2028年中国新零售行业深度调研及投资前景预测报告(全卷)
- python 接口数据驱动_python接口测试实例--数据驱动(程序与数据分离)
- python 两个字典的合并 update
- 2021数据技术嘉年华 • 嘉宾面对面
- iterator遍历_HashMap 的 7 种遍历方式与性能分析!(强烈推荐)
- esp32 io速度_乐鑫科技发布 ESP32-S3 芯片,精准聚焦 AIoT 市场
- Javascript第五章为什么用firstChild获取table中最后一个节点会取到text或者tbody第十一课
- 数学建模之初等模型详解
- matlab能打开mdl文件吗,simulink打开mdl文件的问题
- Python 代码库之Tuple如何append添加元素
- Photoshop的安装教程
- 电气火灾监控系统在地铁供配电系统中的应用
- Unity SetFromToRotation和FromToRotation的区别
- RTX2012概述-1
- LeedCode 717 1比特与2比特字符
- python语义分割数据标签,将数字标签转彩色标签
- openldap 集成 sssd
- 软件“生命”系统进化论——软件以负熵为生
- 共模电感(扼流圈)选型
- 如何解除word文档保护的方法