In CecileProust-Lima/lcmm: Extended Mixed Models Using Latent Classes and Latent Processes

library(lcmm)
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)

 

Background on the model

Joint models are used to analyse simultaneously two related phenomena, the evolution of a variable and the occurence of an event. Joint latent class models (JLCM) consist of a linear mixed model and a proportional hazard model linked by the latent classes. The population is split in several groups, the latent classes, and each class is caracterized by a specific evolution of the dependent variable and an associated risk of event.  

Latent class membership is defined by a discrete random variable $c_{i}$ that equals $g$ if subject $i$ belongs to latent class $g$ ($g$ = 1, ...,$G$). The variable $c_{i}$ is latent; its probability is described using a multinomial logistic model according to covariates $X_{ci}$: $\pi_{ig}= P(c_{i} = g|X_{ci}) = \frac{\exp(\xi_{0g}+X_{ci}\xi_{1g})}{ \sum_{l=1}^{G}\exp(\xi_{0l}+X_{ci}\xi_{1l})}$ where $\xi_{0g}$ is the intercept for class $g$ and $\xi_{1g}$ is the q1-vector of class-specific parameters associated with the q1-vector of time-independent covariates $X_{ci}$. For identifiability, $\xi_{0G} = 0$ and $\xi_{1G} = 0$. When no covariate predicts the latent class membership, this model reduces to a class-specific probability.

For a continuous and Gaussian variable, the trajectories of $Y$ are defined conditionally to the latent class by a linear mixed model. So, conditional on class $g$, the model is defined for subject $i$ at occasion $j$:

  $$Y_{ij}|_{c_{i}=g} = X_{2ij}\beta+X_{3ij}\gamma_{g}+Z_{ij}b_{i}+\epsilon_{ij}$$  

where $X_{2ij}$, $X_{3ij}$ and $Z_{ij}$ are vectors of covariates respectively associated with common fixed effects over classes $\beta$, class-specific fixed effects $\gamma_{g}$ and with individual random effects $b_{i}|_{ci=g}$ called $b_{ig}$ whose distributions are now class-specific. $X_{2}$ and $X_{3}$ can't have common variables.

The proporional hazard model is defined conditionaly on the same class $g$ as :

  $$\lambda(t)|_{c_{i}=g} = \lambda_{0g}(t)\exp(X_{4i}\psi+X_{5i}\eta_g)$$  

where $X_{4i}$ and $X_{5i}$ are vectors of covariates respectively associated with common effects aver classes $\psi$ and class-specific effects $\eta_g$.

