Let say I have a 'kidney catheter' data set. Data are about the recurrence times to infection, at the point of insertion of the catheter, for kidney patients using portable dialysis equipment. Catheters may be removed for reasons other than infection, in which case the observation is censored. Each patient has exactly 2 observations.
If I want to fit a Cox PH model with random effect (here Gamma frailty) using EM algorithm. By using built in survival package R code
coxph(), I can easily do it like this
library(survival) data(kidney) fit<-coxph(Surv(time, status) ~ age + sex + frailty(id, dist='gamma', method='em'), kidney)
But if I want to write a step by step function for EM algorithm, how can I proceed?