A routine showing how to calculate quantiles.

km-quantiles.R

library(survival)
library(Hmisc)

# 1. a weibull random number generator, see http://statwonk.com/weibull.html
rweibull_cens <- function(n, shape, scale) {

  # will happen to see the death time first or censoring?
  rweibull(n, shape = shape, scale = scale) -> a_random_death_time
  rweibull(n, shape = shape, scale = scale) -> a_random_censor_time
  pmin(a_random_censor_time, a_random_death_time) -> observed_time

  (observed_time == a_random_death_time) -> censor
  data.frame(time = observed_time, censor = censor, stringsAsFactors = FALSE)
}

# 2. generate sample data
rweibull_cens(1e3, 1, 500) -> d

# 3. calculate a single median 0.5 quantile (percentile/100)
(survival:::quantile.survfit(
  survfit(Surv(time, censor) ~ 1, data = d),
  probs = 0.5
) -> single_estimate)

# 4. calculate 1,000 of them and visualize a distribution
hist(Hmisc::bootkm(Surv(d$time, d$censor), B = 1e3), breaks = 25, col = "orange",
     xlab = "Time until half convert", main = "Time until half convert (median failure time)")
abline(v = single_estimate, col = "black", lwd = 2, lty = 2)