Weibull Survival Model in R

Question

If I run a Weibull survival model in R with the code

survreg(Surv(t,delta)~expalatory variables, dist="w")

how do I interpret the output of the model? That is, is the form of the model just $1-\exp(e/\lambda)^k)$ with $\lambda$ the scale and $k$ the shape or does it take a different form?

I have found something which says that the output is of the form $$\exp(-\exp(-\alpha_0-\alpha_1-\ldots)^kx^k),$$ where the $\alpha_i$ are the coefficients of the covariates. If so, the output would give me the following parameters:

$$k=k \quad \mbox{and} \quad \lambda=\frac{1}{\exp(-\alpha_0-\alpha_1\ldots)}.$$

Answer 1

Ok so I'm just going to post an answer here using the R help that DWin described. Using the function rweibull in R gives the usual form of the Weibull distribution, with its cumulative function being:

$$F(x)=1-\exp(-\left ( \frac{x}{b}\right )^a)$$

So we will denote the shape parameter of rweibull by $a$ and the scale parameter of rweibull by $b$.

Now the problem is that the output of survreg gives both shape and scale parameters which are not the same as the shape and scale parameters from rweibull. Let us denote the shape parameter from survreg as $a_s$ and the scale parameter of survreg by $b_s$.

Then, from ?survreg we have that:

survreg's scale = 1/(rweibull shape)

survreg's intercept = log(rweibull scale)

So this gives us that:

$$a=\frac{1}{b_s}\quad \mbox{and} \quad b=\exp(a_s)$$

So if we suppose that we run the function survreg with $n$ covariates, then the output will be:

$\alpha_0,\ldots, \alpha_{n-1}$, the coefficents of the covariates and some scale parameter $k$. The Weibull model given in standard form is then given by:

$$F(x)=1-\exp\left (- \left (\frac{x}{\exp(\alpha_0+\alpha_1+\ldots +\alpha_{n-1})} \right ) ^{\frac{1}{k}}\right )$$

Answer 2

I'd like to add an answer with a code example for further clarity.

What we're essentially after is taking the survreg output model and derive from it the survival function. To avoid the common notation confusion I'll actually go ahead and show the code that does that:

fit <- survreg(Surv(time,status) ~ age, data=stanford2) # this is the survreg output model
survreg.lp <- predict(fit, type = "lp")
survreg.scale <- fit$scale

# this is the survival function!
S_t <- function(t, survreg.scale, survreg.lp){
  shape <- 1/survreg.scale 
  scale <- exp(survreg.lp)
  ans <- 1 - pweibull(t, shape = shape, scale = scale)
}

As mentioned by vkehayas R's pweibull parameterisation is:

$$F(x) = 1-exp(-\left(\frac{x}{b}\right)^a$$

where a is the weibull distribution shape and b is the scale.

We then get that a = 1/fit$scale and b = exp(predict(fit, type = "lp"))

We can verify below that the derived survival function

# next let's verify it's correct:
fit <- survreg(Surv(time,status) ~ age, data=stanford2) # this is the survreg output model

# this is the survival function!
S_t <- function(t, survreg.scale, survreg.lp){
  shape <- 1/survreg.scale 
  scale <- exp(survreg.lp)
  ans <- 1 - pweibull(t, shape = shape, scale = scale)
}

new_dat <- data.frame(age = c(0, seq(min(stanford2$age), max(stanford2$age), length.out = 10)))

pct <- seq(0.01, 0.99, 0.01)

surv_curves <- sapply(pct,
                      function(x) predict(fit, type = "quantile", p = 1 - x,
                                          newdata = new_dat))

matplot(y = pct, t(surv_curves), type = "l")

# you can vary the below subject_i variable to see it works for all of them
subject_i <- 1
single_curve <- surv_curves[subject_i, ]
plot(single_curve, pct, type = "l") # this is we know to be true
  
times <- round(seq(1, max(single_curve), length.out = 100))
lp <- predict(fit, newdata = new_dat, type = "lp")[subject_i]
surv <- sapply(times, function(t) S_t(t, survreg.scale = fit$scale, survreg.lp = lp))
lines(times, surv, col = "red", lty = 2) # this is the new S_t function
# They match!

So, to summarize:

a = 1/fit$scale and b = exp(predict(fit, type = "lp"))

Hope this helps. I know I pulled a few hairs before figuring this out.

Answer 3

The help page for ?Weibull says:

The Weibull distribution with shape parameter a and scale parameter b has density given by

f(x) = (a/b) (x/b)^(a-1) exp(- (x/b)^a)

And then the help page for ?survreg says:

# There are multiple ways to parameterize a Weibull distribution. The survreg 
# function embeds it in a general location-scale family, which is a 
# different parameterization than the rweibull function, and often leads
# to confusion.
#   survreg's scale  =    1/(rweibull shape)
#   survreg's intercept = log(rweibull scale)

And there is an elaboration of how this is handled in the Examples section of ?survreg.distributions