Un modelo más potente para conseguirlo, es una metodología relativamente nueva que ha sido propuesta por Wulfsohn y Tsiatis, 1997 (véase también una importante contribución de Henderson, Diggle y Dobson, 2000 ) es el Modelización conjunta de datos longitudinales y de supervivencia .
El modelo consta de dos partes, el submodelo longitudinal:
$$ y_i(t_{ij}) = W_i(t_{ij}) + e_{ij} ~~~, e \sim \text{N}(0,\sigma_e^2)$$
con $ W_i(t_{ij}) = x_i'(t_{ij}) \beta + z_i'(t_{ij}) b_i + u_i \delta $
y el submodelo de supervivencia:
$$ h_i(T_i) = h_0(T_i) \exp(\alpha W_i(T_i) + \phi v_i) $$
con $W_i(T_i) = \beta_{0i} + \beta_{1i} T_i + \delta u_i$
y la estimación se realiza mediante la verosimilitud conjunta. La motivación de este método es exactamente la que usted necesita, para ajustar modelos de supervivencia con covariables dependientes del tiempo.
Hay dos paquetes en CRAN, el JM y el joineR mi preferencia es la primera, principalmente por su documentación. Puedes investigar más a fondo para ver las modificaciones u opciones que puedas necesitar para tu problema.
