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