sampling from specified distribution

Question

In a paper[1] on generating survival data from a Cox proportional hazard model, they specify following cumulative distribution function:

\begin{equation} F(t^0) = \frac{1}{5}\sum_{x=0}^{4} \exp(-\exp(\beta x) t^0). \end{equation}

with $\beta = 0.5\log(2)$

They generate random numbers from this distribution but don't specify how. I don't see a way to invert the distribution and inverse transform sampling.
How would you sample from this kind of distribution?

[1]: Mackenzie, Todd, and Michal Abrahamowicz. "Marginal and hazard ratio specific random data generation: applications to semi-parametric bootstrapping." Statistics and Computing 12.3 (2002): 245-252.

Answer 1

You can think of this as a mixture of exponential distributions with five components. Each mixture component has a probability of $1/5$, and the parameters of the individual components are $\exp(\beta x), x = 0, \dots, 4$.

Sampling from this distribution can be done in two steps, sample code in R. First, randomly select the mixture component to sample from:

x <- sample(5,1) - 1

and second, generate the r.v. from the appropriate component:

t0 <- rexp(1, exp(beta*x))

which can be made much more efficient:

beta <- 0.5*log(2)
N <- 1000
t0 <- rexp(N, exp(beta*(sample(5, N, replace=TRUE)-1)))

We can compare the results to calculated values of the density function as follows:

x <- seq(0.1,7,0.1)
pdf <- rep(0, length(x))
for (i in 1:length(x)) {
  pdf[i] <- mean(dexp(x[i], exp(beta*(0:4))))
}

hist(t0, freq=FALSE)
lines(pdf~x, lwd=2)

which generates the following plot:

which in turn indicates that we likely haven't messed up anywhere.