The correct random slope model for nested data

Question

I'm trying to see how personalities of individuals change with time. The variables in my data are: 1. latency to emerge (response variable in continuous scale) measured for 204 individuals from 14 colonies. Each individual was measured for 6 times (6 trials). The 204 individuals came from 14 different colonies, so individuals (ID) are nested within colonies. Individuals have distinct IDs in the data set such that each ID can belong to only one of the 14 colonies (for example: IDs 1-20 in colony 1, 21-40 in colony 2 and so on). In this scenario, there is an excellent link here which shows that the model output will be the same for (1|colony)+(1|ID) and (1|colony/ID).

I would like to see how the latency to emerge of individuals (IDs) change with time/trials. Remember, individuals are not independent, but nested within colonies. Therefore, I made two models:

m1 <- lmer(Latency~Trial+(1|colony)+(Trial|ID),data=mydata) and

m2 <- lmer(Latency~Trial+(Trial|ID),data=mydata).

I know that model m2 assumes that individuals are independent of colonies (while in reality, individuals are non-independent of their respective colonies). Does model m1 account for this lack of independence? How is model 1 different from model 2 in terms of the random slopes plotted against each individual?

Answer 1

The model you want is probably

Latency~Trial+(Trial|colony)+(Trial|ID)

(or equivalently Latency ~ Trial + (Trial|colony/ID)): the corresponding statistical model for observation $k$ of individual $j$ in colony $i$ (at time $t$) is

$$ \begin{split} y_ijk(t) = & (\beta_0 + \epsilon_{c0,i} + \epsilon_{d0,j}) + \\ & (\beta_1 t + \epsilon_{c1,i} + \epsilon_{d1,j}) t + \epsilon_{r,k} \end{split} $$

where $\{\epsilon_{c0},\epsilon_{c1}\}$ (colony-level random effect on intercept and slope) and $\{\epsilon_{d0},\epsilon_{d1}\}$ (individual-level random effect on intercept and slope) are bivariate Normal deviates, and $\epsilon_r$ (the residual error) is univariate Normal.

Your first model doesn't allow for among-colony variation in slopes; as you state in your question, your second model allows for neither intercept nor slope variation among colonies.

Depending on what your data looks like and what's typical in your field, you might want to consider a log-Normal model (use log(Latency) as your response variable) or a Gamma GLM with a log link ...

Answer 2

Let analyze the random effects in the model.

From the description, the response variable from different colonies are independent. Let check the 2 individuals from the same colony and measurements at time 1 and 6. $Y_{11}, Y_{16}, Y_{21}, Y_{26}$.

The model m2 is: $$Y_{ij} = X\beta + \gamma_{1i} + j\gamma_{2i} + \epsilon_{ij} $$ where $Var(\gamma_{1i}) = \sigma_1^2$, $Var(\gamma_{2i}) = \sigma_2^2$, $Var(\epsilon_{ij}) = \sigma_0^2$, Then the variance-covariance matrix of four $Y$s is $$\Sigma= \left(\begin{matrix}\sigma_1^2+\sigma_2^2+\sigma_0^2&\sigma_1^2+6\sigma_2^2&0&0\\ \sigma_1^2+6\sigma_2^2&\sigma_1^2+36\sigma_2^2+\sigma_0^2&0&0 \\ 0& 0 & \sigma_1^2+\sigma_2^2+\sigma_0^2&\sigma_1^2+6\sigma_2^2\\ 0&0&\sigma_1^2+6\sigma_2^2&\sigma_1^2+36\sigma_2^2+\sigma_0^2\\ \end{matrix} \right)$$

In this model, there are two problems: (1) individual 1 and 2 come from the same colony but covariance is zero, (2) $Y$ at the time 1 and 6 have different variance.

The model m1 is:

$$ Y_{ij} = X\beta + \gamma_{1i} + j\gamma_{2i} + \gamma_3 +\epsilon_{ij} $$ Let $Var(\gamma_3) = \sigma_3^2$.

Then the ariance-covariance matrix of four $Y$s under m1 is adding $\sigma_3^2$ to all of 16 elements of $\Sigma$. So the first problem in m2 is resolved. But the second problem is still there. If you believe that variance of response variable increases along the time, m1 is model that you wanted. If you think it is unreasonable that variance of response variable increases along the time, then the model m3 maybe more reasonable.

$$ Y_{ij} = X\beta + \gamma_{1i} + \gamma_3 +\epsilon_{ij} $$ Under model m3, the ariance-covariance matrix of four $Y$s is: $$\Sigma_3= \left(\begin{matrix}\sigma_1^2+\sigma_3^2+\sigma_0^2&\sigma_1^2+\sigma_3^2&0&0\\ \sigma_1^2+\sigma_3^2&\sigma_1^2+\sigma_3^2+\sigma_0^2&0&0 \\ 0& 0 & \sigma_1^2+\sigma_3^2+\sigma_0^2&\sigma_1^2+\sigma_3^2\\ 0&0&\sigma_1^2+\sigma_3^2&\sigma_1^2+\sigma_3^2+\sigma_0^2\\ \end{matrix} \right)$$
So when you construct the random effect part in the model, need to check which one in more suitable for your situation.

In R, m3 canbe m3<-lmer(Latency~Trial+(1|colony)+(1|ID),data=mydata)

The discussion above is based on the assumption that value of the time (variable trial) are 1,2,3,4,5 and 6 and is continue variable in the model.