I am trying to generate an inverse Weibull distribution using parameters estimated from survreg in R. By this I mean I would like to, for a given probability (which will be a random number in a small simulation model implemented in MS Excel), return the expected time to failure using my parameters. I understand the general form for the inverse Weibull distribution to be:
X=b[-ln(1-rand())]^(1/a)
where a and b are shape and scale parameters respectively and X is the time to failure I want. My problem is in the interpretation of the intercept and covariate parameters from survreg. I have these parameters, the unit of time is days:
Value Std. Error z p
(Intercept) 7.79 0.2288 34.051 0.000
Group 2 -0.139 0.2335 -0.596 0.551
Log(scale) 0.415 0.0279 14.88 0.000
Scale= 1.51
Weibull distribution
Loglik(model)= -8356.7 Loglik(intercept only)= -8356.9
Chisq = 0.37 on 1 degrees of freedom, p= 0.55
Number of Newton-Raphson Iterations: 4
n=1682 (3 observations deleted due to missing values)
I have read in the help files that the coefficients from R are from the "extreme value distribution" but I'm unsure what this really means and how I get 'back to' the standard scale parameter used directly in the formulae. Using b=7.79 and a=1.51 gives nonsensical answers. I really want to be able to generate a time for both the base group and 'Group 2'. I also should note that I did not perform the analysis myself and cannot interrogate the data further.
This is explained in the manual page, ?survreg
(in the "examples" section).
library(survival)
y <- rweibull(1000, shape=2, scale=5)
r <- survreg(Surv(y)~1, dist="weibull")
a <- 1/r$scale # Approximately 2
b <- exp( coef(r) ) # Approximately 5
y2 <- b * ( -ln( 1-runif(1000) ) ) ^(1/a)
y3 <- rweibull(1000, shape=a, scale=5)
# Check graphically that the distributions are the same
plot(sort(y), sort(y2))
abline(0,1)