I have some survival data from different cohorts and can fit Kaplan Meier curves and cox proportional hazard models. I want to be able to simulate data that resembles the data found in say the KM curves. What I would like to know is there any way to simulate data from the KM curves or the cox proportional hazard model?



Remember that the KM curves are estimates of the survival function, ie $\hat S_{km}(t)$. As such, you can simulate draws from this distribution by using the inverse survival function, i.e.

$u_i \sim $ uniform(0,1)

$T_i = S_{km}^{-1}(u_i)$

Then $T_i$ will be distributed according the distribution in the KM curve. Note that the KM curve will be a discrete distribution, while the original data generating process may be continuous. Also, this does not take into account the uncertainty in the KM curves.

Here's how one could do this in R:

y <- rexp(100)
isCen <- y > 2
y[isCen] <- 2
fit <- survfit(Surv(y, !isCen) ~ 1 )
quantile(fit, probs = runif(1) )$quantile

Very technically speaking, a Cox-PH model does not necessarily include an estimated baseline survival distribution, which is required to sample from. However, this can be easily remedied by using the KM curves as the baseline survival estimate. Of course, the following method can be used for any proportional hazards model and baseline survival estimate.

From here, we use the relation that if $\eta_i = X^T \beta$ (i.e. the linear predictor for subject i), the estimated proportional hazards survival curve will be $S(t | \eta_i)$ = $S_o(t)^{\exp(\eta_i)}$ where $S_o(t)$ is the baseline survival curve.

Being extra clever, we can then see that

$S^{-1}(u | \eta_i) = S_o^{-1}(u^{\exp(-\eta_i)})$

So given $\hat S_o^{-1}$ and $\eta_i$, we can take draws via

$u_i \sim$ uniform(0,1)

$T_i = \hat S_o^{-1}(u^{\exp(-\eta_i)})$

  • 1
    $\begingroup$ Awesome, that makes sense. What about simulating from the cox model, is that doable? $\endgroup$ – RustyStatistician Mar 4 '16 at 2:52
  • $\begingroup$ @RustyStatistician: updated it to include the cox model. $\endgroup$ – Cliff AB Mar 4 '16 at 5:04
  • $\begingroup$ Could we talk in a private room or offline? I would like to pick your brain about this question some more for my actual application. $\endgroup$ – RustyStatistician Mar 4 '16 at 17:09
  • $\begingroup$ @RustyStatistician:sure, except I'm not really sure how to the private chat room think on CV. Can you invite me? $\endgroup$ – Cliff AB Mar 4 '16 at 17:22
  • $\begingroup$ How do you calculate $S^{-1}_{km}$? $\endgroup$ – RustyStatistician Mar 4 '16 at 19:44

@CliffAB has provided a very good answer (+1). Let me add a couple of additional possibilities:

  1. You can try various parametric survival models to see if you can find a distribution that is a reasonable fit. Then you can sample from that distribution using standard techniques. (For a simple example of generating censored survival data from a Weibull distribution, see my answer here: How to simulate censored data.)
  2. You can also bootstrap your data, which will take your sample distributions as estimates of the population distributions. This will be more informative if your samples are larger and better behaved. (The results here should be essentially the same as @CliffAB's strategy. That is, this gives you another way of understanding the strategy he suggested.)
  • $\begingroup$ For the record, I actually advise your suggestion (1), provided one can find a parametric distribution that agrees well with the KM curves over directly sampling from the KM curves. This is mainly due to the issues of the KM curves being discrete, where as we often think of the generating process being continuous. $\endgroup$ – Cliff AB Mar 4 '16 at 20:40
  • $\begingroup$ Smooth fits are likely to be just as good or better than K-M estimates if model assumptions are reasonably well satisfied. I wouldn't worry so much about imitating a step function. When only given KM curves (say from published papers) I digitize the curves and fit a regression spline to the log cumulative hazard (- log KM estimates) for simulation. $\endgroup$ – Frank Harrell Mar 4 '16 at 20:45
  • $\begingroup$ @CliffAB, did you suggest that? I'm not seeing it. I can delete this if it is just duplicative. $\endgroup$ – gung - Reinstate Monica Mar 4 '16 at 20:52
  • $\begingroup$ @FrankHarrell, I was thinking along those lines, but it isn't fully worked out in my head. (I was thinking of a lowess fit to the survival curve & using that analogously to CliffAB's answer.) Why not write that up as an official answer? $\endgroup$ – gung - Reinstate Monica Mar 4 '16 at 20:55
  • 1
    $\begingroup$ Frank Harrell's suggestion is also quite good as well. On the other hand, fitting a lowess to the survival curve can be a little complicated; you need to assure that it is constrained to [0,1], non-decreasing, etc. These types of issues also exist for cumulative hazard fits, but more work (and available code) exists for smoothly fitting the cumulative hazards. $\endgroup$ – Cliff AB Mar 4 '16 at 21:05

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.