I have doubts about how to specify the random effects structure in my mixed models.
My experimental design consists of 12 independent enclosures (Encl) with populations of lizards (Subject*_ID).* We applied 2 crossed treatments (Lag: 3 levels and Var: 2 levels). And we repeated the experiment two years (Year), so individuals that survived the first year, were measured again the next year. We analyse the snout vent length (SVL) in mm. Sex (males and females). Individuals were redistributed to different enclosures and treatments the second year, so I include the interaction of Encl:Year in a new column (Encl_Year).
This was my model:
ctrl <- lmeControl(maxIter=200, msMaxIter=200, msVervose=TRUE, opt="optim")
options (contrast=c(factor="contr.sum", ordered="contr.poly"))
model.SVL <- lme(SVL~Lag*Var*Sex*Year, random=list(~1|Subject_ID, ~1|Encl_Year), weight=varIdent(form=~1|Lag*Var*Sex), control=ctrl, data=data)
But I don't know how it would be correct to define random effects. Since it would not be a cross random effect model, because not all levels (Subject_ID) are replicated in the enclosures (Encl:Year), but it is not nested either, because there are repeated individuals in different enclosures. What would be the most correct way to write the model?
Depending on the order:
random=list(~1|Subject_ID, ~1|Encl_Year)
or
random=list(~1|Encl_Year, ~1|Subject_ID)
, the results change quite a lot. I also tried a cross random effect model:
data$Dummy <- factor(1)
data <- groupedData(SVL ~ 1 | Dummy, data)
model.SVL <- lme(SVL~Lag*Var*Sex*Year, random=pdBlocked(list(pdIdent(~ 0 + Subject_ID), pdIdent(~ 0 + Encl_Year))) control=ctrl, weight=varIdent(form=~1|Lag*Var*Sex), data=data)
I should add, that I use the lme function, because there is heteroscedasticity in the residuals that I have corrected with the varident function.
