R语言与概率统计(三) 多元统计分析(下)广义线性回归

广义线性回归

> life<-data.frame(
+   X1=c(2.5, 173, 119, 10, 502, 4, 14.4, 2, 40, 6.6,
+        21.4, 2.8, 2.5, 6, 3.5, 62.2, 10.8, 21.6, 2, 3.4,
+        5.1, 2.4, 1.7, 1.1, 12.8, 1.2, 3.5, 39.7, 62.4, 2.4,
+        34.7, 28.4, 0.9, 30.6, 5.8, 6.1, 2.7, 4.7, 128, 35,
+        2, 8.5, 2, 2, 4.3, 244.8, 4, 5.1, 32, 1.4),
+   X2=rep(c(0, 2, 0, 2, 0, 2, 0, 2, 0, 2, 0, 2, 0, 2, 0, 2,
+            0, 2, 0, 2, 0, 2, 0),
+          c(1, 4, 2, 2, 1, 1, 8, 1, 5, 1, 5, 1, 1, 1, 2, 1,
+            1, 1, 3, 1, 2, 1, 4)),
+   X3=rep(c(0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1),
+          c(6, 1, 3, 1, 3, 1, 1, 5, 1, 3, 7, 1, 1, 3, 1, 1, 2, 9)),
+   Y=rep(c(0,  1,   0,  1), c(15, 10, 15, 10))
+ )
> glm.sol<-glm(Y~X1+X2+X3, family=binomial, data=life)
> summary(glm.sol)

Call:
glm(formula = Y ~ X1 + X2 + X3, family = binomial, data = life)

Deviance Residuals:
    Min       1Q   Median       3Q      Max
-1.6960  -0.5842  -0.2828   0.7436   1.9292  

Coefficients:
             Estimate Std. Error z value Pr(>|z|)
(Intercept) -1.696538   0.658635  -2.576 0.010000 **
X1           0.002326   0.005683   0.409 0.682308
X2          -0.792177   0.487262  -1.626 0.103998
X3           2.830373   0.793406   3.567 0.000361 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 67.301  on 49  degrees of freedom
Residual deviance: 46.567  on 46  degrees of freedom
AIC: 54.567

Number of Fisher Scoring iterations: 5

可见拟合的效果不好

> pre<-predict(glm.sol, data.frame(X1=5,X2=2,X3=0))
> p<-exp(pre)/(1+exp(pre));p#不接受治疗
         1
0.03664087
>
> pre<-predict(glm.sol, data.frame(X1=5,X2=2,X3=1))
> p<-exp(pre)/(1+exp(pre));p#接受治疗
        1
0.3920057
>
> step(glm.sol)
Start:  AIC=54.57
Y ~ X1 + X2 + X3

       Df Deviance    AIC
- X1    1   46.718 52.718
<none>      46.567 54.567
- X2    1   49.502 55.502
- X3    1   63.475 69.475

Step:  AIC=52.72
Y ~ X2 + X3

       Df Deviance    AIC
<none>      46.718 52.718
- X2    1   49.690 53.690
- X3    1   63.504 67.504

Call:  glm(formula = Y ~ X2 + X3, family = binomial, data = life)

Coefficients:
(Intercept)           X2           X3
     -1.642       -0.707        2.784  

Degrees of Freedom: 49 Total (i.e. Null);  47 Residual
Null Deviance:	    67.3
Residual Deviance: 46.72 	AIC: 52.72
> glm.new<-update(glm.sol, .~.-X1)
> summary(glm.new)

Call:
glm(formula = Y ~ X2 + X3, family = binomial, data = life)

Deviance Residuals:
    Min       1Q   Median       3Q      Max
-1.6849  -0.5949  -0.3033   0.7442   1.9073  

Coefficients:
            Estimate Std. Error z value Pr(>|z|)
(Intercept)  -1.6419     0.6381  -2.573 0.010082 *
X2           -0.7070     0.4282  -1.651 0.098750 .
X3            2.7844     0.7797   3.571 0.000355 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 67.301  on 49  degrees of freedom
Residual deviance: 46.718  on 47  degrees of freedom
AIC: 52.718

Number of Fisher Scoring iterations: 5

>
> pre<-predict(glm.sol, data.frame(X1=5,X2=2,X3=0))
> p<-exp(pre)/(1+exp(pre));p#不接受治疗
         1
0.03664087
>
> pre<-predict(glm.sol, data.frame(X1=5,X2=2,X3=1))
> p<-exp(pre)/(1+exp(pre));p#接受治疗
        1
0.3920057
#####再来看一个类似的问题
install.packages(‘AER‘)
data(Affairs,package=‘AER‘)#婚外情数据,包括9个变量,婚外斯通频率,性别,婚龄等。
summary(Affairs)
table(Affairs$affairs)
#我们感兴趣的是是否有过婚外情所以做如下处理
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)
#接下来做逻辑回归
fit.full=glm(ynaffair~.-affairs,data=Affairs,family=binomial())
summary(fit.full)
#除掉较大p值所对应的变量,如性别,是否有孩子、学历和职业在做一次分析
fit.reduced=glm(ynaffair~age+yearsmarried+religiousness+rating,data=Affairs,family=binomial())
summary(fit.reduced)

