r - error message when ploting subjects at risk with survplot -


i following error message when trying plot subjects @ risk along x-axis in survplot:

error in text.default(tt[-1], yy, nri[-1], cex = cex.n.risk, adj = 1) : zero-length 'labels' specified 

any help? rather new survival analysis, , failed find explanation of error. code in general seems fine, except when add n.risk=true option plot, error comes up. hint appreciated. many thanks.

below data code used.

here data

duration <-  structure(list(conflict = c("angola 75-89", "angola 89-91", "angola 92-94",      "azerb (n-k) 89-94", "bosnia 92-95", "cambodia 70-91", "chad 79-79",      "chad 89-96", "chechnya 94-96", "colombia 48-57", "croatia 91-91 (?)",      "croatia 95-95", "domrep 65-65", "el salv 79-91", "georga 89-92",      "georgb 92-94", "guatem 68-96", "india 46-48", "iraq 61-70",      "laos 59-73", "lebanon 58-58", "lebanon 75-89", "liberia 89-93",      "malaysia 48-56", "moldova 92-92", "mozamb 81-92", "nicara 81-89",      "phil 72-96", "rwanda 90-93", "sieleo 91-96", "stafrica 83-91",      "sudan 63-72", "tajik 92-97", "yemen 62-70", "zimbab 72-79",      "guinea-bissau june - november 1998", "liberia 94-96", "papua new guinea 1990 - 2001",      "afghanistan 1978 - 2001", "ethiopia 1961-1993", "indonesia (aceh) 1976 - 2005",      "kenya 2007- 2008", "nepal 1996 - 2006", "somalia 1991 - 2008",      "bangladesh 1997", "burundi 1993-2005", "cote d'ivoire 2002-2007",      "democratic republic of congo 98-03", "northern ireland (68-98)",      "darfur, sudan 2003-2010", "sudan 83-05", "liberia 1999-2003"     ), peacedur = c(2, 17, 58, 175, 157, 206, 7, 117, 34, 322, 43,      157, 520, 204, 192, 171, 144, 0.100000001490116, 48, 25, 199,      230, 12, 626, 89, 195, 232, 148, 8, 6, 204, 141, 138, 357, 348,      122, 40, 23, 0.100000001490116, 60, 40, 8, 24, 0.100000001490116,      133, 28, 22, 69, 128.5, 3, 71, 83), peacefail = c(1, 1, 1, 0,      0, 0, 1, 1, 1, 1, 1, 0, 0, 0, 1, 1, 0, 1, 1, 1, 1, 0, 1, 0, 0,      0, 0, 0, 1, 1, 0, 1, 0, 0, 1, 0, 0, 0, 1, 1, 0, 0, 0, 1, 0, 0,      0, 1, 0, 1, 1, 0), totalps = c(0l, 2l, 3l, 2l, 3l, 2l, 1l, 2l,      4l, 2l, 1l, 1l, 1l, 3l, 2l, 2l, 2l, 2l, 3l, 1l, 2l, 2l, 1l, 2l,      3l, 3l, 4l, 4l, 3l, 3l, 4l, 3l, 2l, 4l, 3l, 2l, 1l, 1l, 2l, 2l,      4l, 1l, 3l, 1l, 3l, 3l, 3l, 2l, 2l, 4l, 3l, 3l), year_end = c(1989l,      1991l, 1994l, 1994l, 1995l, 1991l, 1979l, 1996l, 1996l, 1957l,      1991l, 1995l, 1965l, 1991l, 1992l, 1994l, 1996l, 1948l, 1970l,      1973l, 1958l, 1989l, 1993l, 1956l, 1992l, 1992l, 1989l, 1996l,      1993l, 1996l, 1991l, 1972l, 1997l, 1970l, 1979l, 1998l, 1996l,      2001l, 2001l, 1993l, 2005l, 2008l, 2006l, 2008l, 1997l, 2005l,      2007l, 2003l, 1998l, 2010l, 2005l, 2003l), peacedur.year = c(1,      2, 5, 15, 14, 18, 1, 10, 3, 27, 4, 14, 44, 17, 16, 15, 12, 1,      4, 3, 17, 20, 1, 53, 8, 17, 20, 13, 1, 1, 17, 12, 12, 30, 29,      11, 4, 2, 1, 5, 4, 1, 2, 1, 12, 3, 2, 6, 11, 1, 6, 7), survobj = structure(c(2,      17, 58, 175, 157, 206, 7, 117, 34, 322, 43, 157, 520, 204, 192,      171, 144, 0.100000001490116, 48, 25, 199, 230, 12, 626, 89, 195,      232, 148, 8, 6, 204, 141, 138, 357, 348, 122, 40, 23, 0.100000001490116,      60, 40, 8, 24, 0.100000001490116, 133, 28, 22, 69, 128.5, 3,      71, 83, 1, 1, 1, 0, 0, 0, 1, 1, 1, 1, 1, 0, 0, 0, 1, 1, 0, 1,      1, 1, 1, 0, 1, 0, 0, 0, 0, 0, 1, 1, 0, 1, 0, 0, 1, 0, 0, 0, 1,      1, 0, 0, 0, 1, 0, 0, 0, 1, 0, 1, 1, 0), .dim = c(52l, 2l), .dimnames = list(         null, c("time", "status")), type = "right", class = "surv")), .names = c("conflict",      "peacedur", "peacefail", "totalps", "year_end", "peacedur.year",      "survobj"), row.names = c("1", "2", "3", "4", "5", "6", "7",      "8", "9", "10", "11", "12", "13", "14", "15", "16", "17", "18",      "19", "20", "21", "22", "23", "24", "25", "26", "27", "28", "29",      "30", "31", "32", "33", "34", "35", "36", "37", "38", "39", "40",      "41", "42", "43", "44", "45", "46", "47", "48", "49", "51", "52",          "53"), class = "data.frame") 

