Mice: Weibull regression in censored survival analysis

dellaportas:smith:93 analyse data from grieve:87 on photocarcinogenicity in four groups, each containing 20 mice, who have recorded a survival time and whether they died or were censored at that time. A portion of the data, giving survival times in weeks, are shown below. A * indicates censoring.

The survival distribution is assumed to be Weibull. That is

where is the failure time of an individual with covariate vector and is a vector of unknown regression coefficients. This leads to a baseline hazard function of the form

Setting gives the parameterisation

For censored observations, the survival distribution is a truncated Weibull, with lower bound corresponding to the censoring time. The regression coefficients were assumed a priori to follow independent Normal distributions with zero mean and ``vague'' precision 0.0001. The shape parameter r for the survival distribution was given a Gamma(1, 0.0001) prior, which is slowly decreasing on the positive real line.

Median survival for individuals with covariate vector is given by .

The appropriate graph is shown in Figure 18, using an undirected dashed line to represent a logical range constraint. The BUGS code is given on the next page.

 
Figure 18:   Graphical model for mice example

Model specification for mice example

model mice;

const
  N = 80,    # number of individuals
  M = 4;     # number of treatment groups
var
  t[N],                    # failure time for each mouse
  t.cen[N],                # censoring time for each mouse
  group[N],                # treatment group for each mouse
  mu[N], r,                # Weibull parameters  
  beta[M],                 # log relative risk parameters
  median[M],               # median survival for each group
  irr.control,veh.control, # treatment contrasts
  test.sub,pos.control;  

data t, t.cen, group in "mice.dat";
inits in "mice.in";
{
  for(i in 1:N) {                          # t.cen[i] = 0  
        t[i] ~ dweib(r,mu[i])I(t.cen[i],); # if mouse i fails  
                                           
        mu[i] <- exp(beta[group[i]]);      # relative risk model
  }
  for(j in 1:M) {
        beta[j] ~ dnorm(0.0, 0.0001);      # prior
        median[j] <- pow(log(2) *          # median survival
                       exp(-beta[j]), 1/r);  
  }
  r ~ dgamma(1.0,0.0001);      # slowly decreasing on +ve reals

  irr.control <- beta[1];                  # change 
  veh.control <- beta[2]-beta[1];          # parameterisation           
  test.sub <- beta[3]-beta[1];
  pos.control <- beta[4]-beta[1];
}

We note a number of tricks in setting up this model. First, individuals who are censored are given a missing value in the vector of failure times t, whilst individuals who fail are given a zero in the censoring time vector t.cen (see data file listing below). The truncated Weibull is modelled using I(t.cen[i],) to set a lower bound. Second, we set a parameter beta[j] for each treatment group j, and create a single covariate group taking values 1, 2, 3 or 4 in the data file. Nested subscripts i.e. beta[group[i]] are used to select the required beta[j] to appear in the linear predictor according to the value of group[i] for the ith individual. The contrasts beta[j] with group 1 (the irradiated control) are calculated at the end. Alternatively, we could have included a grand mean term in the relative risk model and constrained beta[1] to be zero.

Data in rectangular format

 12  0 1
 17  0 1
 21  0 1
 25  0 1
 11  0 1
 26  0 1
 .. .. .
 .. .. .
 35  0 1
 NA 40 1
 31  0 1
 36  0 1
 32  0 2
 27  0 2
 23  0 2
 12  0 2
 18  0 2
 NA 40 2
 .. .. .
 .. .. .

Column 1 refers to the failure times t for each mouse, column 2 represents the corresponding censoring times t.cen, and column 3 indicates which treatment group the mouse belongs to.

Analysis

A simple BUGS run took 22 seconds for 1000 iterations after a 500 iteration burn-in. The output was as follows:

These should be compared to the plots shown by dellaportas:smith:93