AIC(fit.full,fit.reduced)#模型比较

#系数解释
exp(coef(fit.reduced))

原文地址:https://www.cnblogs.com/caiyishuai/p/11164543.html

时间: 2024-10-08 18:03:10

R语言与概率统计(三) 多元统计分析(下)广义线性回归的相关文章

R语言与概率统计(三) 多元统计分析

> #############6.2一元线性回归分析 > x<-c(0.10,0.11,0.12,0.13,0.14,0.15,0.16,0.17,0.18,0.20,0.21,0.23) > y<-c(42.0,43.5,45.0,45.5,45.0,47.5,49.0,53.0,50.0,55.0,55.0,60.0) > plot(x~y) > lm.sol<-lm(y ~ x) > summary(lm.sol) Call: lm(formul

R语言结合概率统计的体系分析---数字特征

现在有一个人,如何对这个人怎么识别这个人?那么就对其存在的特征进行提取,比如,提取其身高,其相貌,其年龄,分析这些特征,从而确定了,这个人就是这个人,我们绝不会认错. 同理,对数据进行分析,也是提取出数据的特征,对其特征进行分析,从而确定这些数据所呈现的信息状况,从而确定了这些数据的独特性和唯一性,因为他呈现的信息是唯一的,绝不与别的是相同的. 那么这些特征是什么呢?拥有哪些特征呢?似乎应该是经过无数科学家的总结,终于发现了几个重要的特征,包括数字特征和分布特征,这个数字特征,包括集中位置,分散

R语言与概率统计(六) 主成分分析 因子分析

超高维度分析,N*P的矩阵,N为样本个数,P为指标,N<<P PCA:抓住对y对重要的影响因素 主要有三种:PCA,因子分析,回归方程+惩罚函数(如LASSO) 为了降维,用更少的变量解决问题,如果是二维的,那么就是找到一条线,要使这些点再线上的投影最大,投影最大,就是越分散,就考虑方差最大. 原文地址:https://www.cnblogs.com/caiyishuai/p/11169073.html

R语言与医学统计图形-【14】ggplot2几何对象之直方密度图

ggplot2绘图系统--几何对象之直方图.密度图 1.直方图 参数. geom_histogram(mapping = , data = , stat = 'bin', #统计变换,概率密度为density position = 'stack', binwidth = , #条柱宽度 bins = , #条柱数目,默认30 na.rm = FALSE, show.legend = NA, inherit.aes = TRUE ) 示例. ggplot(diamonds,aes(carat))+

R语言高性能编程(三)

一.使用并行计算加倍提升性能1.数据并行 VS 任务并行实现数据并行的算法scoket 并行性注意并行计算时间并不与执行任务的计算资源数目成正比(计算机核心),amdahl定律:并行代码的速度受限于串行执行的部分,包括并行性带来的开销在非windows系统中,parallel支持分叉集群(交叉法),新的work进程会从父R进程分叉出来,并拷贝数据.好处是不需要显示的创建和销毁集群实现任务并行的算法 2.计算机集群并行执行多个任务只有基于socket的集群可以做到这一点,因为进程不可能被分叉到另外

《概率统计》多元随机变量

楔子 前两篇我们讨论的离散型和连续型随机变量都是单一变量,然而在现实当中,一个试验常常会涉及到多个随机变量.所谓多个随机变量是指在同一个试验结果之下产生的多个随机变量.这些随机变量的取值是由试验结果确定的,因此它们的取值会存在相互关联.这里我们先以离散型随机变量为例,将离散型随机变量的分布列和期望推广到多个随机变量的情况,并且进一步在此基础上讨论多元随机变量条件和独立的重要概念. 好了,此刻我们假设试验中不再只有一个随机变量,而是两个随机变量 X 和 Y,同时描述他们俩的取值概率,我们用什么方式

R语言与医学统计图形-【12】ggplot2几何对象之条图

ggplot2绘图系统--几何对象之条图(包括误差条图) 1.条图 格式: geom_bar(mapping = , data = , stat = 'count', #统计变换默认计数 position = 'stack', #默认堆栈 width = , #条形宽度 binwidth = , na.rm = FALSE, show.legend = , inherit.aes = TRUE) positon: dodge并排 fill堆叠填充标准化为1 stack堆栈 identity不做调

R语言实战读书笔记(七)基本统计分析

summary() sapply(x,fun,options):对数据框或矩阵中的每一个向量进行统计 mean sd:标准差 var:方差 min: max: median: length: range: quantile: vars <- c("mpg", "hp", "wt")head(mtcars[vars]) summary(mtcars[vars]) mystats <- function(x, na.omit = FALS

R语言:常用统计一些方法代码

理论漫衍依赖于若干未知参数时Kolmogorov-Smirnov 检讨ks.test()例一 对一台设备举办寿命检讨,记录十次无妨碍操纵时间,并按从小到大的序次分列如下,用ks检讨要领检讨此设备无妨碍事情时间是否切合rambda=1/1500的指数漫衍呼吁:X<-c(420, 500, 920, 1380, 1510, 1650, 1760, 2100, 2300, 2350)ks.test(X, "pexp", 1/1500)例二 假设从漫衍函数F(x)和G(x)的总体中别离随