### R code from vignette source 'cmpr-DMlate.rnw' ### Encoding: UTF-8 ################################################### ### code chunk number 1: cmpr-DMlate.rnw:2-7 ################################################### options( width=90, # prompt=" ", continue=" ", # Absence of prompts makes it easier for # Students to copy from the final pdf document SweaveHooks=list( fig=function() par(mar=c(3,3,1,1),mgp=c(3,1,0)/1.6,las=1,bty="n") ) ) ################################################### ### code chunk number 2: cmpr-DMlate.rnw:27-34 ################################################### library(survival) library(Epi) library(popEpi) # popEpi::splitMulti returns a data.frame rather than a data.table options("popEpi.datatable" = FALSE) library(tidyverse) clear() ################################################### ### code chunk number 3: cmpr-DMlate.rnw:61-67 ################################################### data(DMlate) # str(DMlate) set.seed(1952) DMlate <- DMlate[sample(1:nrow(DMlate), 2000),] str(DMlate) head(DMlate) ################################################### ### code chunk number 4: cmpr-DMlate.rnw:72-80 ################################################### Ldm <- Lexis(entry = list(per = dodm, age = dodm - dobth, tfd = 0), exit = list(per = dox), exit.status = factor(!is.na(dodth), labels = c("DM","Dead")), data = DMlate) summary(Ldm) ################################################### ### code chunk number 5: cmpr-DMlate.rnw:83-88 ################################################### Cdm <- cutLexis(Ldm, cut = Ldm$dooad, timescale = "per", new.state = "OAD") summary(Cdm) ################################################### ### code chunk number 6: boxCR ################################################### Adm <- subset(Cdm, lex.Cst == "DM") summary(Adm) boxes(Adm, boxpos = TRUE, scale.R = 100, show.BE = TRUE) ################################################### ### code chunk number 7: cmpr-DMlate.rnw:131-140 ################################################### levels(Adm$lex.Xst) m3 <- survfit(Surv(tfd, tfd + lex.dur, lex.Xst) ~ 1, id = lex.id, data = Adm) names(m3) m3$states head(cbind(time = m3$time, m3$pstate)) ################################################### ### code chunk number 8: surv2 ################################################### par( mfrow=c(1,2) ) matplot(m3$time, m3$pstate, type="s", lty=1, lwd=4, col=c("ForestGreen","red","black"), xlim=c(0,15), xaxs="i", ylim=c(0,1), yaxs="i" ) stackedCIF(m3, lwd=3, xlim=c(0,15), xaxs="i", yaxs="i" ) text( rep(12,3), c(0.9,0.3,0.6), levels(Cdm) ) box() ################################################### ### code chunk number 9: surv2 (eval = FALSE) ################################################### ## par( mfrow = c(1,2) ) ## matshade(m3$time, cbind(m3$pstate, ## m3$lower, ## m3$upper)[,c(1,4,7,2,5,8,3,6,9)], ## plot = TRUE, lty = 1, lwd = 4, ## col = c("ForestGreen","red","black"), ## xlim=c(0,15), xaxs="i", ## ylim = c(0,1), yaxs = "i") ## stackedCIF(m3, lwd=3, xlim=c(0,15), xaxs="i", yaxs="i" ) ## text( rep(12,3), c(0.9,0.3,0.6), levels(Cdm) ) ## box() ################################################### ### code chunk number 10: surv3 ################################################### m2 <- survfit(Surv(tfd, tfd + lex.dur, lex.Xst == "OAD" ) ~ 1, data = Adm) M2 <- survfit(Surv(tfd, tfd + lex.dur, lex.Xst == "Dead") ~ 1, data = Adm) par(mfrow = c(1,2)) mat2pol(m3$pstate, c(2,3,1), x = m3$time, col = c("red", "black", "transparent"), xlim=c(0,15), xaxs="i", yaxs = "i", xlab = "time since DM", ylab = "" ) lines(m2$time, 1 - m2$surv, lwd = 3, col = "red" ) mat2pol(m3$pstate, c(3,2,1), x = m3$time, yaxs = "i", col = c("black","red","transparent"), xlim=c(0,15), xaxs="i", yaxs = "i", xlab = "time since DM", ylab = "" ) lines(M2$time, 1 - M2$surv, lwd = 3, col = "black" ) ################################################### ### code chunk number 11: cmpr-DMlate.rnw:313-316 ################################################### Sdm <- splitMulti(Adm, tfd = seq(0,20,0.1) ) summary(Adm) summary(Sdm) ################################################### ### code chunk number 12: cmpr-DMlate.rnw:324-328 ################################################### round(cbind( with(subset(Sdm, lex.Xst == "OAD" ), quantile(tfd + lex.dur, 0:10/10)), with(subset(Sdm, lex.Xst == "Dead"), quantile(tfd + lex.dur, 0:10/10))), 3) ################################################### ### code chunk number 13: cmpr-DMlate.rnw:332-336 ################################################### okn <- c(0,0.5,3,6) dkn <- c(0,2.0,5,9) OAD.glm <- glm.Lexis(Sdm, ~ Ns(tfd, knots = okn), from = "DM", to = "OAD" ) Dead.glm <- glm.Lexis(Sdm, ~ Ns(tfd, knots = dkn), from = "DM", to = "Dead") ################################################### ### code chunk number 14: cmpr-DMlate.rnw:342-346 ################################################### int <- 0.01 nd <- data.frame(tfd = seq(0, 15, int)) l.glm <- ci.pred( OAD.glm, nd) m.glm <- ci.pred(Dead.glm, nd) ################################################### ### code chunk number 15: OAD-mort ################################################### matshade(nd$tfd, cbind(l.glm, m.glm) * 100, plot = TRUE, log = "y", ylim = c(2, 20), col = rep(c("red","black"), 2), lwd = 3) ################################################### ### code chunk number 16: int-ill ################################################### t <- seq(0, 15, 0.01) nd <- data.frame(tfd = t) mu <- ci.pred(Dead.glm, nd)[,1] head(cbind(t, mu)) plot(t, mu, type="l", lwd = 3, xlim = c(0, 7), xaxs = "i", ylim = c(0, max(mu)), yaxs = "i") polygon(t[c(1:501,501:1)], c(mu[1:501], rep(0, 501)), col = "gray", border = "transparent") ################################################### ### code chunk number 17: cmpr-DMlate.rnw:404-407 ################################################### mid <- function(x) x[-1] - diff(x) / 2 (x <- c(1:5, 7, 10)) mid(x) ################################################### ### code chunk number 18: cmpr-DMlate.rnw:415-416 ################################################### sum(diff(t[1:501]) * mid(mu[1:501])) ################################################### ### code chunk number 19: cmpr-DMlate.rnw:430-432 ################################################### Mu <- c(0, cumsum(diff(t) * mid(mu))) head(cbind(t, Mu)) ################################################### ### code chunk number 20: cmpr-DMlate.rnw:480-484 ################################################### cR <- ci.Crisk(mods = list(OAD = OAD.glm, Dead = Dead.glm), nd = nd) str(cR) ################################################### ### code chunk number 21: crisk ################################################### matshade(as.numeric(dimnames(cR$Crisk)[[1]]), cbind(cR$Crisk[,1,], cR$Crisk[,2,], cR$Crisk[,3,]), plot = TRUE, lwd = 2, col = c("limegreen","red","black")) ################################################### ### code chunk number 22: cmpr-DMlate.rnw:503-504 ################################################### mat2pol(cR$Crisk[,3:1,1], col = c("forestgreen","red","black")[3:1]) ################################################### ### code chunk number 23: cmpr-DMlate.rnw:513-518 ################################################### matshade(as.numeric(dimnames(cR$Srisk)[[1]]), cbind(cR$Srisk[,1,], cR$Srisk[,2,]), plot = TRUE, lwd = 2, col = c("black","red"), ylim = 0:1, yaxs = "i") ################################################### ### code chunk number 24: cmpr-DMlate.rnw:553-555 ################################################### str(cR$Stime) round(cR$Stime[1:3*500+1,"Surv",], 1) ################################################### ### code chunk number 25: cmpr-DMlate.rnw:560-562 ################################################### (mY <- matrix(rep(1:3 * 5, 3), 3, 3)) round(100 * cR$Stime[1:3*500+1,"Surv",] / mY, 1) ################################################### ### code chunk number 26: cmpr-DMlate.rnw:567-574 ################################################### time <- as.numeric(dimnames(cR$Stime)[[1]]) / 100 matshade(time, cR$Stime[,"Surv",] / cbind(time, time, time) * 100, plot=TRUE, ylim = 0:1*100, yaxs = "i", xaxs = "i")