CecileProust-Lima / lcmm

R package lcmm
https://CecileProust-Lima.github.io/lcmm/
55 stars 13 forks source link

Matching the mathematical model with the R code for lcmm #31

Closed MaryWaterloo closed 4 years ago

MaryWaterloo commented 5 years ago

Dear Dr. Proust-Lima,

I am running your lcmm package on a sample of 500 older adults with up to 12 longitudinal assessments that created a 4-level ordinal outcome. I am running the model without any covariates. Based on my understanding, I am using a version of formula number 8 in your 2017 article (Proust-Lima, C., Philipps, V., & Liquet, B. (2017). Estimation of extended mixed models using latent classes and latent processes: the R package lcmm. Journal of Statistical Software, 78(2), 1-56. doi:10.18637/jss.v078.i02) Ʌ_i (t)|(c_i=g)=X_L1i (t)^T β+X_L2i (t)^T α_g+Z_i (t_ij)^T u_ig+ω_i (t_ij)

My best model in R is a model with random intercept and quadratic for the effect of age (I think because my outcome has only 4 levels and the maximum follow-up period is 12 assessments, adding random slope did not improve my model. I mean he variablilty in slope is not enough that can be modeled) m <- lcmm(outcome4 ~ poly (age, degree = 2, raw = TRUE), random = ~ 1, mixture = ~ poly (age, degree = 2, raw = TRUE), subject = "ID", data = long, link = "thresholds",ng=4)

I have two questions that I greatly appreciate your advise on them:

Question 1: Did I introduce my model to R correctly? I want to fit an lcmm model with random intercept and quadratic for the effect of age.

Question 2: How does my lcmm model in R match with formula (8) in your article? which parameters did I estimate and which ones I did not estimate?

Kind regards, Maryam

VivianePhilipps commented 5 years ago

Hello Maryam,

you are right, the lcmm call you use will effectively estimate a model with quadratic trend with age and random intercept.

For each class, you have 3 fixed effects : intercept, degree 1 polynomial and degree 2 polynomial. So the X_L2i(t) vector contains these 3 values. The corresponding paramters are (α0_1, α1_1, α2_1) for class 1, (α0_2, α1_2, α2_2) for class 2,(α0_3, α1_3, α2_3) for class 3, (α0_3, α1_3, α2_3) for class 4. All parameters are estimated, excepted α0_1 which is fixed to 0. But you don't need to specify it, it is internally done in the function.

X_L1i(t) should contain all variables appearing in the 'fixed' arguments of the lcmm function (ie outcome4 ~ poly (age, degree = 2, raw = TRUE)) but not appearing in the 'mixture argument'. You don't have such covariate, so you don't have any X_l1i(t) vector.

Z_i(t) contains the random effects, and as you only have an intercept, Z_i(t) = 1. And you estimate the variance of u_i.

ω_i (t) is fir an additional correlation (with cor argument). You do not mentioned this option, so you do not have any ω_i (t).

In summary, you have 12 parameters to estimate for this part of the model, 11 fixed effects and 1 variance.

Hope this clarify the question,

Viviane