the creation of survival object

library(survival) library(rms) duration$survobj <- with(duration, surv(peacedur, peacefail==1)) 

the fit + converstion npsurv

km.duration.totalps <- survfit(survobj ~ totalps, data = duration, conf.type = "log-log") class(km.duration.totalps) <- c(class(km.duration.totalps), "npsurv") 

the plot:

survplot(km.duration.totalps,           xlab="duration in months", ylab="survival prob",          conf="none",          label.curves = true,                              time.inc=12,          levels.only  = false,                             n.risk=true) 

totalps=0 has 1 event. n.risk doesn't that.

running survplot(km.duration.totalps[-1], ...) works setting n.risk = false

library(survival) library(rms) duration$survobj <- with(duration, surv(peacedur, peacefail==1))  km.duration.totalps <- survfit(survobj ~ totalps, data = duration, conf.type = "log-log") class(km.duration.totalps) <- c(class(km.duration.totalps), "npsurv")  summary(km.duration.totalps)  # call: survfit(formula = survobj ~ totalps, data = duration, conf.type = "log-log") #  # totalps=0  # time       n.risk      n.event     survival      std.err lower 95% ci upper 95% ci  # 2            1            1            0          nan           na           na  #  # totalps=1  # time n.risk n.event survival std.err lower 95% ci upper 95% ci # 0.1     10       1    0.900  0.0949       0.4730        0.985 # 7.0      9       1    0.800  0.1265       0.4087        0.946 # ...    par(mfrow = c(2,1)) survplot(km.duration.totalps,          xlab="duration in months", ylab="survival prob",          conf="none",          label.curves = true,                              time.inc=12,          levels.only  = false,                             n.risk=false)  survplot(km.duration.totalps[-1],          xlab="duration in months", ylab="survival prob",          conf="none",          label.curves = true,                              time.inc=12,          levels.only  = false,                             n.risk=true) 

enter image description here

digging deeper root of problem, line in rms:::survplot.npsurv

text(tt[-1], yy, nri[-1], cex = cex.n.risk, adj = 1) 

your nri totalps=0 vector of length 1, r attempting plot

text(1, 1, integer(0)) 

try , same error. remedy that, either use solutions above (plotting totalps=0 not interesting anyway since straight line), or can edit source code below, inserting if statement near end of code trick if (length(nri) > 1)

