
##Test2: Root diameter/ root systems of known dimensions in laboratory_Figure4

library (reshape)

addr_input = "D://"
addr_output="D://"
##

morpDF = read.table(file = paste (addr_input, "data.Figure 4.CSV", sep = ""), header = T,  sep = ";", dec=".")
morpDF$MethodName = factor (ifelse (morpDF$Method == "C", "Camera",
                                  ifelse (morpDF$Method == "F", "Flatbed",
                                          ifelse (morpDF$Method == "H", "Handheld",
                                                  ifelse (morpDF$Method == "S", "Smartphone",  "Error")))))


# --------------------------------
# COMPUTE THE RELATIVE DIAMETER

morpDF <- morpDF[!(morpDF$MeanDiameter_mm == 0 & morpDF$Method == "F"),] # Remove the data where method = 1 and RER = 0
morpDF$RootName <- factor(paste(morpDF$Rhizotron, morpDF$RootName, sep="-"))
totalLevel = levels(factor(morpDF$RootName))
for (level in totalLevel) {
  if(length(morpDF$MeanDiameter_mm[morpDF$RootName == level & morpDF$Method == "F"]) > 0){
    morpDF$diam_rel[morpDF$RootName == level] <- morpDF$MeanDiameter_mm[morpDF$RootName == level] / 
                                                    morpDF$MeanDiameter_mm[morpDF$RootName == level & morpDF$Method == "F"]
  }
}
morpDF <- morpDF[!is.na(morpDF$diam_rel),]

# --------------------------------

nameFig = "boxplot_mean_diam_amap_GL.fig4.pdf"

pdf (file = paste(addr_output, nameFig, sep =""), width = 7, height = 6 )

  x0 = morpDF$MethodName
  y0 = morpDF$diam_rel
  morpDF$id <- as.numeric(morpDF$MethodName)
  
  boxplot ( y0 ~ x0, type ="n", xlab = NA, ylab = NA,
            outline=T, boxwex=0.5)
  
  totalLevel = levels(factor(morpDF$RootName))
  for (level in totalLevel) {
    # level = "R0_S2"
    dataI = morpDF[morpDF$RootName == level,]
    points (x = dataI$id,
            y = dataI$diam_rel,
            col="#91765C60", pch = 21, bg = "#FF0D1760", cex = 1)
  }
  
  boxplot ( y0 ~ x0, type ="n", xlab = NA, ylab = NA, col="#FFFFFF50",
            outline=T, boxwex=0.5, add=T)
  
  abline(h=c(1.1, 0.9), lty=2)
  
  mtext (side = 1, "Methods", cex = 1.2, line = 3)
  mtext (side = 2, "Relative root mean diameter (mm)", cex = 1.2, line = 2.5)
  
dev.off()


## ANOVA AND NORMALITY TEST

lm1 = lm (formula = diam_rel ~ MethodName, data = morpDF)
summary(lm1)
fit <- (aov (lm1))
summary(fit)
TukeyHSD(fit)

shapiro.test(morpDF$diam_rel)
shapiro.test(morpDF$diam_rel ^ 2)
shapiro.test(morpDF$diam_rel ^ 0.5)
plot (density(morpDF$diam_rel), col = 3 )
hist (density(morpDF$diam_rel), col = 3 )
hist (morpDF$diam_rel )
par (mfrow = c(2,2) )
plot (density(morpDF$diam_rel), col = 3 )
plot (density(morpDF$diam_rel ^ 0.5), col = 3 )
plot (density(log(morpDF$diam_rel + 1) ), col = 3 )


lm4 = lm (formula = diam_rel ~ MethodName + as.factor(RootName), data = morpDF)
summary (lm4)
lm4 = lm (formula = diam_rel ~ MethodName * as.factor(RootName), data = morpDF)
summary (lm4)


