WLC Data analysis - Respiration

RScript written by Erinne Stirling 2018

R version 3.5.1 (2018-07-02) -- "Feather Spray"
Copyright (C) 2018 The R Foundation for Statistical Computing
Platform: x86_64-apple-darwin15.6.0 (64-bit)

#Create the environment
library(readr)#for reading data
library(dplyr)#for manipulating data
library(ggplot2)#for plotting
library(agricolae)#for stastics

#Call the data
WLC<-read_csv("2018 STOTEN WLC Rscript respiration")%>%#Read the file
  transmute(TIME, REALTIME, CHAR, ID, RESP,#units in mg/gOC
            LITTER=ifelse(LITTER=="BLK", "Control",
                          ifelse(LITTER=="UL", "NFA", "FA")))

#Test for interactions
anova(lm(RESP~CHAR*LITTER, data=WLC))

WLC1<-WLC%>%
  filter(TIME==1)
anova(lm(RESP~CHAR*LITTER, data=WLC1))#***Litter
HSD.test(aov(RESP~LITTER, data=WLC1), c("LITTER"), console=TRUE, alpha=0.05)
#---------------
WLC2<-WLC%>%
  filter(TIME==2)
anova(lm(RESP~CHAR*LITTER, data=WLC2))#***Litter
HSD.test(aov(RESP~LITTER, data=WLC2), c("LITTER"), console=TRUE, alpha=0.05)

#---------------
WLC3<-WLC%>%
  filter(TIME==3)
anova(lm(RESP~CHAR*LITTER, data=WLC3))#***Litter
HSD.test(aov(RESP~LITTER, data=WLC3), c("LITTER"), console=TRUE, alpha=0.05)

#---------------
WLC4<-WLC%>%
  filter(TIME==4)
anova(lm(RESP~CHAR*LITTER, data=WLC4))#***Char:Litter
HSD.test(aov(RESP~CHAR*LITTER, data=WLC3), c("CHAR","LITTER"), console=TRUE, alpha=0.05)

#---------------
WLC5<-WLC%>%
  filter(TIME==5)
anova(lm(RESP~CHAR*LITTER, data=WLC5))#***Char:Litter
HSD.test(aov(RESP~CHAR*LITTER, data=WLC5), c("CHAR","LITTER"), console=TRUE, alpha=0.05)

#---------------
WLC6<-WLC%>%
  filter(TIME==6)
anova(lm(RESP~CHAR*LITTER, data=WLC6))#**Char:Litter
HSD.test(aov(RESP~CHAR*LITTER, data=WLC6), c("CHAR","LITTER"), console=TRUE, alpha=0.05)

#Adjust data for plotting
RESP_sd<-aggregate(WLC$RESP,
                   list(Litter = WLC$LITTER,
                        Char = WLC$CHAR, Time=WLC$TIME), FUN=sd)
RESP_mean<-aggregate(WLC$RESP,
                     list(Litter = WLC$LITTER,
                          Char = WLC$CHAR, Time=WLC$TIME), FUN=mean)
RESP_mean<-RESP_mean%>%
  mutate(sd=RESP_sd$x, Time=c("00","00","00","00","00","00",
                              "02","02","02","02","02","02",
                              "07","07","07","07","07","07",
                              "13","13","13","13","13","13",
                              "65","65","65","65","65","65",
                              "80","80","80","80","80","80",
                              "94","94","94","94","94","94"))%>%
  mutate(X=x, ymax=X+sd,ymin=X-sd,
         Litter=as.factor(Litter))%>%
  mutate(ymax=ifelse(x<sd,ifelse(X<0,ymax==0,ymax),ymax),
         ymin=ifelse(x<sd,ifelse(X>0,ymin==0,ymin),ymin))

levels(RESP_mean$Litter)[c(1,3,2)]
RESP_mean<-RESP_mean[order(RESP_mean$Char),]

#Plotting
ggplot(RESP_mean, aes(RESP_mean$Time, RESP_mean$X)) +
  facet_grid(Litter~Char) +
  geom_bar(stat="identity") + 
  geom_errorbar(aes(ymax=RESP_mean$ymax, ymin=RESP_mean$ymin), width=0.2) +
  xlab("Time (days)") +
  ylab(expression(paste("RESP (mgC/gOC)"))) +  scale_y_continuous(limits=c(0,320),breaks=c(0,50,100,150,200,250,300),labels=c(0,50,100,150,200,250,300)) +
  scale_x_discrete(labels=c("0","2","7","13","65","80","94")) +
  theme(panel.background=element_rect(fill="white"), panel.border = element_rect(linetype = "solid", fill = NA),
        strip.text.y=element_text(angle=0, size=12))