I'm trying to analyze a data set but I'm having trouble finding the right analysis. I'm looking at performance (Work) of a material.
Unfortunately my experimental design turned out to be unbalanced. I had 16 individuals set up in one of two treatments (A or B) (between treatments). After treatment we collected material and switch the individuals to the other treatment and collected material again. Individuals were switch about 6 times, and I was not able to collect material each time.
After collection I measured performance (Work) at four different conditions (a, b, c or d)(within treatment).
I want to see if Treatment had an effect on the performance of the material in any of the conditions. So, I ran a linear mixed model on the log transformed Work data using R. The original Work data did not met ANOVA assumptions. My response variable is Log of Work, my factors are between treatment and within treatment, and I used individual and order as my random effects.
I set zero contrast to factors, including random effects. Is that OK? Should I use a different analysis for my type of data?
Thank you for your help.
Here is the data file: https://docs.google.com/spreadsheets/d/1NnYWbvpjoCLAQJNaGnFDEf_DHdDPtRgjHSZdSe6ChVI/edit?usp=sharing
Here is my code:
library(lme4) # for glmer
library(multcomp) # for glht
library(emmeans) # for emm
# read in a data file
Mchnqs = read.csv("Example.csv")
colnames(Mchnqs)[1] <- gsub('^...','',colnames(Mchnqs)[1])
View(Mchnqs)
Mchnqs$Individual = factor(Mchnqs$Individual) # convert to nominal factor
Mchnqs$Order = factor(Mchnqs$Order) # convert to nominal factor
Mchnqs$Within.Treatment = factor(Mchnqs$Within.Treatment) # convert to nominal factor
Mchnqs$Betweenn.Treatment = factor(Mchnqs$Betweenn.Treatment) # convert to nominal factor
summary(Mchnqs)
with(Mchnqs, interaction.plot(Betweenn.Treatment, Within.Treatment, Log.Work)) # for convenience
boxplot(Log.Work ~ Betweenn.Treatment*Within.Treatment, data=Mchnqs, xlab="Between Treatment", ylab="Work (Log uJ)") # boxplots
# set sum-to-zero contrasts
contrasts(Mchnqs$Betweenn.Treatment) <- "contr.sum"
contrasts(Mchnqs$Within.Treatment) <- "contr.sum"
contrasts(Mchnqs$Order) <- "contr.sum"
contrasts(Mchnqs$Individual) <- "contr.sum"
# LMM test on Log.Work
# Betweenn.Treatment:Within.Treatment are
# fixed effects. Individual and Order are random effects.
m = lmer(Log.Work ~ (Betweenn.Treatment * Within.Treatment) + (1|Order) + (1|Individual), data=Mchnqs)
Anova(m, type=3)
qqnorm(residuals(m)); qqline(residuals(m)) # plot residuals
# perform post hoc pairwise comparisons with holm's bonferroni correction
summary(glht(m, emm(pairwise ~ Betweenn.Treatment * Within.Treatment)), test=adjusted(type="holm"))