As I understand, you have a Between x Within x Within x Within design.

In attachment an R script for such a design. It includes a parametric ANOVA, which might be helpfull to select the models you want to test. It also includes the Shapiro test and figures of the interactions.

The ANOVA is from the 'ez' package. ges = generalized eta squared effect size

Here's an idea: With ARTool in R, if you treat the anova result as an object --- e.g. Result = anova(model) ---, you can extract the residuals. With the actual y values and the residuals, you can calculate various forms of "accuracy" for the model, like Ephron's pseudo r-square, MAE, MAPE, RMSE, and so on. You could compare these across models. It wouldn't be a test per se, but might be a reasonable way to compare models.

Below is an example of extracting Efron's pseudo R-squared from an ARTool object. I added a function to my rcompanion package to do so. You can also request other statistics like MAPE and NRMSE.

Because the ARTool model object has neither the actual y values nor the predicted values, you have to manually pass the actual y values to the function.

If I understand correctly, the random effects in the model call don't affect the residuals returned in the model object. That is, they affect the p values, but not the residuals. The upshot is that the Efron pseudo r-squared will not be affected by random terms in the model.

# # #

### Example from: https://rcompanion.org/handbook/F_16.html

Location = as.factor(c(rep("Olympia" , 6), rep("Ventura", 6),

rep("Northampton", 6), rep("Burlington", 6)))

Tribe = as.factor(c(rep(c("Jedi", "Sith"), 12)))

Midichlorians = c(10, 4, 12, 5, 15, 4, 15, 9, 15, 11, 18, 12,

8, 13, 8, 15, 10, 17, 22, 22, 20, 22, 20, 25)

Data = data.frame(Tribe, Location, Midichlorians)

library(ARTool)

model = art(Midichlorians ~ Tribe + Location + Tribe:Location, data = Data)

anova(model)

library(rcompanion)

efronRSquared(actual=Data$Midichlorians, residual=model$residuals)