Weibull survival model with random intercept

We can define a multilevel survival model as follows,

$$h(t_{ij}) = h_{0}(t_{ij}) \exp(\mathbf{X}_{ij}\mathbf{\beta} + u_{i})$$

where

$$h_{0}(t_ij) = \lambda \gamma t_{ij} ^ {\gamma-1}$$

and

$$u_{i} \sim N(0,\sigma^{2})$$

We have our vector of baseline covariates $\mathbf{X}_{ij}$ and associated conditional log hazard ratios $\mathbf{\beta}$ (conditional on the frailty, $u_{i}$). We can load an example simulated dataset from the stmixed package

use http://fmwww.bc.edu/repec/bocode/s/stmixed_example1

which represents a multi-centre trial scenario, with 100 centres and each centre recuiting 60 patients, resulting in 6000 observations. Two covariates were collected, a binary covariate $X_{1}$ (coded 0/1), and a continuous covariate, $X_{2}$, within the range [0,1].

Let’s inspect the first few rows of our dataset,

 list centre stime event x1 x2 in 1/5, noobs

  +-------------------------------------------+
  | centre      stime   event   x1         x2 |
  |-------------------------------------------|
  |      1   .2714585       1    1   .4306063 |
  |      1   .7353575       1    0   .5285968 |
  |      1   .8213144       1    0   .1442798 |
  |      1    .758931       1    1   .8612115 |
  |      1   .9547321       0    0    .501694 |
  +-------------------------------------------+

where centre denotes the centre number for each patient, stime contains our survival time, event our indicator variable for those that died (1) or were censored (0), x1 our binary covariate, and x2 our continuous covariate. We can fit our model as follows

megenreg (stime x1 x2 M1[centre], family(weibull, failure(event)))

Getting starting values:

Fitting full model:

Iteration 0:   log likelihood = -2090.3163  
Iteration 1:   log likelihood = -1967.0645  
Iteration 2:   log likelihood = -1963.2032  
Iteration 3:   log likelihood = -1963.2003  
Iteration 4:   log likelihood = -1963.2003  

Log likelihood = -1963.2003                     Number of obs     =      6,000

------------------------------------------------------------------------------
             |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
ln_lambda1   |
          x1 |   1.008781   .0363988    27.71   0.000     .9374409    1.080122
          x2 |  -.9623624   .0603866   -15.94   0.000    -1.080718   -.8440067
       _cons |   .1779824   .0976107     1.82   0.068     -.013331    .3692958
-------------+----------------------------------------------------------------
ln_gamma1    |
       _cons |   .6896406    .013228    52.13   0.000     .6637142     .715567
-------------+----------------------------------------------------------------
lns1_1       |
       _cons |  -.1084075   .0749276    -1.45   0.148    -.2552628    .0384478
------------------------------------------------------------------------------

which gives us an estimate of $\sigma =$ 0.897 (95% CI: 0.775, 1.039), showing substantial heterogeneity between centres.

comments powered by Disqus