Hi all
Please could you help me to figure out how to simulate competing risks survival data?
I am aware that I can achieve this using survsim command but I would like to generate this from first principles.
The idea is to simulate a HR for treatment for the event of interest and then have a HR for treatment for the competing risk. I want to understand how much the HR from the cox regression is biased compared to the HR from stcrreg. In the code, I tried to assume minimal censoring but this can be modified later.
I have tried to do this with the following code but it does not seem correct. Could you help me figure out the correct code?
From my understanding this is a key paper explaining the algorithm for the competing risks simulation which I am struggling to understand
Beyersmann, Jan, Aurelien Latouche, Anika Buchholz, and Martin Schumacher. "Simulating competing risks data in survival analysis." Statistics in medicine 28, no. 6 (2009): 956-971.
Code:
Thanks
A
Please could you help me to figure out how to simulate competing risks survival data?
I am aware that I can achieve this using survsim command but I would like to generate this from first principles.
The idea is to simulate a HR for treatment for the event of interest and then have a HR for treatment for the competing risk. I want to understand how much the HR from the cox regression is biased compared to the HR from stcrreg. In the code, I tried to assume minimal censoring but this can be modified later.
I have tried to do this with the following code but it does not seem correct. Could you help me figure out the correct code?
From my understanding this is a key paper explaining the algorithm for the competing risks simulation which I am struggling to understand
Beyersmann, Jan, Aurelien Latouche, Anika Buchholz, and Martin Schumacher. "Simulating competing risks data in survival analysis." Statistics in medicine 28, no. 6 (2009): 956-971.
Code:
Code:
clear set obs 3000 gen treatment = runiform()>0.5 gen double u = runiform() // 2 represents the competing risk and 1 represents the event of interest gen risk_events = runiform()>0.95 // competing risk, 50/50 chance recode risk_events 0=1 1=2 tab risk_events // the hazard for the event of interest local mediantime = 12 local b_cons = -ln(0.5)/`mediantime' display `b_cons' gen xb1 = ln(`b_cons') + ln(0.70)*treatment gen lambda1 = exp(xb1) // the hazard for the competing risk local mediantime = 7 local b_cons = -ln(0.5)/`mediantime' display `b_cons' gen xb2 = ln(`b_cons') + ln(1)*treatment gen lambda2 = exp(xb2) gen t = -log(u)/lambda1 if risk_events==1 replace t = -log(u)/lambda2 if risk_events==2 // censoring time gen cens=runiform(1000,2000) // generate a non-informative censoring indicator gen event = 1 replace event = 0 if t>cens gen t0 = 0 replace t0 = t if event==1 gen state1 = 1 replace state1 = risk_events if event==1 replace state1 = 0 if event==0 tab event state1, col gen t1 = t replace t1 = cens if event==0 stset t1, fail(state1=1) sts graph, by(treatment) stcox treatment stcrreg treatment, compete(state1==2)
A
Comment