R语言计算logistic回归C指数,R进行Logistic回归学习笔记

library(stats)

help(infert)                     # Description of data

infert

str(infert)                      # Check type of variables

summary(infert)                  # Statistical summary

## Model1 Develop a simple logistic regression:

model1

## Model output

summary(model1)                    # Output summary information

confint(model1)                    # Output 95% CI for the coefficients

exp(coef(model1))                  # Output OR (exponentiated coefficients)

exp(confint(model1))               # 95% CI for exponentiated coefficients

predict(model1, type="risk")       # predicted values

residuals(model1, type="deviance") # residuals

## Model2 Develop a logistic regression adjusted for other potential confounders:

model2

family = binomial())

summary(model2)

## Model3 Develop a conditional logistic regression

library(survival)

model3

summary(model3)

## Model4 Conditional logistic regression

## Transfer the numerous variable into category variable

model4

data = infert)

summary(model4)

## Model5 Conditional logistic regression

## Conduct a subgroup analysis

## "subset" is not aviable for "clogit"

## create the subset first

str(infert$education)

infert1

model5

data = infert1)

summary(model5)

## Model6 General logistic regression

## Change the reference category for category variable

## Note: Argument "contrasts" not matched in clogit

infert$spontaneous

as.factor(infert$spontaneous)

infert$induced

as.factor(infert$induced)

infert$spontaneous

relevel(infert$spontaneous, ref = "2")

infert$induced

relevel(infert$induced, ref = "1")

model6

data = infert)

summary(model6)

## Model7 General logistic regression

## Performs stepwise model selection by AIC

## Note: Model selection in R is based on AIC but not P value

model7

summary(model7)

model8

summary(model8)

model9

summary(model9)

## Compare the AUC of two ROC curves

infert$pred6

infert$pred8

library("pROC")

rocobj1

infert$pred6,

percent=TRUE,ci=TRUE,col="#1c61b6")

rocobj2

infert$pred8,

add=TRUE,percent=TRUE,ci=TRUE,col="#008600")

testobj

text(50, 50, labels=paste("p-value =",

format.pval(testobj$p.value)), adj=c(0, .5))

legend("bottomright", legend=c("Model5", "Model8"),

col=c("#1c61b6", "#008600"),lwd=2,ncol=2)