# Probabilistic estimates of Kaplan-Meier curves using S(t) and std.err for each time point

I am working on a survival analysis where I am fitting different parametric models to survival data. Varying the models that I have fit to the data is straightforward, via Cholesky decompositions- but I'd like to generate probabilistic Kaplan-Meier curves using the "survival" and "std.err" parameters shown in the data below (I am working in R).

Call: survfit(formula = Surv(futime, fustat) ~ 1, data = PF)

time n.risk n.event survival std.err lower 95% CI upper 95% CI
0    221       0    1.000  0.0000       1.0000        1.000
1    149      52    0.764  0.0286       0.7104        0.823
2     50      24    0.587  0.0396       0.5147        0.670
3     33       3    0.541  0.0445       0.4608        0.636
4     26       1    0.522  0.0469       0.4377        0.623
5     16       4    0.428  0.0577       0.3288        0.558
6      9       3    0.336  0.0653       0.2300        0.492
7      8       1    0.299  0.0679       0.1917        0.467
8      6       2    0.224  0.0685       0.1233        0.408
9      2       2    0.140  0.0642       0.0571        0.344
10      2       0    0.140  0.0642       0.0571        0.344
11      2       0    0.140  0.0642       0.0571        0.344


One idea that I've had is to use the estimated survival probability and variance for each time point to sample from a normal distribution for each time point. There is an obvious problem with this - I could end up sampling a higher survival probability for time 3 than for the time before it, which is illogical.

If I enforce crude "rules" (i.e. force the next sampled S(t) to be equal to or less than the previous S(t) estimate) I will end up systematically underestimating S(t).

I feel as though I am barking up the wrong tree and I'm pretty stumped. Has anyone come across a solution for this before? I'd be grateful for any insights.

Thanks,

constarc

• Could you please edit your question to say more about why you would "like to vary the Kaplan-Meier data itself" and to say more about just what you mean by that phrase? – EdM Oct 8 '19 at 14:37
• I've just updated my question - hopefully it should be a bit clearer. – constarc Oct 8 '19 at 15:19

One trick would be to sample from the distribution of the hazard function $$h(t)$$ itself rather than from the integrated hazard function that underlies the survival function $$S(t)$$. For a Kaplan-Meier analysis this would only be done at discrete observed event times. At event time $$i$$ the estimated hazard $$\hat h_i$$ is $$d_i/n_i$$, where $$d_i$$ is the number of events at that time and $$n_i$$ is the number still at risk. The variance of $$\hat h_i$$ at time $$i$$ is then $$h_i (1-h_i)/n_i$$. (See the derivation of Statistical considerations on the Wikipedia Kaplan-Meier page.) You would start for example with the 221 cases initially at risk and adjust both the numbers of events and the numbers at risk at each event time accordingly (although I haven't thought enough about this to say how to deal properly with censoring in this situation).