I've been working on simulating survival experiments, and came across what seems to be an odd-looking artifact of Kaplan-Meier analysis which I have been struggling to explain in a concise manner mathematically. I've managed to convince myself that it's not an artifact of my simulation process itself. Here's an overview of what I'm doing:
- Generate survival times for a population of size n, randomly. Random values on (0,1) are drawn from a uniform distribution, and the survival times are then taken based on those values from a logistic distribution with a given median.
- Iterate over this population to simulate a survival experiment. Starting with the n members of the population, on each "day" (iteration) they are "scored" and the ones that are set to "die" on that day are recorded as such and not considered in future iterations.
- Continue iterating until the remaining population has 0 members.
- Repeat steps 1-3 10,000 times.
- Repeat steps 1-4 for a range of starting population sizes, all based on the same generating distribution.
- Fit each individual "experiment trial" with Kaplan-Meier (R survival package), and retrieve the median, mean, and 95% quantile survival estimates
- Summarize across the 10k "experiment trial" results for a given population size (i.e. s.d. of medians)
This is implemented in R, and I've run it both with and without ensuring reproducibility of random numbers, and with two different RNGs and saw the same patterns.
I see an odd oscillatory pattern in the summarized quantile survival estimates (median or 95% survival). I do not see such a pattern with the summaries of mean survival.
It gets more pronounced with looking at single-increment increases in population. Here's the same thing, but with 500 "experiment trials" (instead of 10k) and individual increments in population size from 5 - 300:
At that scale, the oscillation in the s.d. of the medians isn't really obvious, but both the mean and the s.d. of the 95% survival just looks wacky. The period of the repeating pattern in the mean of the 95% values across the increases in population seems to be 20. 20 is the median of the generating distribution I used, and I'm not sure if that's related.
I also see a pattern similar to the alternating one in the s.d. of the medians, in the medians of the medians when I show all the medians from each "experiment trial" as bars in a barplot. Since these are alternating, and the population size sequence is 5, 10, 15 etc, presumably that particular behavior is just because I'm going between even and odd population sizes, though I haven't explored that particular question in much depth. That alternating pattern also appears when I calculate s.d. or SEM or mean-squared error on the same data.
I don't have a reproducible example yet, but if anyone else out there has done this kind of large-scale simulation of survival and seen a similar pattern please let me know. I wasn't able to find anything relevant in the literature, based on the search terms I came up with at least.