# Convergence of Latent Class Linear Mixed Model with LCMM in R

I am seeking help getting my latent class linear mixed models to converge using the LCMM package in R.

I am studying common patterns of heroin and cocaine use over the life course of adults who have injected drugs. My data come from an active cohort study which has collected nearly 30 years of data on the drug use behaviors of adults in Baltimore, Maryland. The goal is to identify common trajectories of use over the life course.

To do this, I want to use latent class linear mixed models that model the conditional probability of past six-month use (a binary outcome) as a function of participant age. The outcome is modeled with a threshold (probit) link function. I included linear, quadratic, and cubic age terms that can vary by latent class, as well as a few other covariates (e.g., gender, race, cohort recruitment wave, decade) that are held constant across cohorts. I included a random intercept to account for clustering of observations within participants. I fit models with as few as g=1 to as many as g=6 latent classes.

Here is the code I used to estimate a 1-class model (i.e. a simple multilevlel model):

fit.1class <- lcmm(fixed = heroin ~ age.std + age2.std + age3.std + female + black + cohort + visit.era,
random = ~1,
subject = "newid",
data = alive.complete,
verbose = TRUE)


Here is the code I used to estimate an 3-class model, using the 1-class model for starting values (note that, in the code below, i = 3).

fit <- lcmm(fixed = heroin ~ age.std + age2.std + age3.std + female + black + cohort + visit.era,
random = ~1,
mixture = ~ age.std + age2.std + age3.std,
ng = i,
subject = "newid",
data = alive.complete,
maxiter = 500,
B = fit.1class,
verbose = TRUE)


This 3-class model did not converge in 500 iterations. Further, my simple understanding is that the convergence algorithm is "getting stuck" somehow, since the second derivative convergence criteria = 1, and this did not change after running for 500 more iterations (not shown).

One possible problem is that one of the latent classes that is being estimated is essentially 100% probability of use at every age. This is actually a reasonable class -- there are some participants who never stop using -- but it seems like this could cause problems for parameter estimation since the predicted model values are at the boundary.

Any suggestions for how to get these models with g>1 classes to converge? Here are some ideas I am considering, in no particular order:

• Fit a g>1-class model without random effects, and then use that model to generate starting values (instead of a one-class model without random effects to generate starting values).
• Allow the random intercept variance to vary across classes.
• Use natural smoothing splines (as generated by ns() in the "splines" package) instead of simple linear, quadratic, and cubic terms, to force the function to be linear outside the range of the data.
• Drop the quadratic and cubic terms entirely.
• Impose a reasonable constraint on some of the model parameters. For example, I could fix the conditional probability of use at age 20 to be 95% in all classes, a reasonable constraint since participants must have injected drugs to enroll in the study.
• Reconsider whether this modeling approach is really appropriate for my data. Maybe latent classes are just not a good way to model drug use over the lifecourse?

Any suggestions you have would be greatly appreciated.

The LCMM package in R has a huge bug which cause the model not to converge. If you use large numbers (e.g., years) as time points, the model do not converge. I think changing the time points to small numbers will fix the problem, without needing to change the convergence parameters.

You have highlighted one important issue: estimates on the boundaries of the parameter space can cause convergence issues. From the documentation:

the program may not converge and reach the maximum number of iterations fixed at 100. In this case, the user should check that parameter estimates at the last iteration are not on the boundaries of the parameter space. If the parameters are on the boundaries of the parameter space, the identifiability of the model is critical. This may happen especially with splines parameters that may be too close to 0 (lower boundary) or classmb parameters that are too high or low (perfect classification). When identifiability of some parameters is suspected, the program can be run again from the former estimates by fixing the suspected parameters to their value with option posfix. This usually solves the problem

If this does not help, the next thing I would do is to simplify the model by dropping the non-linear terms.

You can use the lcmm::gridsearch function to use random starting values. This should help, if the optimiser is getting stuck at a local minimum.

More generally, you can control the optimiser by specifying convl, the threshold for the convergence criterion based on the log-likelihood stability, convG, threshold for the convergence criterion based on the derivatives, maxiter, the maximum number of optimiser iterations

To expand on Robert's answer, use the grid search option to avoid local maxima (note that it's maxima, not minima because it's a maximum likelihood estimator), such as:

fit <- gridsearch(lcmm(fixed = heroin ~ age.std + age2.std + age3.std + female +
black + cohort + visit.era,
random = ~1,
mixture = ~ age.std + age2.std + age3.std,
ng = i,
subject = "newid",
data = alive.complete,
maxiter = 500,
B = fit.1class,
verbose = TRUE),
rep=20, maxiter=20)


Also, do you have the time dimension in the equation? For example, if you have measurements over 30 years, you would want to specify time (assuming your data is in long format), for example:

lcmm(fixed = heroin ~ time + covariates...)


And to allow for quadratic or cubic shapes in your time dimension, use:

lcmm(fixed = heroin ~ poly(time, degree = 3) + covariates...)