how I can simulate survival data with shared frailties?

[Andrew Xavier]

Hi All,

I am trying to understand a bit more about the mechanics of shared frailty survival models and I would like to know how I can simulate survival data with shared frailties.

Let's say that for each subject in the study, I have multiple time to event values which means that there the time to event values are dependent.

Below is an example of a STATA dataset and code that can be used as a starting point.


Ideally, I would like to understand how I can simulate data with following parameters

. streg female, dist(exponential) shared(patient) frailty(gamma)

Gamma shared frailty                                Number of obs     =     76
Group variable: patient                             Number of groups  =     38
                                                    Obs per group:   
No. of subjects =    76                                           min =      2
No. of failures =    58                                           avg =      2
Time at risk    = 7,424                                           max =      2
                                                    LR chi2(1)        =  12.07
Log likelihood = -99.658307                         Prob > chi2       = 0.0005

------------------------------------------------------------------------------
          _t | Haz. ratio   Std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
      female |   .2402041   .0959946    -3.57   0.000     .1097509    .5257175
       _cons |   .0320384   .0115257    -9.56   0.000      .015829    .0648469
-------------+----------------------------------------------------------------
    /lntheta |  -1.081838   .4993981    -2.17   0.030     -2.06064   -.1030357
-------------+----------------------------------------------------------------
       theta |   .3389719   .1692819                      .1273724    .9020947
------------------------------------------------------------------------------

webuse catheter, clear

stset time, fail(infect)

streg female, dist(exponential) shared(patient) frailty(gamma)
stcox female, shared(patient) 

stcox, shared(patient) estimate
predict nu_cox, effects

streg, dist(exponential) shared(patient) frailty(gamma)
predict nu_alpha1, alpha1

Many thanks

Andrew

[Andrew Xavier]

This is an example provided by STATA:

**** begin example

drop _all
*** simulate weibull data (PH)
*** no censoring
set seed 1357

*** number of obs
local size 2000

*** set parameters
** shape: p
** scale: lambda = exp(xb)
local p = 1.7
local b_age = .15
local b_drug = -2
local b_cons = -6

set obs `size'


** create independent variables

gen lnage = .1*rnormal()+4
gen age = exp(lnage)
gen drug = runiform()<.5
gen group=mod(_n,30)

* theta:4
gen alpha= rgamma(0.25,4)

bysort group: replace alpha=alpha[1]

*** linear form
gen xb = `b_drug'*drug + `b_age'*age + `b_cons'
gen lambda = exp(xb)*alpha

*** dependent variable
gen t = rweibullph(`p', lambda)

**** fit model
stset t
streg drug age, distrib(weibull) frailty(gamma) nohr shared(group)

stcox drug age, shared(group)

*These two are similar. They should become closer with larger observation number and group number, but running the Cox model with shared frailty is time-consuming.

predict a, effects
gen a1=log(alpha)

**** end example