# Linear mixed effects model on physiological parameter as function of time

I have 1000 individuals with measurements of plasmatic HDL levels over the last 15 years. The number of measurements per individual ranges from 1 to 15, and were not taken in the same dates. I would like to compare the HDL levels (quantitative trait) over time in the individuals and make a ranking of those individuals whose HDL levels are consistently the highest over time. I want to adjust for sex and Tanner (categorical variable related to puberty development, it can be 1,2,3,4 or 5).

I have a data.frame (data) with 5 columns: ID (individual ID), HDL, Age (at the time of measurement), Sex and Tanner. This is a data frame with a toy example:

ID  HDL   Age   Sex  Tanner
1   50    12.3  M    2
1   52.1  15.4  M    4
1   55.3  17.1  M    5
2   45    12.1  M    1
2   46.3  13.1  M    1
3   60    14.3  F    3
3   55    16.2  F    5
... ...   ...   ...  ...


This is the model I figured out:

library(lme4)
model = lmer(HDL ~ Age + Sex + (1|Tanner) + (1|ID), data)


Age and Sex are my fixed variables and I expect to adjust for variation between Tanner. Individuals of the same Tanner should have more similar HDL values for reasons other than Age and Sex.

After running the regression, I extract the intercept and slope from each individual using coef(data)\$ID. Then, I multiply each intercept by the corresponding slope and get the final value I will use for the ranking.

Is my reasoning correct?

• You could start by searching for "r mixed models repeated measures", I'm sure you'd get a few ideas... – Dominic Comtois Aug 1 '17 at 2:33
• "Is my reasoning correct" strikes me as a software-neutral & on topic question. I'm voting to reopen. – gung - Reinstate Monica Aug 3 '17 at 16:15

## 1 Answer

You are on the right track, but you have to consider carefully the hidden assumptions in your model and whether you need to use a more complicated model to achieve your goals.

Treating ID as a random effect in this study with 1000 individuals makes a lot of sense. You have modeled this as (1|ID), however, so you are only allowing for intercepts to differ among individuals. That would seem to be consistent with your interest in "ranking ... individuals whose HDL levels are consistently the highest over time."

You say, however, that you wish to "extract the intercept and slope from each individual." (Emphasis added.) If that is the case, then you need to specify the slope (presumably for Age here) that you wish to allow to differ among individuals. Analyzing slopes might pose difficulties, however, with only 1 observation for some individuals.

Treating Tanner as a random effect, with only 5 levels, is not so obviously correct. If this is an ordered categorical variable for which you expect an ordered relation to HDL, you probably should model it as an ordinal predictor variable instead.

Specifying Age as a linear predictor makes sense if you think that HDL changes strictly linearly with age, other things being equal. I suppose that is possible, but it seems a bit unlikely. Often it makes more sense to model age in a more flexible way, for example with splines.

Finally, although your use of coef() to extract random effects from the model might work OK in this case, interpreting its output in more complicated models with nesting can be troublesome, as for example in this question. With the R lme4 package, the ranef() function might provide more intuitive results. It might make sense to look at both functions and make sure that you understand the types of results that they return.

• Thanks so much for your reply, @EdM. Most of the slopes for individual subjects are positive regarding HDL against Age, but a substantial number is negative, and the slopes vary a lot. I will then specify the slope as (Age | ID). I have removed subjects with <3 observations to avoid the problem of number of measurements = 1. Assuming HDL increases linearly with Age, would you say the following model is now correct?. model2 = lmer(HDL ~ Age + Sex + Tanner + (Age|ID), data). Thanks ! – Lucas Aug 4 '17 at 18:30
• I cannot give you the green check since I have less than 15 points (I registered here a couple of days ago) – Lucas Aug 4 '17 at 18:32
• @Lucas the specifications of mixed models like yours may include implicit assumptions about correlations between intercepts and slopes. See here and here for examples. As for your model being "correct," that depends on what you are trying to accomplish and how well your model represents the data. To me, intercepts alone would seem most related to being "consistently the highest over time." You will have to evaluate carefully how well your model fits your data, as for any regression. – EdM Aug 4 '17 at 19:00