How to generate simulation paths for Weibull distribution in R?

Question

I'm trying to generate simulation paths across time in R for survival and hazard functions using the Weibull distribution. Here are the steps I'm taking so far. I need direction because I'm not able to generate simulation paths with the below:

First: Here's an example of a single basic hazard curve generated in R which works fine (code followed by image):

# function for hazard:
weibHaz <- function(x, shape, scale) dweibull(x, shape=shape,
              scale=scale)/pweibull(x, shape=shape, scale=scale,lower.tail=F)
plot hazard:
curve(weibHaz(x, shape=1.5, scale=1/0.03), from=0, to=80, ylab='Hazard', xlab='Time', col='red')

Second: Generate random numbers from Weibull distribution and store them in simulated_data, using same parameters as in the above hazard curve:

set.seed(123) # set a seed for reproducibility
simulated_data <- rweibull(n = 20, shape=1.5, scale=1/0.03) 

Third: plot the simulation results (code followed by image):

# ts() creates a time series object from the simulated data
# frequency = 1 indicates the time series data is not periodic
ts_data <- ts(simulated_data, start = 1, frequency = 1)
plot(ts_data, main = "Simulated Time Series Data", xlab = "Time", ylab = "Data")

The third step above obviously shows incorrect results. I need to see 20 curves over 80 time periods. The below image shows the type of simulation paths I've generated before using other data and using autoregressive time-series models. This is the form of output data and plot I'm looking for for the hazard curve generated with the Weibull distribution. Any recommendations for doing this correctly?

Answer 1

Your simulated_data are just 20 draws from a single fixed Weibull distribution having the specified shape and scale values. They are just samples of individual survival times drawn from that single distribution. This shows the associated survival function:

 plot(1:70,1-pweibull(1:70,shape=1.5, scale=1/0.03),type="l",bty="n",xlab="Time",ylab="Survival probability")

If the Weibull distribution is fixed, then the survival function is also fixed and there would be no variability in survival curves.

If you want to evaluate uncertainty in survival-curve estimates based on modeled data, then you want to draw from distributions of the shape and scale values consistent with the variance-covariance matrix of their modeled estimates, for example with the mvrnorm() function of the R MASS package, and display the set of survival curves in the way indicated above.