##目的是构建P=P(Y=1 |X)与影响因素x之间的关系，出现阳性概率的结果
rm(list = ls())
library(MPV)
knitr::kable(head(p13.2))
log_fix <- glm(y~x,family = "binomial",data = p13.2)
summary(log_fix)

#结果解读
exp(coef(log_fix)[2])#结果=1.000201，即x每增加一个单位，y相应的增加1.00020倍
exp(confint(log_fix)[2,])##95%CI,若pvalue<0.05，则95%不包括1，反之包括1

φ <- log_fix$deviance/log_fix$df.residual

log_fix <- glm(y~x,family = quasibinomial(),data = p13.2)
summary(log_fix)

rm(list = ls())
library(foreign)
library(nnet)
ml <- read.dta("http://stats.idre.ucla.edu/stat/data/hsbdemo.dta")
head(ml)
ml$prog2 <- relevel(ml$prog,ref = "academic")
test <- multinom(prog2~ses + write,data = ml)
summary(test)
##计算pvalue
z <- summary(test)$coefficients/summary(test)$standard.errors
p <- (1 - pnorm(abs(z),0,1))*2#双侧，所以*2
p
exp(coef(test))#得到OR值

#因变量有多个，且为有序变量
require(foreign)
require(MASS)
dat <- read.dta("http://stats.idre.ucla.edu/stat/data/ologit.dta")
head(dat)
#查看数据
ftable(xtabs(~public+apply+pared,data = dat))m <- polr(apply~pared+public+gpa,data = dat,Hess = TRUE,method = "logistic")
summary(m)
ctable <- coef(summary(m))
(ctable <- coef(summary(m)))#计算p值
p <- pnorm(abs(ctable[,"t value"]),lower.tail = FALSE)*2
(ctable <- cbind(ctable,"p value"=p))(ci <- confint(m))#95%CIexp(cbind(OR=coef(m),ci))#OR值##解读：pared↑ 1 unit → OR=2.85（class1→class2；class2→class3）
##OR c1vsCk>1 = OR c2vsCk>2

#1：M case control  1个病例；M m个对照组，m<=4
library(faraway)
data(amlxray,package = "faraway")
head(amlxray)
str(amlxray)
library(survivv)
cmod <- clogit(disease~Sex+Mray+Fray+CnRay+strata(ID),data = amlxray)##strate函数中的ID如果自
#自己数据中没有，需要加入进去。summary(cmod)
##对于条件logistic，重点是需要知道什么时候需要采用条件logistic回归。