# Data We use the paquid sample included in the package. Please refer to the introduction vignette for more details about these data. We consider here only the subjects at risk of dementia at the begining of the study : wzxhzdk:1 We also create some variables that will be used in the example : wzxhzdk:2 # First use of Jointlcmm function We model jointly the trajectory of normMMSE and time to dementia. As a JLCM is estimated for a fixed number of latent classes, we begin by specifying the model with 1 latent class. ## 1. Linear mixed model for normMMSE trajectory We begin by specifying the linear mixed model for normMMSE. We will consider the same specification as in the hlme vignette, that is a quadratic trajectory with age adjusted for CEP. wzxhzdk:3 ## 2. Survival model for dementia diagnosis Joinlcmm assumes a parametric baseline risk function. We thus need to determine the family of baseline risks. To do so, we will use the jointlcmm function in which the longitudinal part will be the same as in hlme. In the application, the risk of dementia is described according to age so we have a problem of delayed entry. The program handles it by specifying age_init in the Surv object. We try different families of baseline risks (Weibull, Splines, piecewise constant) and we systematically adjust on CEP and male. wzxhzdk:4 The Weibull model gives the best fit. wzxhzdk:5 It is possible to change the parameterization of the survival model (with log instead of +/-sqrt) wzxhzdk:6 ## 3. Estimation with different numbers of latent classes Once the specification of the model under G=1 is done (one class means independent models), we can estimate the model with more than one class. This is a cumbersome step in the analysis since the estimation has to be replicated for different numbers of latent classes AND various initial values to avoid the convergence toward a local maximum. ### 3.1 Model with two latent classes wzxhzdk:7 We could choose to consider proportional hazards given the fit obtained but we continue with class-specific baseline risk functions to allow for more flexibility. Note that the results in terms of classification seem to be very close. wzxhzdk:8 The estimates of the model are in the summary. wzxhzdk:9 From this first model, we can look at different output functions available in the package to evaluate the quality of fit of the model. wzxhzdk:10 wzxhzdk:11 wzxhzdk:12 And the predictions of the model: wzxhzdk:13 The model obtained with the first call of jointlcmm is not necessarily the maximum likelihood estimator for 2 classes. The model must be refitted with other initial values. There are different possibilities in the package: - random departure from the asymptotic distribution of the estimates under G=1 - initial values chosen by the user - a grid search with replicates R times the random departures with a maximu of M iterations of the algorithm each time. The program finishes the estimation with the departure which gave the best log likelihood after the M iterations. This is what is recommended with latent class models to ensure the convergence toward the global maximum. The grid search can take a lot of time as replicating R model estimation. I recommend to use 100 random departures and if possible between 30-50 iterations. Here, we will illustrate the procedure with less replicates and iterations to reduce the processing time. wzxhzdk:14 In this example, we always converge to the same maximum whatever the departure but the latent classes might be exchanged. wzxhzdk:15 ### 3.2 Model with more than 2 latent classes Estimation of the models with 3 and 4 classes from default values (based on G=1 model estimates). wzxhzdk:16 The program did not converge properly after 100 iterations (the maximum number of iterations specified by default). This is explained by the fact that the program converged to an non optimal point (local maximum). Indeed, when looking at the results (loglikelihood and estimates), the model mj3 converged toward the 1 class model and the model mj4 converged toward the 2 class model. wzxhzdk:17 Such a problem may happen and the user should be careful with this issue. We can use other initial values to search for the global maximum in another part of the parameter space: wzxhzdk:18 wzxhzdk:19 From the summary, we choose the 4 class model. wzxhzdk:20 The score test statistics rejects the conditional independence at the 5% level but the shape of the curve of ST according to the number of class (with the asymptot) tells us that we will probably never reach the significance level. ## 4. Analysis of the results with the 4 class model ## 4.1. Summary of the estimation The summary gives all the important information wzxhzdk:21 ## 4.2. Evaluation of the classification A major element of joint latent class models is the posterior classification, and the discriminatory power of this classification. wzxhzdk:22 ## 4.3 Graphs for the fit We can assess the fit of the model by comparing the predictions to the observations by time intervals. wzxhzdk:23 wzxhzdk:24 The residuals of the longitudinal part: wzxhzdk:25 ## 4.4. Graph of predicted trajectories according to a profile of covariates wzxhzdk:26 ## 4.5. Graph of predicted cumulative incidence according to a profile of covariates wzxhzdk:27 wzxhzdk:28 # C. To go further ... ## 1. Individual dynamic prediction The joint latent class model can be used to provide individual dynamic prediction of the event from the observed repeated measures of the marker. This is usually done for a "new" subject (not in the estimation data) but for the class, we focus on subject 72: wzxhzdk:29 When the objective is to provide dynamic predictions, the predictive power of the model should be specifically checked using appropriate techniques. The package has an internal measure, the EPOCE which quantifies the pronostic information. For measures such as the AUC or Brier Score, other packages can be used from the dynamic predictions obtained with Jointlcmm, for instance timeROC (Blanche et al., Biometrics 2015). Here is an example of the use of EPOCE: wzxhzdk:30 wzxhzdk:31 Here the model with 3 latent classes seems to have a better predictive power than the 4 class model. ## 2. Competing risks Jointlcmm function can account for competing risks with the same structure of call. The only difference is in the definition of the time to event which is the minimum time between all the causes of event and the censoring. The indicator of event also indicates 0 for censoring or k for cause k. Suppose we have further information in the paquidS sample, namely a variable Age_CR that includes the first event between dementia and death, and Indic_CR that indicates the cause of event, dementia or death before dementia. In the following we give some examples of joint latent class models with competing risks. wzxhzdk:32 In the regression, we can consider specific effects with cause() : wzxhzdk:33 Note that different baseline risk functions can be considered for the two events: wzxhzdk:34 We can now use the same technique of estimation with G>1 for this model with competing risks ...

CecileProust-Lima/lcmm documentation built on March 3, 2024, 5:23 p.m.