# NOT RUN {
## consider pbc dataset in survival package as an example
library(survival)
dat = pbc[1:312,]
dat$sex = ifelse(dat$sex=='f', 1, 0)
## predciting the event of 'death'
time = dat$time
event = ifelse(dat$status==2, 1, 0)
## standard prediction model: age, bilirubin, and albumin
z.std = as.matrix(subset(dat, select = c(age, bili, albumin)))
## new prediction model: age, bilirubin, albumin, and protime
z.new = as.matrix(subset(dat, select = c(age, bili, albumin, protime)))
## fitting cox models
mstd = coxph(Surv(time, event) ~ ., data.frame(time, event, z.std), x=TRUE)
mnew = coxph(Surv(time, event) ~ ., data.frame(time, event, z.new), x=TRUE)
## Calculation of the risk category NRI at 2000 days
nricens(mdl.std = mstd, mdl.new = mnew, t0 = 2000, cut = c(0.2, 0.4),
niter = 10)
## Next, consider binary prediction models
library(survival)
dat = pbc[1:312,]
dat$sex = ifelse(dat$sex=='f', 1, 0)
## subjects censored before 2000 days are excluded
dat = dat[ dat$time > 2000 | (dat$time < 2000 & dat$status == 2), ]
## predciting the event of 'death' before 2000 days
event = ifelse(dat$time < 2000 & dat$status == 2, 1, 0)
## standard prediction model: age, bilirubin, and albumin
z.std = as.matrix(subset(dat, select = c(age, bili, albumin)))
## new prediction model: age, bilirubin, albumin, and protime
z.new = as.matrix(subset(dat, select = c(age, bili, albumin, protime)))
## glm fit (logistic model)
mstd = glm(event ~ ., binomial(logit), data.frame(event, z.std), x=TRUE)
mnew = glm(event ~ ., binomial(logit), data.frame(event, z.new), x=TRUE)
## Calculation of risk difference NRI
nribin(mdl.std = mstd, mdl.new = mnew, cut = 0.02, niter = 0,
updown = 'diff')
# }
Run the code above in your browser using DataLab