timedep function for more details.
Joint Frailty model
Fit a joint either with gamma or log-normal frailty model for recurrent and terminal events using a penalized likelihood estimation on the hazard function or a parametric estimation. Right-censored data is allowed. Left-truncated data and stratified analysis are not possible. Joint frailty models allow studying, jointly, survival processes of recurrent and terminal events, by considering the terminal event as an informative censoring.
There is two kinds of joint frailty models that can be fitted with frailtyPenal :
- The first one (Rondeau et al. 2007) includes a common frailty term to the individuals $(\omega_i)$ for the two rates which will take into account the heterogeneity in the data, associated with unobserved covariates. The frailty term acts differently for the two rates ( $\omega_i$ for the recurrent rate and $\omega_i^{\alpha}$ for the death rate). The covariates could be different for the recurrent rate and death rate.
For the $j^{th}$ recurrence $(j=1,...,n_i)$ and the $i^{th}$ subject $(i=1,...,G)$, the joint gamma frailty model for recurrent event hazard function $r_{ij}(.)$ and death rate $\lambda_i(.)$ is :
$$\left{
\begin{array}{ll}
r_{ij}(t|\omega_i)=\omega_ir_0(t)\exp(\bold{\beta_1^{'}Z_i(t)}) & \mbox{(Recurrent)} \
\lambda_i(t|\omega_i)=\omega_i^{\alpha}\lambda_0(t)\exp(\bold{\beta_2^{'}Z_i(t)}) & \mbox{(Death)} \
\end{array}
\right.$$
where $r_0(t)$ (resp. $\lambda_0(t)$) is the recurrent (resp. terminal) event baseline hazard function, $\bold{\beta_1}$ (resp. $\bold{\beta_2}$) the regression coefficient vector, $\bold{Z_i(t)}$ the covariate vector. The random effects of frailties $\omega_i\sim\bold{\Gamma}(\frac{1}{\theta},\frac{1}{\theta})$ and are iid.
The joint log-normal frailty model will be :
$$\left{
\begin{array}{ll}
r_{ij}(t|\eta_i)=r_0(t)\exp(\eta_i+\bold{\beta_1^{'}Z_i(t)}) & \mbox{(Recurrent)} \
\lambda_i(t|\eta_i)=\lambda_0(t)\exp(\alpha \eta_i+\bold{\beta_2^{'}Z_i(t)}) & \mbox{(Death)} \
\end{array}
\right.$$
where $$\eta_i\sim N(0,\sigma^2)$$
- The second one (Rondeau et al. 2011) is quite similar but the frailty term is common to the individuals from a same group. This model is useful for joint modelling two clustered survival outcomes. In this case, j is for the subject and i for the cluster.
$$\left{
\begin{array}{ll}
r_{ij}(t|u_i)=u_ir_0(t)\exp(\bold{\beta_1^{'}Z_{ij}(t)}) & \mbox{(Time to event)} \
\lambda_{ij}(t|u_i)=u_i^{\alpha}\lambda_0(t)\exp(\bold{\beta_2^{'}Z_{ij}(t)}) & \mbox{(Death)} \
\end{array}
\right.$$
In case of a log-normal distribution of the frailties, we will have :
$$\left{
\begin{array}{ll}
r_{ij}(t|v_i)=r_0(t)\exp(v_i+\bold{\beta_1^{'}Z_{ij}(t)}) & \mbox{(Time to event)} \
\lambda_{ij}(t|v_i)=\lambda_0(t)\exp(\alpha v_i+\bold{\beta_2^{'}Z_{ij}(t)}) & \mbox{(Death)} \
\end{array}
\right.$$
where $$v_i\sim N(0,\sigma^2)$$
This joint frailty model can also be applied to clustered recurrent events and a terminal event (example on "readmission" data below).
From now on, you can also consider time-varying effects covariates in your model, see timedep function for more details.
Nested Frailty model
Fit a nested frailty model using a Penalized Likelihood on the hazard function or using a parametric estimation. Nested frailty models allow survival studies for hierarchically clustered data by including two iid gamma random effects.
Left truncated and censored data are allowed. Stratification analysis is allowed(maximum of strata=2).
The hazard function conditional on the two frailties $v_i$ and $w_{ij}$ for the $k^{th}$ individual of the $j^{th}$ subgroup of the $i^{th}$ group is :
$$\left{
\begin{array}{ll}
\lambda_{ijk}(t|v_i,w_{ij})=v_iw_{ij}\lambda_0(t)exp(\bold{\beta^{'}X_{ijk}}) \
v_i\sim\Gamma\left(\frac{1}{\alpha},\frac{1}{\alpha}\right) \hspace{0.05cm}i.i.d. \hspace{0.2cm} \bold{E}(v_i)=1 \hspace{0.2cm}\bold{Var}(v_i)=\alpha \
w_{ij}\sim\Gamma\left(\frac{1}{\eta},\frac{1}{\eta}\right)\hspace{0.05cm}i.i.d. \hspace{0.2cm} \bold{E}(w_{ij})=1 \hspace{0.2cm} \bold{Var}(w_{ij})=\eta
\end{array}
\right.$$
where $\lambda_0(t)$ is the baseline hazard function, $X_{ijk}$
denotes the covariate vector and $\beta$ the corresponding vector of regression parameters.frailtyPenal(formula, formula.terminalEvent, data, Frailty = FALSE,
joint = FALSE, recurrentAG = FALSE, cross.validation =
FALSE, n.knots, kappa1, kappa2, maxit = 350, hazard =
"Splines", nb.int1, nb.int2, RandDist = "Gamma",
betaknots = 1, betaorder = 3, B)SurvIC, subcluster, terminal, num.id, timedep### COX proportional hazard model (SHARED without frailties) ###
### estimated with penalized likelihood ###
data(kidney)
frailtyPenal(Surv(time,status)~sex+age,
n.knots=12,kappa1=10000,data=kidney,Frailty=FALSE)
### Shared Frailty model ###
frailtyPenal(Surv(time,status)~cluster(id)+sex+age,
n.knots=12,kappa1=10000,data=kidney,Frailty=TRUE)
# with an initialisation of regression coefficients
frailtyPenal(Surv(time,status)~cluster(id)+sex+age,
n.knots=12,kappa1=10000,data=kidney,Frailty=TRUE,B=c(-1.44,0))
# with truncated data
data(dataNested)
frailtyPenal(Surv(t1,t2,event) ~ cluster(group),
data=dataNested,n.knots=10,kappa1=10000,
Frailty=TRUE,cross.validation=TRUE,recurrentAG=F)
# stratified analysis
data(readmission)
frailtyPenal(Surv(time,event)~cluster(id)+dukes+strata(sex),
n.knots=10,kappa1=10000,kappa2=10000,data=readmission,
Frailty=TRUE)
# recurrentAG=TRUE
frailtyPenal(Surv(t.start,t.stop,event)~cluster(id)+sex+dukes+
charlson,data=readmission,Frailty=TRUE,
n.knots=6,kappa1=100000,recurrentAG=TRUE)
# cross.validation=TRUE
frailtyPenal(Surv(t.start,t.stop,event)~cluster(id)+sex+dukes+
charlson,data=readmission, Frailty=TRUE,n.knots=6,
kappa1=5000,recurrentAG=TRUE,cross.validation=TRUE)
# log-normal distribution
frailtyPenal(Surv(t.start,t.stop,event)~sex+dukes+charlson+cluster(id),
data=readmission, Frailty=TRUE,n.knots=6,kappa1=5000,
recurrentAG=TRUE,RandDist="LogN")
### Joint Frailty model (recurrent and terminal events) ###
data(readmission)
# Gap-time
modJoint.gap <- frailtyPenal(Surv(time,event)~cluster(id)+sex+dukes+charlson+
terminal(death),formula.terminalEvent=~sex+dukes+charlson,
data=readmission,n.knots=14,kappa1=9.55e+9,kappa2=1.41e+12,
Frailty=TRUE,joint=TRUE,recurrentAG=FALSE)
# Calendar time
modJoint.calendar <- frailtyPenal(Surv(t.start,t.stop,event)~cluster(id)+sex
+dukes+charlson+terminal(death),formula.terminalEvent=~sex
+dukes+charlson,data=readmission,n.knots=10,kappa1=9.55e9,
kappa2=1.41e12,Frailty=TRUE,joint=TRUE,recurrentAG=TRUE)
# log-normal distribution
modJoint.log <- frailtyPenal(Surv(t.start,t.stop,event)~cluster(id)+sex
+dukes+charlson+terminal(death),formula.terminalEvent=~sex
+dukes+charlson,data=readmission,n.knots=10,kappa1=9.55e9,
kappa2=1.41e12,Frailty=TRUE,joint=TRUE,recurrentAG=TRUE,
RandDist="LogN")
### Joint Frailty model for clustered data ###
# here is a generated cluster (31 clusters of 13 subjects)
readmission <- transform(readmission,group=id%%31+1)
joi.clus <- frailtyPenal(Surv(t.start,t.stop,event)~cluster(group)+
num.id(id)+dukes+charlson+sex+chemo+terminal(death),
formula.terminalEvent=~dukes+charlson+sex+chemo,
data=readmission,Frailty=TRUE,joint=TRUE,recurrentAG=TRUE,
n.knots=10,kappa1=2.11e+08,kappa2=9.53e+11)
### Nested Frailty model ###
data(dataNested)
modClu <- frailtyPenal(Surv(t1,t2,event)~cluster(group)+
subcluster(subgroup)+cov1+cov2,Frailty=TRUE,
data=dataNested,n.knots=8,kappa1=50000)
modClu.str <- frailtyPenal(Surv(t1,t2,event)~cluster(group)+
subcluster(subgroup)+cov1+strata(cov2),Frailty=TRUE,
data=dataNested,n.knots=8,kappa1=50000,kappa2=50000)Run the code above in your browser using DataLab