tt<-read.csv("6.06.csv")
library(gtsummary)
library(flextable)
library(pROC)
library(mice)


str(tt)

tt$sex<-as.factor(tt$sex)
tt$smoking.history<-as.factor(tt$smoking.history)
tt$drinking.history<-as.factor(tt$drinking.history)
tt$hemorrhagic.stroke<-as.factor(tt$hemorrhagic.stroke)
tt$ischemic.stroke<-as.factor(tt$ischemic.stroke)
tt$GCS.group<-as.factor(tt$GCS.group)
tt$X3FNENS<-as.factor(tt$WFNS3)
tt$X2FNENS<-as.factor(tt$WFNS2)
tt$aneurysm.locatuion<-as.factor(tt$aneurysm.location)
tt$treatment <-as.factor(tt$treatment )
tt$hypertension<-as.factor(tt$hypertension)
tt$multiple.aneurysm<-as.factor(tt$multiple.aneurysm)
tt$diabetes<-as.factor(tt$diabetes)
tt$hyperlipidemia <-as.factor(tt$hyperlipidemia )
tt$coronary.heart.disease<-as.factor(tt$coronary.heart.disease)
tt$a62<-as.factor(tt$a62)
tt$Hospitalization.day<-as.numeric(tt$Hospitalization.day)
tt$aneurysm.size<-as.numeric(tt$aneurysm.size)



Table1<-tbl_summary(tt, by = y, missing = 'no',
)  %>% 
  add_overall() %>%
  add_p(pvalue_fun = ~style_pvalue(.x, digits = 3))%>%
  separate_p_footnotes() %>% 
  modify_caption("**Table 1. Patient Characteristics**") %>% 
  bold_labels()%>%
  as_flex_table() %>% # 导出 Word
  flextable::save_as_docx(path = 'Table1.docx')
Table1

shapiro.test(tt$age)
t.test(age~y,data=tt ,mu=0,var.equal=F)
shapiro.test(tt$Hospitalization.day)
shapiro.test(tt$aneurysm.size)
shapiro.test(tt$mFisher)
wilcox.test(tt$WFNS,mu=65,alternative="less")
wilcox.test(tt$mFisher,mu=65,alternative="less")



Table.univariate.logistc.regression<-tbl_uvregression(tt, method= glm, y = y, method.args = list(family = binomial),exponentiate = TRUE,
)  %>% 
  add_global_p() %>%
  as_flex_table() %>% # 导出 Word
  flextable::save_as_docx(path = 'Table.univariate logistc regression.docx')
Table.univariate.logistc.regression


model.step <- glm(y ~  WFNS+aneurysm.location+ a62+ GCS.group, data = tt,
                  family = binomial)
model.step1 <- step(model.step,direction='backward')
summary(model.step1)

model.step <- glm(y ~  NLR1  + WFNS+aneurysm.location+ a62+ GCS.group, data = tt,
                  family = binomial)
model.step2 <- step(model.step,direction='backward')
summary(model.step2)

model.step <- glm(y ~  NLR2  + WFNS+aneurysm.location+ a62+ GCS.group, data = tt,
                  family = binomial)
model.step3 <- step(model.step,direction='backward')
summary(model.step3)

model.step <- glm(y ~  NLR3  + WFNS+aneurysm.location+ a62+ GCS.group, data = tt,
                  family = binomial)
model.step4 <- step(model.step,direction='backward')
summary(model.step4)


model1<- glm(y ~  WFNS+aneurysm.location+ a62, data = tt,
             family = binomial)
model1_table<-tbl_regression(model1, exponentiate = TRUE) 

model2<- glm(y ~  NLR1+WFNS+aneurysm.location+ a62, data = tt,
             family = binomial)
model2_table<-tbl_regression(model2, exponentiate = TRUE)

model3<- glm(y ~ NLR2+ WFNS+aneurysm.location+ a62, data = tt,
             family = binomial)
model3_table<-tbl_regression(model3, exponentiate = TRUE)

