第5章 广义及一般性模型
Logistic 模型
library(readxl)
d5.1=read_excel('C:/Users/12290/Desktop/多元统计表格/mvstats5.xlsx','d5.1')
head(d5.1)
## # A tibble: 6 × 4
## y x1 x2 x3
## <dbl> <dbl> <dbl> <dbl>
## 1 1 1 17 1
## 2 0 1 44 0
## 3 0 1 48 1
## 4 0 1 55 0
## 5 1 1 75 1
## 6 1 0 35 0
logit.glm<-glm(y~x1+x2+x3,family=binomial,data=d5.1)
summary(logit.glm)
##
## Call:
## glm(formula = y ~ x1 + x2 + x3, family = binomial, data = d5.1)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.5636 -0.9131 -0.7892 0.9637 1.6000
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.597610 0.894831 0.668 0.5042
## x1 -1.496084 0.704861 -2.123 0.0338 *
## x2 -0.001595 0.016758 -0.095 0.9242
## x3 0.315865 0.701093 0.451 0.6523
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 62.183 on 44 degrees of freedom
## Residual deviance: 57.026 on 41 degrees of freedom
## AIC: 65.026
##
## Number of Fisher Scoring iterations: 4
logit.step<-step(logit.glm,direction="both")
## Start: AIC=65.03
## y ~ x1 + x2 + x3
##
## Df Deviance AIC
## - x2 1 57.035 63.035
## - x3 1 57.232 63.232
## <none> 57.026 65.026
## - x1 1 61.936 67.936
##
## Step: AIC=63.03
## y ~ x1 + x3
##
## Df Deviance AIC
## - x3 1 57.241 61.241
## <none> 57.035 63.035
## + x2 1 57.026 65.026
## - x1 1 61.991 65.991
##
## Step: AIC=61.24
## y ~ x1
##
## Df Deviance AIC
## <none> 57.241 61.241
## + x3 1 57.035 63.035
## + x2 1 57.232 63.232
## - x1 1 62.183 64.183
summary(logit.step)
##
## Call:
## glm(formula = y ~ x1, family = binomial, data = d5.1)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.4490 -0.8782 -0.8782 0.9282 1.5096
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.6190 0.4688 1.320 0.1867
## x1 -1.3728 0.6353 -2.161 0.0307 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 62.183 on 44 degrees of freedom
## Residual deviance: 57.241 on 43 degrees of freedom
## AIC: 61.241
##
## Number of Fisher Scoring iterations: 4
pre1<-predict(logit.step,data.frame(x1=1))
p1<-exp(pre1)/(1+exp(pre1))
pre2<-predict(logit.step,data.frame(x1=0))
p2<-exp(pre2)/(1+exp(pre2))
c(p1,p2)
## 1 1
## 0.32 0.65
对数线性模型
d5.2=read_excel('C:/Users/12290/Desktop/多元统计表格/mvstats5.xlsx','d5.2')
head(d5.2)
## # A tibble: 6 × 3
## y x1 x2
## <dbl> <dbl> <dbl>
## 1 53 1 1
## 2 434 2 1
## 3 111 3 1
## 4 38 1 2
## 5 108 2 2
## 6 48 3 2
log.glm<-glm(y~x1 + x2,family=poisson(link=log),data=d5.2)
summary(log.glm)
##
## Call:
## glm(formula = y ~ x1 + x2, family = poisson(link = log), data = d5.2)
##
## Deviance Residuals:
## 1 2 3 4 5 6
## -10.784 14.444 -8.468 -2.620 4.960 -3.142
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 6.15687 0.14196 43.371 < 2e-16 ***
## x1 0.12915 0.04370 2.955 0.00312 **
## x2 -1.12573 0.08262 -13.625 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 662.84 on 5 degrees of freedom
## Residual deviance: 437.97 on 3 degrees of freedom
## AIC: 481.96
##
## Number of Fisher Scoring iterations: 5
完全随机设计模型
d5.3=read_excel('C:/Users/12290/Desktop/多元统计表格/mvstats5.xlsx','d5.3')
head(d5.3)
## # A tibble: 6 × 2
## Y A
## <dbl> <dbl>
## 1 2.36 1
## 2 2.38 1
## 3 2.48 1
## 4 2.45 1
## 5 2.47 1
## 6 2.43 1
anova(lm(Y~factor(A),data=d5.3))
## Analysis of Variance Table
##
## Response: Y
## Df Sum Sq Mean Sq F value Pr(>F)
## factor(A) 2 0.122233 0.061117 40.534 8.94e-07 ***
## Residuals 15 0.022617 0.001508
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
d5.4=read_excel('C:/Users/12290/Desktop/多元统计表格/mvstats5.xlsx','d5.4')
head(d5.4)
## # A tibble: 6 × 3
## Y A B
## <dbl> <dbl> <dbl>
## 1 582 1 1
## 2 491 2 1
## 3 601 3 1
## 4 758 4 1
## 5 562 1 2
## 6 541 2 2
anova(lm(Y~factor(A)+factor(B),data=d5.4))
## Analysis of Variance Table
##
## Response: Y
## Df Sum Sq Mean Sq F value Pr(>F)
## factor(A) 3 15759 5253 0.4306 0.7387
## factor(B) 2 22385 11192 0.9174 0.4491
## Residuals 6 73198 12200
随机区组设计模型
d5.5=read_excel('C:/Users/12290/Desktop/多元统计表格/mvstats5.xlsx','d5.5')
head(d5.5)
## # A tibble: 6 × 3
## Y A B
## <dbl> <dbl> <dbl>
## 1 52 1 1
## 2 48 1 1
## 3 44 1 1
## 4 44 1 1
## 5 84 1 2
## 6 88 1 2
anova(lm(Y~factor(A)+factor(B),data=d5.5))
## Analysis of Variance Table
##
## Response: Y
## Df Sum Sq Mean Sq F value Pr(>F)
## factor(A) 1 1600 1600.00 14.804 0.002017 **
## factor(B) 1 2500 2500.00 23.132 0.000341 ***
## Residuals 13 1405 108.08
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
正交设计模型
d5.6=read_excel('C:/Users/12290/Desktop/多元统计表格/mvstats5.xlsx','d5.6')
d5.6
## # A tibble: 8 × 5
## A B C D Y
## <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 1 1 1 1 86
## 2 1 1 2 2 95
## 3 1 2 1 2 91
## 4 1 2 2 1 94
## 5 2 1 1 2 91
## 6 2 1 2 1 96
## 7 2 2 1 1 83
## 8 2 2 2 2 88
anova(lm(Y~factor(A)+factor(B),data=d5.6))
## Analysis of Variance Table
##
## Response: Y
## Df Sum Sq Mean Sq F value Pr(>F)
## factor(A) 1 8 8 0.3333 0.5887
## factor(B) 1 18 18 0.7500 0.4261
## Residuals 5 120 24