Nrthugu

Probably an easy one.

I have data with fixed and random effects I'd like to fit a mixed effects model to:

set.seed(1)



df <- data.frame(group = c(rep("A",40),rep("B",40)),

                 treatment = rep(c(rep("T",20),rep("CT",20)),2),

                 class = c(rep("AT1",10),rep("ACT1",10),rep("AT2",10),rep("ACT2",10),rep("BT1",10),rep("BCT1",10),rep("BT2",10),rep("BCT2",10)),

                 value = rnorm(80),

                 stringsAsFactors = F)



df$group <- factor(df$group, levels = c("A","B"))

df$treatment <- factor(df$treatment, levels = c("CT","T"))

The fixed effects are group and treatment and the random effect is class, which to my understanding is nested within the group and treatment combinations.

The model I want to fit is:

value ~ group*treatment

Where the effect of interest if the group:treatment interaction.
Of course I want to account for class as a random effect, but I can't seem to find what the syntax for that is. I tried:
(1|group*treatment/class) and (1|group:treatment/class) but both give an error.

Defining a group:treatment column in df:

df <- df %>% dplyr::mutate(group_treatment = paste0(group,"_",treatment))

And fitting:

fit <- lmer(value ~ group*treatment + (1|group_treatment/class), data = df)

Does seem to work, but I'm wondering if that's the only way or whether there's a more explicit syntax for such cases of random effect nesting.

Any idea?