model4<- glm(y ~  NLR3+WFNS+aneurysm.location+ a62, data = tt,
             family = binomial)
model4_table<-tbl_regression(model4, exponentiate = TRUE)

Table2 <-tbl_merge( tbls = list(model1_table, model2_table, model3_table, model4_table), tab_spanner = c("Model 1", "Model 2", "Model 3", "Model 4"))
Table2 %>%
  as_flex_table() %>% # 导出 Word 
  flextable::save_as_docx(path = 'Table 2.docx') 

model.auc1<- glm(y ~  WFNS+a62, data = tt,
             family = binomial)
model.auc2<- glm(y ~  NLR1+WFNS+ a62, data = tt,
             family = binomial)
model.auc3<- glm(y ~ NLR2+ WFNS+ a62, data = tt,
             family = binomial)
model.auc4<- glm(y ~  NLR3+WFNS+ a62, data = tt,
             family = binomial)



tt$prob1<-predict(model.auc1,newdata = tt,type = 'response')
tt$prob2<-predict(model.auc2,newdata = tt,type = 'response')
tt$prob3<-predict(model.auc3,newdata = tt,type = 'response')
tt$prob4<-predict(model.auc4,newdata = tt,type = 'response')
roc1<-roc(tt$y, tt$prob1)
roc2<-roc(tt$y, tt$prob2)
roc3<-roc(tt$y, tt$prob3)
roc4<-roc(tt$y, tt$prob4)
plot(roc1,col="red",legacy.axes=T,)
plot(roc2,add=TRUE, col="blue")
plot(roc3,add=TRUE, col="green")
plot(roc4,add=TRUE, col="black")
round(auc(roc1),3)##AUC
round(ci(roc1),3)##95%CI
round(auc(roc2),3)##AUC
round(ci(roc2),3)##95%CI
round(auc(roc3),3)##AUC
round(ci(roc3),3)##95%CI
round(auc(roc4),3)##AUC
round(ci(roc4),3)##95%CI
legend(0.6,0.35,
       legend=c("Model 1-auc: 0.669(0.615-0.784)","Model 2-auc: 0.734(0.649-0.820)","Model 3-auc: 0.778(0.693-0.864)","Model 4-auc: 0.739(0.651-0.826)"),
       col=c("red","blue","green","black"),lty=1, text.width =0.47,text.font = 0.5,cex =0.7)


tt<-read.csv("6.06.csv")
library(mice)
str(tt)

tt1<- mice(tt, #数据集
           method = "pmm", #采用pmm
           m=5, # 5次插补
           printFlag = FALSE #不显示历史记录
           ,seed=859929351)


fit1.glm<-with(tt1, glm(y ~WFNS+aneurysm.location+ a62, 
                        family = binomial))
glm.pool1<-pool(fit1.glm)
summary(glm.pool1)


fit2.glm<-with(tt1, glm(y ~  NLR1  + WFNS+aneurysm.location+ a62, 
                        family = binomial))
glm.pool2<-pool(fit2.glm)
summary(glm.pool2) 

fit3.glm<-with(tt1, glm(y ~  NLR2  + WFNS+aneurysm.location+ a62, 
                        family = binomial))

glm.pool3<-pool(fit3.glm)
summary(glm.pool3)

fit4.glm<-with(tt1, glm(y ~  NLR3  + WFNS+aneurysm.location+ a62, 
                        family = binomial))

glm.pool4<-pool(fit4.glm)
summary(glm.pool4)





poolTable<-function(fit1.pool){
  
  Table<-summary(glm.pool1)
  
  OR<-round(exp(Table$estimate),3)
  
  LCI<-round(exp(Table$estimate-1.96*Table$std.error),3)
  
  UCI<-round(exp(Table$estimate+1.96*Table$std.error),3)
  
  β<-round(Table$estimate,3)
  
  S.E<-round(Table$std.error,3)
  
  P<-round(Table$p.value,4)
  
  library(reshape)
  
  library(tidyr)
  
  Table<-rename(Table,c(term="Variable"))
  
  Table<-cbind(Table[1],β,S.E,OR,LCI,UCI,P)
  
  Table$k1<-"(";Table$k2<-"-";Table$k3<-")"
  
  Table<-Table[,c("Variable","β","S.E","OR","k1","LCI","k2","UCI","k3","P")]
  
  Table<-unite(Table,"OR(95%CI)",c(OR,k1,LCI,k2,UCI,k3),sep = "",remove = T)
  
  return(Table)
  
}