so can use new function full table/plot without errors (i wouldn't, though, cause can see, label covers @ risk table)

survplot2(km.duration.totalps,           xlab="duration in months", ylab="survival prob",           conf="none",           label.curves = true,                               time.inc=12,           levels.only  = false,                              n.risk=true) 

enter image description here

code:

survplot2 <- function (fit, xlim, ylim, xlab, ylab, time.inc, conf = c("bands", "bars", "diffbands", "none"),           add = false, label.curves = true,            abbrev.label = false, levels.only = false, lty, lwd = par("lwd"),            col = 1, col.fill = gray(seq(0.95, 0.75, length = 5)), loglog = false,            fun, n.risk = false, logt = false, dots = false, dotsize = 0.003,            grid = null, srt.n.risk = 0, sep.n.risk = 0.056, adj.n.risk = 1,            y.n.risk, cex.n.risk = 0.6, pr = false, ...) {   conf <- match.arg(conf)   polyg <- ordgridfun(grid = false)$polygon   conf.int <- fit$conf.int   if (!length(conf.int) | conf == "none")      conf.int <- 0   opar <- par(c("mar", "xpd"))   on.exit(par(opar))   fit.orig <- fit   units <- fit$units   if (!length(units))      units <- "day"   maxtime <- fit$maxtime   if (!length(maxtime))      maxtime <- max(fit$time)   mintime <- min(fit$time, 0)   pret <- pretty(c(mintime, maxtime))   maxtime <- max(pret)   mintime <- min(pret)   if (missing(time.inc)) {     time.inc <- switch(units, day = 30, month = 1, year = 1,                         (maxtime - mintime)/10)     if (time.inc > maxtime)        time.inc <- (maxtime - mintime)/10   }   if (n.risk && !length(fit$n.risk)) {     n.risk <- false     warning("fit not have number @ risk\nis parametric model\nn.risk set f")   }   trans <- loglog | !missing(fun)   if (missing(ylab)) {     if (loglog)        ylab <- "log(-log survival probability)"     else if (trans)        ylab <- ""     else ylab <- "survival probability"   }   if (loglog)      fun <- function(w) logb(-logb(ifelse(w == 0 | w == 1,                                           na, w)))   else if (!trans)      fun <- function(w) w   if (missing(xlab)) {     if (logt)        xlab <- paste("log survival time in ", units, "s",                      sep = "")     else xlab <- if (units == " ")        ""     else paste(units, "s", sep = "")   }   if (missing(xlim))      xlim <- if (logt)        logb(c(maxtime/100, maxtime))   else c(mintime, maxtime)   if (trans) {     fit$surv <- fun(fit$surv)     fit$surv[is.infinite(fit$surv)] <- na     if (conf.int > 0) {       fit$lower <- fun(fit$lower)       fit$upper <- fun(fit$upper)       fit$lower[is.infinite(fit$lower)] <- na       fit$upper[is.infinite(fit$upper)] <- na       if (missing(ylim))          ylim <- range(c(fit$lower, fit$upper), na.rm = true)     }     else if (missing(ylim))        ylim <- range(fit$surv, na.rm = true)   }   else if (missing(ylim))      ylim <- c(0, 1)   if (length(grid)) {     dots <- false     if (is.logical(grid))        grid <- if (grid)          gray(0.8)     else null   }   if (logt | trans) {     dots <- false     grid <- null   }   olev <- slev <- names(fit$strata)   if (levels.only)      slev <- gsub(".*=", "", slev)   sleva <- if (abbrev.label)      abbreviate(slev)   else slev   ns <- length(slev)   slevp <- ns > 0   labelc <- is.list(label.curves) || label.curves   if (!slevp)      labelc <- false   ns <- max(ns, 1)   y <- 1:ns   stemp <- if (ns == 1)      rep(1, length(fit$time))   else rep(1:ns, fit$strata)   if (n.risk | (conf.int > 0 & conf == "bars")) {     stime <- seq(mintime, maxtime, time.inc)     v <- summary(fit, times = stime, print.it = false)     vs <- if (ns > 1)        as.character(v$strata)   }   xd <- xlim[2] - xlim[1]   yd <- ylim[2] - ylim[1]   if (n.risk && !add) {     mar <- opar$mar     if (mar[4] < 4) {       mar[4] <- mar[4] + 2       par(mar = mar)     }   }   lty <- if (missing(lty))      seq(ns + 1)[-2]   else rep(lty, length = ns)   lwd <- rep(lwd, length = ns)   col <- rep(col, length = ns)   if (labelc || conf == "bands")      curves <- vector("list", ns)   tim <- srv <- list()   par(xpd = na)   (i in 1:ns) {     st <- stemp ==     time <- fit$time[st]     surv <- fit$surv[st]     if (logt)        time <- logb(time)     s <- !is.na(time) & (time >= xlim[1])     if (i == 1 & !add) {       plot(time, surv, xlab = xlab, xlim = xlim, ylab = ylab,             ylim = ylim, type = "n", axes = false)       mgp.axis(1, @ = if (logt)          pretty(xlim)         else seq(xlim[1], max(pretty(xlim)), time.inc), labels = true)       mgp.axis(2, @ = pretty(ylim))       if (dots || length(grid)) {         xlm <- pretty(xlim)         xlm <- c(xlm[1], xlm[length(xlm)])         xp <- seq(xlm[1], xlm[2], = time.inc)         yd <- ylim[2] - ylim[1]         if (yd <= 0.1)            yi <- 0.01         else if (yd <= 0.2)            yi <- 0.025         else if (yd <= 0.4)            yi <- 0.05         else yi <- 0.1         yp <- seq(ylim[2], ylim[1] + if (n.risk && missing(y.n.risk))            yi           else 0, = -yi)         if (dots)            (tt in xp) symbols(rep(tt, length(yp)),                                   yp, circles = rep(dotsize, length(yp)), inches = dotsize,                                   add = true)         else abline(h = yp, v = xp, col = grid, xpd = false)       }     }     tim <- time[s]     srv <- surv[s]     if (conf.int > 0 && conf == "bands") {       blower <- fit$lower[st][s]       bupper <- fit$upper[st][s]     }     if (max(tim) > xlim[2]) {       srvl <- srv[tim <= xlim[2] + 1e-06]       s.last <- srvl[length(srvl)]       k <- tim < xlim[2]       tim <- c(tim[k], xlim[2])       srv <- c(srv[k], s.last)       if (conf.int > 0 && conf == "bands") {         low.last <- blower[time <= xlim[2] + 1e-06]         low.last <- low.last[length(low.last)]         up.last <- bupper[time <= xlim[2] + 1e-06]         up.last <- up.last[length(up.last)]         blower <- c(blower[k], low.last)         bupper <- c(bupper[k], up.last)       }     }     if (logt) {       if (conf %nin% c("bands", "diffbands"))          lines(tim, srv, type = "s", lty = lty[i], col = col[i],                lwd = lwd[i])       if (labelc || conf %in% c("bands", "diffbands"))          curves[[i]] <- list(tim, srv)     }     else {       xxx <- c(mintime, tim)       yyy <- c(fun(1), srv)       if (conf %nin% c("bands", "diffbands"))          lines(xxx, yyy, type = "s", lty = lty[i], col = col[i],                lwd = lwd[i])       if (labelc || conf %in% c("bands", "diffbands"))          curves[[i]] <- list(xxx, yyy)     }     if (pr) {       zest <- rbind(time[s], surv[s])       dimnames(zest) <- list(c("time", "survival"), rep("",                                                          sum(s)))       if (slevp)          cat("\nestimates ", slev[i], "\n\n")       print(zest, digits = 3)     }     if (conf.int > 0) {       if (conf == "bands") {         if (logt)            polyg(x = c(tim, max(tim), rev(tim)), y = c(blower,                                                        rev(bupper), max(bupper)), col = col.fill[i],                  type = "s")         else polyg(x = c(max(tim), tim, rev(c(tim, max(tim)))),                     y = c(fun(1), blower, rev(c(fun(1), bupper))),                     col = col.fill[i], type = "s")       }       else if (conf == "diffbands")          survdiffplot(fit.orig, conf = conf, fun = fun)       else {         j <- if (ns == 1)            true         else vs == olev[i]         tt <- v$time[j]         ss <- v$surv[j]         lower <- v$lower[j]         upper <- v$upper[j]         if (logt)            tt <- logb(ifelse(tt == 0, na, tt))         tt <- tt + xd * (i - 1) * 0.01         errbar(tt, ss, upper, lower, add = true, lty = lty[i],                 col = col[i])       }     }     if (n.risk) {       j <- if (ns == 1)          true       else vs == olev[i]       tt <- v$time[j]       nrisk <- v$n.risk[j]       tt[1] <- xlim[1]       if (missing(y.n.risk))          y.n.risk <- ylim[1]       yy <- y.n.risk + yd * (ns - i) * sep.n.risk       nri <- nrisk       nri[tt > xlim[2]] <- na       text(tt[1], yy, nri[1], cex = cex.n.risk, adj = adj.n.risk,             srt = srt.n.risk)         ## add condition here       if (length(nri) > 1)         text(tt[-1], yy, nri[-1], cex = cex.n.risk, adj = 1)         if (slevp)          text(xlim[2] + xd * 0.025, yy, adj = 0, sleva[i],               cex = cex.n.risk)     }   }   if (conf %in% c("bands", "diffbands"))      (i in 1:ns) lines(curves[[i]][[1]], curves[[i]][[2]],                            lty = lty[i], lwd = lwd[i], col = col[i], type = "s")   if (labelc)      labcurve(curves, sleva, type = "s", lty = lty, lwd = lwd,               opts = label.curves, col. = col)   invisible(slev) } 

Comments

Popular posts from this blog

javascript - RequestAnimationFrame not working when exiting fullscreen switching space on Safari -

Python ctypes access violation with const pointer arguments -