基于 R的 广义线性模型分析

只是一些杂乱的笔记,没有注明出处。

简单的逻辑回归模型 logistic regression model,适用于因变量是二分变量的分析。
回归的本质是建立一个模型用来预测,而逻辑回归的独特性在于,预测的结果是只能有两种,true or false。
在R里面做逻辑回归也很简单,只需要构造好数据集,然后用glm函数(广义线性模型(generalized linear model))建模即可,预测用predict函数。

library(DAAG)
head(anesthetic)
cdplot(factor(nomove)~conc,data=anesthetic,main='条件密度图',ylab='病人移动',xlab='麻醉剂量')
#原始数据建模
# glm()是建立广义线性模型的函数。逻辑回归要采用的就是广义线性模型.
anes1=glm(nomove~conc,family=binomial(link='logit'),data=anesthetic,weights= nomove)
#注意这里的weights参数是必须的,因为R无法识别这个占比所基于的基数是多少。
summary(anes1)
#汇总数据
anestot=aggregate(anesthetic[,c('move','nomove')],by=list(conc=anesthetic$conc),FUN=sum)
anestot$conc=as.numeric(as.character(anestot$conc))
anestot$total=apply(anestot[,c('move','nomove')],1,sum)
anestot$prop=anestot$nomove/anestot$total
#汇总数据建模1
anes2=glm(cbind(nomove,move)~conc,family=binomial(link='logit'),data=anestot)
#汇总数据建模2
anes3=glm(prop~conc,family=binomial(link='logit'),weights=total,data=anestot)
# 根据logistic模型,我们可以使用predict函数来预测结果,下面根据上述模型来绘图
x=seq(from=0,to=3,length.out=30)
y=predict(anes1,data.frame(conc=x),type='response')
plot(prop~conc,pch=16,col='red',data=anestot,xlim=c(0.5,3),main='Logistic回归曲线图',ylab='病人静止概率',xlab='麻醉剂量')
lines(y~x,lty=2,col='blue')

可以用两种形式的数据来做分析:

  1. 汇总数据。一种方法是将两种结果的向量合并做为因变量,如anes2模型。
    另一种是将比率做为因变量,总量做为权重进行建模,如anes3模型。这两种建模结果是一样的。
  2. 原始数据,因变量0,1编码。

模型输出的结果,采用summary()来查看,可以看到具体的模型信息,回归系数,系数是否显著等。(需要注意的是,Logistic为非线性模型,回归系数是通过极大似然估计方法计算所得)
summary() 展示拟合模型的详细结果
如果某个变量有多个不同的水平,你想查看一下这个变量总体是否有统计学意义,可以做一下联合检验,用anova()来查看分析偏差表。
anova() 生成一个拟合模型的方差分析表,或者比较两个或更多拟合模型的方差分析表
结果可以这么理解,随着变量从第一个到最后一个逐个加入模型,模型的课解释方差发生变化,通过看p值,是否通过显著性检验。
说明由上述这些变量组成的模型是有意义的,并且是正确的。
anova(object = model2,test = "Chisq")

可以看到这个数据集是关于申请学校是否被录取的,根据学生的GRE成绩,GPA和排名来预测该学生是否被录取。
其中GRE成绩是连续性的变量,学生可以考取任意正常分数。
而GPA也是连续性的变量,任意正常GPA均可。
最后的排名虽然也是连续性变量,但是一般前几名才有资格申请,所以这里把它当做因子,只考虑前四名!
而我们想做这个逻辑回归分析的目的也很简单,就是想根据学生的成绩排名,绩点信息,托福或者GRE成绩来预测它被录取的概率是多少!
结果解读:

> mydata$rank <- factor(mydata$rank)
> mylogit <- glm(admit ~ gre + gpa + rank, data = mydata, family = "binomial")
> summary(mylogit) 
Call:
glm(formula = admit ~ gre + gpa + rank, family = "binomial", data = mydata)
Deviance Residuals: 
     Min       1Q   Median       3Q      Max  
 -1.6268  -0.8662  -0.6388   1.1490   2.0790  
Coefficients:
              Estimate Std. Error z value Pr(>|z|)    
(Intercept) -3.989979   1.139951  -3.500 0.000465 ***
gre          0.002264   0.001094   2.070 0.038465 *  
gpa          0.804038   0.331819   2.423 0.015388 *  
rank2       -0.675443   0.316490  -2.134 0.032829 *  
rank3       -1.340204   0.345306  -3.881 0.000104 ***
rank4       -1.551464   0.417832  -3.713 0.000205 ***
 ---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 
(Dispersion parameter for binomial family taken to be 1) 
    Null deviance: 499.98  on 399  degrees of freedom
Residual deviance: 458.52  on 394  degrees of freedom
AIC: 470.52
Number of Fisher Scoring iterations: 4

根据对这个模型的summary结果可知:
GRE成绩每增加1分,被录取的优势对数(log odds)增加0.002
而GPA每增加1单位,被录取的优势对数(log odds)增加0.804,不过一般GPA相差都是零点几。
最后第二名的同学比第一名同学在其它同等条件下被录取的优势对数(log odds)小了0.675,看来排名非常重要啊!!!
这里必须解释一下这个优势对数(log odds)是什么意思了,如果预测这个学生被录取的概率是p,那么优势对数(log odds)就是log2(p/(1-p)),一般是以自然对数为底。

1

2

3

4

零偏差和剩余偏差之间的差异越大越好。通过分析表格,我们可以看到每次添加一个变量时出现偏差的情况。同样,增加Pclass,Sex and Age可以显着减少残余偏差。这里的大p值表示没有变量的模型或多或少地解释了相同的变化量。最终你想看到的是一个显着的下降和偏差AIC。

glmer()模型下,需要用Anova 大写,Anova(model,type=3)来查看卡方值和p值,这样输出结果和jamovi是相同的。

> Anova(mymodel,type = 3)
Analysis of Deviance Table (Type III Wald chisquare tests)

Response: punishornot
                        Chisq Df Pr(>Chisq)    
(Intercept)            6.2574  1  0.0123674 *  
power                  1.3830  1  0.2395917    
member                36.7568  2  1.043e-08 ***
playera              948.5700  6  < 2.2e-16 ***
power:member          15.0146  2  0.0005491 ***
power:playera         22.7321  6  0.0008914 ***
member:playera        26.1596 12  0.0101894 *  
power:member:playera  20.0121 12  0.0668580 .  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 200,667评论 5 472
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 84,361评论 2 377
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 147,700评论 0 333
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 54,027评论 1 272
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 62,988评论 5 361
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 48,230评论 1 277
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 37,705评论 3 393
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 36,366评论 0 255
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 40,496评论 1 294
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 35,405评论 2 317
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 37,453评论 1 329
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 33,126评论 3 315
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 38,725评论 3 303
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 29,803评论 0 19
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 31,015评论 1 255
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 42,514评论 2 346
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 42,111评论 2 341

推荐阅读更多精彩内容