## 
Q1<-poolTable(glm.pool1)
write.csv(Q1,file ="Tab1.csv")


poolTable<-function(fit2.pool){
  
  Table<-summary(glm.pool2)
  
  OR<-round(exp(Table$estimate),3)
  
  LCI<-round(exp(Table$estimate-1.96*Table$std.error),3)
  
  UCI<-round(exp(Table$estimate+1.96*Table$std.error),3)
  
  β<-round(Table$estimate,3)
  
  S.E<-round(Table$std.error,3)
  
  P<-round(Table$p.value,4)
  
  library(reshape)
  
  library(tidyr)
  
  Table<-rename(Table,c(term="Variable"))
  
  Table<-cbind(Table[1],β,S.E,OR,LCI,UCI,P)
  
  Table$k1<-"(";Table$k2<-"-";Table$k3<-")"
  
  Table<-Table[,c("Variable","β","S.E","OR","k1","LCI","k2","UCI","k3","P")]
  
  Table<-unite(Table,"OR(95%CI)",c(OR,k1,LCI,k2,UCI,k3),sep = "",remove = T)
  
  return(Table)
  
}

## 
Q2<-poolTable(glm.pool2)
write.csv(Q2,file ="Tab2.csv")

poolTable<-function(fit3.pool){
  
  Table<-summary(glm.pool3)
  
  OR<-round(exp(Table$estimate),3)
  
  LCI<-round(exp(Table$estimate-1.96*Table$std.error),3)
  
  UCI<-round(exp(Table$estimate+1.96*Table$std.error),3)
  
  β<-round(Table$estimate,3)
  
  S.E<-round(Table$std.error,3)
  
  P<-round(Table$p.value,4)
  
  library(reshape)
  
  library(tidyr)
  
  Table<-rename(Table,c(term="Variable"))
  
  Table<-cbind(Table[1],β,S.E,OR,LCI,UCI,P)
  
  Table$k1<-"(";Table$k2<-"-";Table$k3<-")"
  
  Table<-Table[,c("Variable","β","S.E","OR","k1","LCI","k2","UCI","k3","P")]
  
  Table<-unite(Table,"OR(95%CI)",c(OR,k1,LCI,k2,UCI,k3),sep = "",remove = T)
  
  return(Table)
  
}

## 
Q3<-poolTable(glm.pool3)
write.csv(Q3,file ="Tab3.csv")

poolTable<-function(fit4.pool){
  
  Table<-summary(glm.pool4)
  
  OR<-round(exp(Table$estimate),3)
  
  LCI<-round(exp(Table$estimate-1.96*Table$std.error),3)
  
  UCI<-round(exp(Table$estimate+1.96*Table$std.error),3)
  
  β<-round(Table$estimate,3)
  
  S.E<-round(Table$std.error,3)
  
  P<-round(Table$p.value,4)
  
  library(reshape)
  
  library(tidyr)
  
  Table<-rename(Table,c(term="Variable"))
  
  Table<-cbind(Table[1],β,S.E,OR,LCI,UCI,P)
  
  Table$k1<-"(";Table$k2<-"-";Table$k3<-")"
  
  Table<-Table[,c("Variable","β","S.E","OR","k1","LCI","k2","UCI","k3","P")]
  
  Table<-unite(Table,"OR(95%CI)",c(OR,k1,LCI,k2,UCI,k3),sep = "",remove = T)
  
  return(Table)
  
}

## 
Q4<-poolTable(glm.pool4)
write.csv(Q4,file ="Tab4.csv")