library(tidyr)
library(nlme)
library(multcomp)
#Read in data for lateral root counts. Ensure that it is saved as a comma-separated values file (csv).
dat<-read.table("example_data_for_script.csv",header=TRUE,sep=",")
#Convert plant ID# into a factor.
dat[,2]<-as.factor(dat[,2])
#Drop column for air gap counts.
dat<-dat[,-c(4)]
#Calculate lateral root densities by dividing count in each category by length.
dat1<-cbind(dat[,1:2],dat[,3]/dat[,5],dat[,4]/dat[,5])
colnames(dat1)[3:4]<-c("Block","Plate")
#Create two new columns, "Side" and "Density", based on block and plate density columns.
dat1<-gather(dat1,Side,Density,c(Block,Plate))
#Separate "Condition" column into separate columns for the block and plate sides.
dat1<-separate(dat1,col=Condition, into=c("Block_Condition","Plate_Condition"), sep="/")
dat1[,1]<-as.character(dat1[,1])
dat1[,2]<-as.character(dat1[,2])
#Create two new columns for the stimulus presented on the current side of the root (cis)
#and the stimulus on the other side (trans) for each row.
cis_Condition<-vector()
trans_Condition<-vector()
for(i in 1:nrow(dat1)){
ifelse(dat1[i,"Side"]=="Block",{
cis_Condition<-c(cis_Condition,dat1[i,"Block_Condition"])
trans_Condition<-c(trans_Condition,dat1[i,"Plate_Condition"])},{
cis_Condition<-c(cis_Condition,dat1[i,"Plate_Condition"])
trans_Condition<-c(trans_Condition,dat1[i,"Block_Condition"])})
}
dat2<-cbind(cis_Condition,trans_Condition,dat1[,3:5])
#Perform ANOVA on linear mixed-effects model. Side of the root, cis-condition, and
#trans-condition are all treated as fixed effects, and plant ID is treated as a random
#effect to account for repeated measures on each plant.
lme.out<-lme(Density~Side*cis_Condition*trans_Condition,random=~1|Plant,data=dat2)
anova(lme.out)
#A significant interaction between cis-condition and trans-condition suggests that the
#two sides of the root are behaving non-autonomously. If this is the case, perform post-hoc
#analysis using Tukey's HSD test to compare densities within each combination of cis-
#and trans-conditions.
#Create a column for the interaction term between the cis- and trans-conditions.
#Levels of the interaction term are denoted as "cis-condition.trans-condition".
dat2$Interaction<-interaction(dat2$cis_Condition,dat2$trans_Condition)
#Create a new linear mixed-effects model based on the interaction term.
lme.out.1<-lme(Density~Interaction, random=~1|Plant,data=dat2)
#Perform all pairwise comparisons between the different levels of the interaction term using
#Tukey's HSD test.
summary(glht(lme.out.1, linfct=mcp(Interaction = "Tukey")))
#Assuming "A" and "B" are two levels of conditions, a significant difference between "A.A"
#and "A.B" suggests that the response to condition A depends on the condition on the opposite
#side of the root. Note that significant differences between "A.X" and "B.X" are to be expected
#if cis-condition itself has a significant effect on lateral root density.
1