This is pretty complicated, so I hope you already know how to simulate mixed models and survival functions. You have privately comunicated with me that you have my book Simulating Data with SAS, so study the relevant sections in the book to help carry out steps 4) and 5) below.
The model is explained in Section 2.1-2.3. There are actually 12 models in this paper (Table 2). I'll describe the main ideas.
1) Choose a model in Table 2.
2) Choose parameters from the "Joint analysis" column of Table 3: fixed effects (beta) coefficients, covariance matrix (Sigma), and random coefficients (gamma)
3) For the m subjects, generate U, which is an mx2 matrix that contains m MVN(0, Sigma) observations. Since you posted this question to the IML forum, you can use
U = RandNormal(m, {0 0}, Sigma); /* random intercepts and slopes */
The i_th row of U is used for the random intercepts and slopes for the i_th patient.
4) If you aren't using a time-varying survival function (such as the exponential hazard function), use Eqn 2 or the Cox PH model to draw a survival time for the i_th patient.
5) Use the mixed model (Eqn 1) to simulate the correlated longitudinal data for the i_th patient. Be sure to stop when the measurement time s_{ij} is greater than the survival time, and mark that patient as censored.
Good luck!
I do not cover time-varying survival functions in my book, but after you are able to simulate the case of constant hazard you can post your program and maybe someone can help you modify it to handle the time-varying case.
