Help with setting up a random effects model

I am slowly becoming familiar with mixed models. I am facing, however, some difficulties in defining a model for a new data set I am working on.

The dataset contains the reaction times RT for 50 subjects when they performed 3 tasks (the three tasks are A, B, and C). The task order was randomized between subjects. That is, subject 1 performed A, B, C; subject 2 performed B, A, C, and so forth. Each subject did the task only once, that is I have 3 RT measurements for each subject.

The data set could look like:

RT Subject Task Order
1.5 1 A I
1.1 1 B II
1.8 1 C III
1.0 2 A II
1.2 2 B I
1.5 2 C III

I would like to model the distribution of reaction times with a mixed model. I would like to know if there is any difference in reaction time for the task type, and the overall variability between subjects. One model could be:

RT ~ task + (task | subject)

I believe, however, that the order of the task was presented may have an effect on the reaction time. Therefore, by including the categorical variable order, I hope I would get some hints if there is a learning effect. The model would like:

RT ~ task * order + (task * order | subject)

Would this model be correct? Is this the way one would define a model to take into account the learning effect?

Is this question related to Stan or one of its interfaces?

It is a “general” question really. Once I understand how to set up a model, I would implement it in STAN or brms. I am sorry if this kind of question should not be ask in this forum.

If you are trying to get this model run in brms or rstanarm that it qualifies as a valid question I guess.

The model in general looks sensible, but you will likely not be able to fit it as you have just 3 RTs per person (if I understood you correctly), but 3 * 3 = 9 varying effects per person.

In fact, I tried to fit it in brms with this formula and with the lognormal family. The model was very slow to run and got many divergences and high r_hat values. What would be your suggestion? To remove some random effects (e.g., remove the interactions)?

You seem to have too few data to reliably fit more than 2 varying effects, right now you have 9.

Ok, that’s good information anyway, thank you.

You counted 3 \cdot 3=9 varying effects per person because both task and order have 3 levels? Thus, my design matrix X for the fixed effects would have 9 columns:

TaskA
TaskB
TaskC
OrderII
OrderIII
TaskB:OrderII
TaskC:OrderII
TaskB:OrderIII
TaskC:OrderIII

Then, the Z matrix for the random effects would have a 9 coefficients for each subject as well (i.e., 9 \cdot 50 = 450 total columns).

Given the measurements I have, you suggest that I can afford only 2 random effects per subject, which are too few even for the simplest model with only the task type as variable, because I have 3 coefficients to estimate.

TaskA
TaskB
TaskC

Do you have any suggestions on how to get around this issue?

Run more trials per person, trial and order. Otherwise, you have to use
RT ~ task * order + (1 | subject)

1 Like

Thank you, I really appreciate your insights.