Logistic回归实例2

2744557306

# logistic回归
- g# r  V4 o9 N2 h! g$ n6 D实际上线性最小二乘回归和Logistic回归都是广义线性模型的一个特例。当随机变量Y服从高斯分布，那么得到的是线性最小二乘回归，当随机变量服从伯努利分布，则得到的是Logistic回归。

8 U( W1 C% z8 S! y3 a( dR软件提供了拟合计算广义线性模型的函数glm()，其命令格式如下：fitted.model <- glm(formula, family=family.generator, data=data.frame) 其中，formula是拟合公式；family是分布族，即前面讲到的广义线性模型的种类，如正态分布、Poisson分布、二项分布等。


有了上面这些分布族和连接函数，我们就可以完成相应的广义线性模型的拟合问题。

1）正态分布 正态分布族的使用方法： fm <- glm(formula, family=gaussian(link=identity), data=data.frame) 其中，link=identity可以不写，因为正态分布的连接函数缺省值是恒等（identity）。事实上，整个参数family=gaussian也可以不写，因为分布族的缺省值就是正态分布。 注意：正态分布的广义线性模型实际上与线性模型是相同的，也就是 fm <- glm(formula, family=gaussian, data=data.frame) 与线性模型 fm <- lm(formula, data=data.frame)有完全相同的结果，但效率却低得多。

, t. Y. T' O% c+ H- b/ L3 _2）二项分布
+ J0 G. j" M: A# U* U  I

, k) r- X7 Y8 B, V  {logistic回归模型是一个非线性回归模型，自变量可以是连续变量，也可以是分类变量，或哑变量。但可以使用线性回归模型对参数进行估计，所以Logistic回归模型属于广义线性模型。
  e8 c+ k/ e8 U9 v9 x2 w( d
. E- C7 M# e8 `3 jLogistic回归模型的公式为： fm <- glm(formula, family=binomial(link=logit), data=data.frame) 其中，link=logit可以不写，因为logit是二项分布族连接函数的缺省状态。
  m; J% C6 O; I) x- ^# a0 e实例一、Norell实验，高压电线对牲畜的影响

1. #1、加载数据
   ) Y7 z/ u, J' d+ h
2. norell<-data.frame( x=0:5, n=rep(70,6), success=c(0,9,21,47,60,63) )$ h0 o: w' G& ~( m: U8 {
   
3. norell$Ymat<-cbind(norell$success, norell$n-norell$success)  
   # R4 f, o! Z0 `2 O
4. $ S& c9 b  c& F/ k; ]
   
5. #2、建模
   / _$ a- ^) p. A4 Q+ G5 q
6. glm.sol <- glm(Ymat ~ x, family=binomial, data=norell)5 ]; Y0 t: H. o: d
   
7. 
   & i$ f' s4 z3 F2 g6 t- u; N
8. #3、模型评估
   0 y  l+ A& G5 x# N& N9 x  K
9. summary(glm.sol)

复制代码

1. ##
   ! \6 J( i- h$ _1 p; z1 X0 G& ]7 l
2. ## Call:$ _0 ]3 F2 h: c# W
   
3. ## glm(formula = Ymat ~ x, family = binomial, data = norell)( t- R; ]2 a* Q: M
   
4. ##
   ) W\" z1 c' b: m\" f* g1 k+ `  s
5. ## Deviance Residuals: ) C3 H  V( g0 W# e1 h: |/ e* |
   
6. ##       1        2        3        4        5        6  
   & n. q5 G. f2 _1 t9 s& Z
7. ## -2.2507   0.3892  -0.1466   1.1080   0.3234  -1.6679  / H. ?! X) o! L$ O- B3 c) b
   
8. ##
   # G0 Q% a6 ^6 Y5 j9 s6 J
9. ## Coefficients:
   1 K\" ]8 z# V, H2 |+ U2 A
10. ##             Estimate Std. Error z value Pr(>|z|)   
    9 L4 d\" y6 s. l\" v7 b6 P  Q- Q
11. ## (Intercept)  -3.3010     0.3238  -10.20   <2e-16 ***/ D4 w7 [2 H1 @( w
    
12. ## x             1.2459     0.1119   11.13   <2e-16 ***. _! V) r$ p$ w* R; {+ i! [* ]
    
13. ## ---# r* o$ p; J' e: v% Z. P6 m
    
14. ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1: s1 g! ?9 Y9 Y, Y9 r& B1 X9 Z, B\" a
    
15. ##
    # g  \: C2 [- u# t6 q$ j
16. ## (Dispersion parameter for binomial family taken to be 1)
    7 s+ d- u8 b# h% }: A- Z& J& ?
17. ##
    ( J, {5 n' M) }3 k% @7 p! Z
18. ##     Null deviance: 250.4866  on 5  degrees of freedom
    0 R/ u; I6 `0 C) K( k
19. ## Residual deviance:   9.3526  on 4  degrees of freedom
    / H4 Q9 h% R8 B( D: @7 P
20. ## AIC: 34.093
    # q1 W\" F' H2 C$ r
21. ##
    % \& O4 p, {; S9 Q4 D
22. ## Number of Fisher Scoring iterations: 4

复制代码

1. #与线性回归模型相同，在得到回归模型后，可以作预测：电流强度为3.5毫安时，有响应的牛的概率
   # ~3 O$ l* n6 O6 u3 c
2. 
   & [9 I/ a8 T. j9 N! ~: C' I7 p4 @
3. #4、预测
   ' z% t+ Q' e, Y% m
4. pre <- predict(glm.sol, data.frame(x=3.5))
   ( Q$ D* i5 ~' z9 P* l
5. (p <- exp(pre)/(1+exp(pre)))

复制代码

1. ##        1
   2 {+ U/ B, R6 h8 ]! ^! F
2. ## 0.742642

复制代码

1. #求有50%的牛响应时的电流强度：当P=0.5时，ln(P/(1-P))=0,所以X=-b0/b13 d0 s# O' M4 q$ y: m
   
2. glm.sol$coefficients

复制代码

1. ## (Intercept)           x
   - b1 @+ z/ V\" t2 y$ ]  C: k$ {
2. ##   -3.301035    1.245937

复制代码

1. (X <- -glm.sol$coefficients[[1]]/glm.sol$coefficients[[2]])

复制代码

1. ## [1] 2.649439

复制代码

1. #5、画出响应比例与logistic回归曲线：6 I. l( g+ ^3 ^/ }2 T7 j( ?
   
2. d <- seq(0, 5, length=100)) U7 m& z* R5 R  A) S3 t
   
3. pre <- predict(glm.sol, data.frame(x=d))
   & |, m) l( r\" C) w# I! K5 n
4. p <- exp(pre)/(1+exp(pre))- E\" o* g8 |8 @# a! u+ i
   
5. norell$y <- norell$success/norell$n
   9 F+ ~, y0 Y& ^# `; \* f3 Y
6. plot(norell$x, norell$y)
   + d& q8 p* b2 W# Q$ K, T
7. lines(d, p)

复制代码

1. #其中，d是给出曲线横坐标的点，pre是计算预测值，p是相应的预测概率。用plot函数和lines给出散点图和对应的预测曲线。

复制代码

; P) c! }9 p% _) q