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 `task`

s (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.