SimpleToast 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)
library(dplyr)
library(ggplot2)
library(RColorBrewer)
library(agricolae)
library(reshape2)
library(nls2)

palette<-brewer.pal(n=5, name="RdYlBu")

###Prepare the Data
TOAST<- read.csv("CSV/CSV toasted resp_2.csv")%>%#'SimpleToast litter respiration.csv'
  mutate(TYPE=ifelse(SOIL=="NUB", "Native",
                     ifelse(SOIL=="NHB", "Native", "Pine")),
         STATUS=ifelse(SOIL=="NUB", "Unburnt", 
                       ifelse(SOIL=="PUB", "Unburnt", "Burnt")), 
         LITTER=ifelse(LITTER=="BLK", "Control",
                       ifelse(LITTER=="LNU", "N UnSc",
                              ifelse(LITTER=="LNS", "N Sc", 
                                     ifelse(LITTER=="LPU", "P UnSc", "P Sc")))))%>%
    select(TYPE, STATUS, LITTER, TIME, CUM_RESP)#Load the data

#Decay Curves
TOAST_blank<-TOAST%>%
  filter(LITTER=="Control")%>%
  group_by(STATUS,TYPE, TIME)

plot(TOAST_blank$TIME,TOAST_blank$CUM_RESP)

##----
fit_1<-TOAST_blank%>%
  filter(STATUS=="Burnt", TYPE=="Native")
fit_a<-nls(CUM_RESP~C*(1-exp(k*TIME)),data=fit_1,
                     algorithm="port",
                     start=c(C=10,k=-1),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
a<-summary(fit_a)
fit_1$pred <- predict(fit_a)
plot(fit_1$TIME,fit_1$CUM_RESP)
lines(fit_1$TIME,fit_1$pred, col="blue",lty=2)

##----
fit_2<-TOAST_blank%>%
  filter(STATUS=="Burnt", TYPE=="Pine")
fit_b<-nls(CUM_RESP~C*(1-exp(k*TIME)),data=fit_2,
         algorithm="port",
         start=c(C=10,k=-1),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
b<-summary(fit_b)
fit_2$pred <- predict(fit_b)
plot(fit_2$TIME,fit_2$CUM_RESP)
lines(fit_2$TIME,fit_2$pred, col="blue",lty=2)

##----
fit_3<-TOAST_blank%>%
  filter(STATUS=="Unburnt", TYPE=="Native")
fit_c<-nls(CUM_RESP~C*(1-exp(k*TIME)),data=fit_3,
           algorithm="port",
           start=c(C=10,k=-1),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
c<-summary(fit_c)
fit_3$pred <- predict(fit_c)
plot(fit_3$TIME,fit_3$CUM_RESP)
lines(fit_3$TIME,fit_3$pred, col="blue",lty=2)

##----
fit_4<-TOAST_blank%>%
  filter(STATUS=="Unburnt", TYPE=="Pine")
fit_d<-nls2(CUM_RESP~C*(1-exp(k*TIME)),data=fit_4,
           algorithm="brute-force",
           start=c(C=4,k=-0.3),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
d<-summary(fit_d)
fit_4$pred <- predict(fit_d)
plot(fit_4$TIME,fit_4$CUM_RESP)
lines(fit_4$TIME,fit_4$pred, col="blue",lty=2)


TOAST_blank_ave<-data.frame(id=c(rep("Blank",4)),
  C=c(a$parameters[1,1],b$parameters[1,1],c$parameters[1,1],d$parameters[1,1]),
  k=c(a$parameters[2,1],b$parameters[2,1],c$parameters[2,1],d$parameters[2,1]))

##--------------------Pine Scorched litter
TOAST_pSc<-TOAST%>%
  filter(LITTER=="P Sc")%>%
  group_by(STATUS,TYPE, TIME)%>%
  mutate(CUM_RESP=1.3*CUM_RESP)

##----
fit_1<-TOAST_pSc%>%
  filter(STATUS=="Burnt", TYPE=="Native")
fit_a<-nls(CUM_RESP~C*(1-exp(k*TIME)),data=fit_1,
           algorithm="port",
           start=c(C=10,k=-1),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
a<-summary(fit_a)
fit_1$pred <- predict(fit_a)
plot(fit_1$TIME,fit_1$CUM_RESP)
lines(fit_1$TIME,fit_1$pred, col="blue",lty=2)

##----
fit_2<-TOAST_pSc%>%
  filter(STATUS=="Burnt", TYPE=="Pine")
fit_b<-nls(CUM_RESP~C*(1-exp(k*TIME)),data=fit_2,
           algorithm="port",
           start=c(C=10,k=-1),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
b<-summary(fit_b)
fit_2$pred <- predict(fit_b)
plot(fit_2$TIME,fit_2$CUM_RESP)
lines(fit_2$TIME,fit_2$pred, col="blue",lty=2)

##----
fit_3<-TOAST_pSc%>%
  filter(STATUS=="Unburnt", TYPE=="Native")
fit_c<-nls(CUM_RESP~C*(1-exp(k*TIME)),data=fit_3,
           algorithm="port",
           start=c(C=10,k=-1),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
c<-summary(fit_c)
fit_3$pred <- predict(fit_c)
plot(fit_3$TIME,fit_3$CUM_RESP)
lines(fit_3$TIME,fit_3$pred, col="blue",lty=2)

##----
fit_4<-TOAST_pSc%>%
  filter(STATUS=="Unburnt", TYPE=="Pine")
fit_d<-nls(CUM_RESP~C*(1-exp(k*TIME)),data=fit_4,
           algorithm="port",
           start=c(C=10,k=-1),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
d<-summary(fit_d)
fit_4$pred <- predict(fit_d)
plot(fit_4$TIME,fit_4$CUM_RESP)
lines(fit_4$TIME,fit_4$pred, col="blue",lty=2)


TOAST_pSc_ave<-data.frame(id=c(rep("P Sc",4)),
  C=c(a$parameters[1,1],b$parameters[1,1],c$parameters[1,1],d$parameters[1,1]),
  k=c(a$parameters[2,1],b$parameters[2,1],c$parameters[2,1],d$parameters[2,1]))

##--------------------Pine UnScorched litter
TOAST_pUnSc<-TOAST%>%
  filter(LITTER=="P UnSc")%>%
  group_by(STATUS,TYPE, TIME)%>%
  mutate(CUM_RESP=1.3*CUM_RESP)

##----
fit_1<-TOAST_pUnSc%>%
  filter(STATUS=="Burnt", TYPE=="Native")
fit_a<-nls2(CUM_RESP~C*(1-exp(k*TIME)),data=fit_1,
           algorithm="brute-force",
           start=c(C=12,k=-0.2),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
a<-summary(fit_a)
fit_1$pred <- predict(fit_a)
plot(fit_1$TIME,fit_1$CUM_RESP)
lines(fit_1$TIME,fit_1$pred, col="blue",lty=2)

##----
fit_2<-TOAST_pUnSc%>%
  filter(STATUS=="Burnt", TYPE=="Pine")
fit_b<-nls(CUM_RESP~C*(1-exp(k*TIME)),data=fit_2,
           algorithm="port",
           start=c(C=10,k=-1),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
b<-summary(fit_b)
fit_2$pred <- predict(fit_b)
plot(fit_2$TIME,fit_2$CUM_RESP)
lines(fit_2$TIME,fit_2$pred, col="blue",lty=2)

##----
fit_3<-TOAST_pUnSc%>%
  filter(STATUS=="Unburnt", TYPE=="Native")
fit_c<-nls(CUM_RESP~C*(1-exp(k*TIME)),data=fit_3,
           algorithm="port",
           start=c(C=10,k=-1),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
c<-summary(fit_c)
fit_3$pred <- predict(fit_c)
plot(fit_3$TIME,fit_3$CUM_RESP)
lines(fit_3$TIME,fit_3$pred, col="blue",lty=2)

##----
fit_4<-TOAST_pUnSc%>%
  filter(STATUS=="Unburnt", TYPE=="Pine")
fit_d<-nls(CUM_RESP~C*(1-exp(k*TIME)),data=fit_4,
           algorithm="port",
           start=c(C=10,k=-1),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
d<-summary(fit_d)
fit_4$pred <- predict(fit_d)
plot(fit_4$TIME,fit_4$CUM_RESP)
lines(fit_4$TIME,fit_4$pred, col="blue",lty=2)


TOAST_pUnSc_ave<-data.frame(id=c(rep("P UnSc",4)),
                          C=c(a$parameters[1,1],b$parameters[1,1],c$parameters[1,1],d$parameters[1,1]),
                          k=c(a$parameters[2,1],b$parameters[2,1],c$parameters[2,1],d$parameters[2,1]))

##--------------------Native Scorched litter
TOAST_nSc<-TOAST%>%
  filter(LITTER=="N Sc")%>%
  group_by(STATUS,TYPE, TIME)%>%
  mutate(CUM_RESP=1.3*CUM_RESP)

##----
fit_1<-TOAST_nSc%>%
  filter(STATUS=="Burnt", TYPE=="Native")
fit_a<-nls(CUM_RESP~C*(1-exp(k*TIME)),data=fit_1,
           algorithm="port",
           start=c(C=10,k=-1),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
a<-summary(fit_a)
fit_1$pred <- predict(fit_a)
plot(fit_1$TIME,fit_1$CUM_RESP)
lines(fit_1$TIME,fit_1$pred, col="blue",lty=2)

##----
fit_2<-TOAST_nSc%>%
  filter(STATUS=="Burnt", TYPE=="Pine")
fit_b<-nls2(CUM_RESP~C*(1-exp(k*TIME)),data=fit_2,
           algorithm="brute-force",
           start=c(C=12,k=-0.35),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
b<-summary(fit_b)
fit_2$pred <- predict(fit_b)
plot(fit_2$TIME,fit_2$CUM_RESP)
lines(fit_2$TIME,fit_2$pred, col="blue",lty=2)

##----
fit_3<-TOAST_nSc%>%
  filter(STATUS=="Unburnt", TYPE=="Native")
fit_c<-nls(CUM_RESP~C*(1-exp(k*TIME)),data=fit_3,
           algorithm="port",
           start=c(C=10,k=-1),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
c<-summary(fit_c)
fit_3$pred <- predict(fit_c)
plot(fit_3$TIME,fit_3$CUM_RESP)
lines(fit_3$TIME,fit_3$pred, col="blue",lty=2)

##----
fit_4<-TOAST_nSc%>%
  filter(STATUS=="Unburnt", TYPE=="Pine")
fit_d<-nls(CUM_RESP~C*(1-exp(k*TIME)),data=fit_4,
           algorithm="port",
           start=c(C=10,k=-1),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
d<-summary(fit_d)
fit_4$pred <- predict(fit_d)
plot(fit_4$TIME,fit_4$CUM_RESP)
lines(fit_4$TIME,fit_4$pred, col="blue",lty=2)


TOAST_nSc_ave<-data.frame(id=c(rep("N Sc",4)),
                          C=c(a$parameters[1,1],b$parameters[1,1],c$parameters[1,1],d$parameters[1,1]),
                          k=c(a$parameters[2,1],b$parameters[2,1],c$parameters[2,1],d$parameters[2,1]))

##--------------------Native Uncorched litter
TOAST_nUnSc<-TOAST%>%
  filter(LITTER=="N UnSc")%>%
  group_by(STATUS,TYPE, TIME)%>%
  mutate(CUM_RESP=1.3*CUM_RESP)

##----
fit_1<-TOAST_nUnSc%>%
  filter(STATUS=="Burnt", TYPE=="Native")
fit_a<-nls(CUM_RESP~C*(1-exp(k*TIME)),data=fit_1,
           algorithm="port",
           start=c(C=10,k=-1),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
a<-summary(fit_a)
fit_1$pred <- predict(fit_a)
plot(fit_1$TIME,fit_1$CUM_RESP)
lines(fit_1$TIME,fit_1$pred, col="blue",lty=2)

##----
fit_2<-TOAST_nUnSc%>%
  filter(STATUS=="Burnt", TYPE=="Pine")
fit_b<-nls(CUM_RESP~C*(1-exp(k*TIME)),data=fit_2,
           algorithm="port",
           start=c(C=10,k=-1),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
b<-summary(fit_b)
fit_2$pred <- predict(fit_b)
plot(fit_2$TIME,fit_2$CUM_RESP)
lines(fit_2$TIME,fit_2$pred, col="blue",lty=2)

##----
fit_3<-TOAST_nUnSc%>%
  filter(STATUS=="Unburnt", TYPE=="Native")
fit_c<-nls(CUM_RESP~C*(1-exp(k*TIME)),data=fit_3,
           algorithm="port",
           start=c(C=10,k=-1),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
c<-summary(fit_c)
fit_3$pred <- predict(fit_c)
plot(fit_3$TIME,fit_3$CUM_RESP)
lines(fit_3$TIME,fit_3$pred, col="blue",lty=2)

##----
fit_4<-TOAST_nUnSc%>%
  filter(STATUS=="Unburnt", TYPE=="Pine")
fit_d<-nls(CUM_RESP~C*(1-exp(k*TIME)),data=fit_4,
           algorithm="port",
           start=c(C=10,k=-1),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
d<-summary(fit_d)
fit_4$pred <- predict(fit_d)
plot(fit_4$TIME,fit_4$CUM_RESP)
lines(fit_4$TIME,fit_4$pred, col="blue",lty=2)

TOAST_nUnSc_ave<-data.frame(id=c(rep("N UnSc",4)),
                            C=c(a$parameters[1,1],b$parameters[1,1],c$parameters[1,1],d$parameters[1,1]),
                            k=c(a$parameters[2,1],b$parameters[2,1],c$parameters[2,1],d$parameters[2,1]))

Curves<-rbind(TOAST_blank_ave,TOAST_pSc_ave,TOAST_pUnSc_ave,TOAST_nSc_ave,TOAST_nUnSc_ave)

##-----Stats

stats<-Curves%>%
  group_by(id)

anova(aov(C~id, data=stats)) #**
HSD.test(aov(C~id, data=stats), c("id"),
         console=TRUE,alpha=0.05)

aggregate(C~id, data=stats, FUN="sd")

anova(aov(k~id, data=stats)) #ns
HSD.test(aov(k~id, data=stats), c("id"),
         console=TRUE,alpha=0.05)

aggregate(k~id, data=stats, FUN="sd")

####------------Plot
fun_control <- function(x) {
  3.673312*(1-exp(-0.2615971*x))
}

fun_NUnSc <- function(x) {
  15.870470*(1-exp(-0.2869384*x))
}

fun_NSc <- function(x) {
  22.838834*(1-exp(-0.2462890*x))
}

fun_PUnSc <- function(x) {
  16.454169*(1-exp(-0.2105133*x))
}

fun_PSc <- function(x) {
  35.715913*(1-exp(-0.2884479*x))
}

TOAST_lines<-ggplot(data.frame(x=c(0,8)),aes(x=x))+
      stat_function(fun=fun_control, colour="grey")+
      stat_function(fun=fun_NUnSc, colour="lightblue")+
      stat_function(fun=fun_NSc, colour="blue")+
      stat_function(fun=fun_PUnSc, colour="orange")+
      stat_function(fun=fun_PSc, colour="red")+
  scale_y_continuous(limits=c(0,40))+
  ggsave("TOAST_a5.svg", width=200, height=100, units="mm")

TOAST%>%
  ggplot(aes(x=TIME,y=CUM_RESP))+
  geom_point(aes(colour=LITTER))+
  scale_colour_manual(values=c("grey", "blue", "lightblue", "red", "orange"))+
  xlab("Time (days)") +
  ylab(expression(paste("Cumulative Respiration (",mu,"gCO2-C/gOC)"))) +
  scale_y_continuous(limits=c(0,40))+
  ggsave("TOAST_b5.svg", width=200, height=100, units="mm")

###Prepare the Data: Extracts
TOAST<- read.csv("CSV/CSV toasted resp ext_3.csv")%>%#'SimpleToast extract respiration.csv'
  mutate(TYPE=ifelse(SOIL=="NUB", "Native",
                     ifelse(SOIL=="NHB", "Native", "Pine")),
         STATUS=ifelse(SOIL=="NUB", "Unburnt", 
                       ifelse(SOIL=="PUB", "Unburnt", "Burnt")), 
         LITTER=ifelse(LITTER=="BLK", "Control",
                       ifelse(LITTER=="LNU", "N UnSc",
                              ifelse(LITTER=="LNS", "N Sc", 
                                     ifelse(LITTER=="LPU", "P UnSc", "P Sc")))),
         CUM_RESP=CUM_RESP_ADJ)%>%
  #filter(BATCH==2)%>%
  select(TYPE, STATUS, LITTER, TIME, CUM_RESP)#Load the data

#Decay Curves
TOAST_blank<-TOAST%>%
  filter(LITTER=="Control")%>%
  group_by(STATUS,TYPE, TIME)#%>%
#summarise_at(vars(CUM_RESP), funs(mean))

plot(TOAST_blank$TIME,TOAST_blank$CUM_RESP)
plot(TOAST_nUnSc$TIME,TOAST_nUnSc$CUM_RESP)
plot(TOAST_pUnSc$TIME,TOAST_pUnSc$CUM_RESP)
plot(TOAST_nSc$TIME,TOAST_nSc$CUM_RESP)
plot(TOAST_pSc$TIME,TOAST_pSc$CUM_RESP)

##----
fit_1<-TOAST_blank%>%
  filter(STATUS=="Burnt", TYPE=="Native")
fit_a<-nls(CUM_RESP~C*(1-exp(k*TIME)),data=fit_1,
           algorithm="port",
           start=c(C=10,k=-1),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
a<-summary(fit_a)
fit_1$pred <- predict(fit_a)
plot(fit_1$TIME,fit_1$CUM_RESP)
lines(fit_1$TIME,fit_1$pred, col="blue",lty=2)

##----
fit_2<-TOAST_blank%>%
  filter(STATUS=="Burnt", TYPE=="Pine")
fit_b<-nls(CUM_RESP~C*(1-exp(k*TIME)),data=fit_2,
           algorithm="port",
           start=c(C=10,k=-1),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
b<-summary(fit_b)
fit_2$pred <- predict(fit_b)
plot(fit_2$TIME,fit_2$CUM_RESP)
lines(fit_2$TIME,fit_2$pred, col="blue",lty=2)

##----
fit_3<-TOAST_blank%>%
  filter(STATUS=="Unburnt", TYPE=="Native")
fit_c<-nls2(CUM_RESP~C*(1-exp(k*TIME)),data=fit_3,
           algorithm="port",
           start=c(C=4,k=-0.3),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
c<-summary(fit_c)
fit_3$pred <- predict(fit_c)
plot(fit_3$TIME,fit_3$CUM_RESP)
lines(fit_3$TIME,fit_3$pred, col="blue",lty=2)

##----
fit_4<-TOAST_blank%>%
  filter(STATUS=="Unburnt", TYPE=="Pine")
fit_d<-nls2(CUM_RESP~C*(1-exp(k*TIME)),data=fit_4,
            algorithm="port",
            start=c(C=4,k=-0.3),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
d<-summary(fit_d)
fit_4$pred <- predict(fit_d)
plot(fit_4$TIME,fit_4$CUM_RESP)
lines(fit_4$TIME,fit_4$pred, col="blue",lty=2)


TOAST_blank_ave<-data.frame(id=c(rep("Blank",4)),
                            C=c(a$parameters[1,1],b$parameters[1,1],c$parameters[1,1],d$parameters[1,1]),
                            k=c(a$parameters[2,1],b$parameters[2,1],c$parameters[2,1],d$parameters[2,1]))

##--------------------Pine Scorched litter
TOAST_pSc<-TOAST%>%
  filter(LITTER=="P Sc")%>%
  group_by(STATUS,TYPE, TIME)

##----
fit_1<-TOAST_pSc%>%
  filter(STATUS=="Burnt", TYPE=="Native")
fit_a<-nls(CUM_RESP~C*(1-exp(k*TIME)),data=fit_1,
           algorithm="port",
           start=c(C=10,k=-1),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
a<-summary(fit_a)
fit_1$pred <- predict(fit_a)
plot(fit_1$TIME,fit_1$CUM_RESP)
lines(fit_1$TIME,fit_1$pred, col="blue",lty=2)

##----
fit_2<-TOAST_pSc%>%
  filter(STATUS=="Burnt", TYPE=="Pine")
fit_b<-nls(CUM_RESP~C*(1-exp(k*TIME)),data=fit_2,
           algorithm="port",
           start=c(C=10,k=-1),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
b<-summary(fit_b)
fit_2$pred <- predict(fit_b)
plot(fit_2$TIME,fit_2$CUM_RESP)
lines(fit_2$TIME,fit_2$pred, col="blue",lty=2)

##----
fit_3<-TOAST_pSc%>%
  filter(STATUS=="Unburnt", TYPE=="Native")
fit_c<-nls(CUM_RESP~C*(1-exp(k*TIME)),data=fit_3,
           algorithm="port",
           start=c(C=10,k=-1),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
c<-summary(fit_c)
fit_3$pred <- predict(fit_c)
plot(fit_3$TIME,fit_3$CUM_RESP)
lines(fit_3$TIME,fit_3$pred, col="blue",lty=2)

##----
fit_4<-TOAST_pSc%>%
  filter(STATUS=="Unburnt", TYPE=="Pine")
fit_d<-nls(CUM_RESP~C*(1-exp(k*TIME)),data=fit_4,
           algorithm="port",
           start=c(C=10,k=-1),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
d<-summary(fit_d)
fit_4$pred <- predict(fit_d)
plot(fit_4$TIME,fit_4$CUM_RESP)
lines(fit_4$TIME,fit_4$pred, col="blue",lty=2)


TOAST_pSc_ave<-data.frame(id=c(rep("P Sc",4)),
                          C=c(a$parameters[1,1],b$parameters[1,1],c$parameters[1,1],d$parameters[1,1]),
                          k=c(a$parameters[2,1],b$parameters[2,1],c$parameters[2,1],d$parameters[2,1]))

##--------------------Pine UnScorched litter
TOAST_pUnSc<-TOAST%>%
  filter(LITTER=="P UnSc")%>%
  group_by(STATUS,TYPE, TIME)

##----
fit_1<-TOAST_pUnSc%>%
  filter(STATUS=="Burnt", TYPE=="Native")
fit_a<-nls2(CUM_RESP~C*(1-exp(k*TIME)),data=fit_1,
            algorithm="port",
            start=c(C=12,k=-0.2),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
a<-summary(fit_a)
fit_1$pred <- predict(fit_a)
plot(fit_1$TIME,fit_1$CUM_RESP)
lines(fit_1$TIME,fit_1$pred, col="blue",lty=2)

##----
fit_2<-TOAST_pUnSc%>%
  filter(STATUS=="Burnt", TYPE=="Pine")
fit_b<-nls(CUM_RESP~C*(1-exp(k*TIME)),data=fit_2,
           algorithm="port",
           start=c(C=10,k=-1),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
b<-summary(fit_b)
fit_2$pred <- predict(fit_b)
plot(fit_2$TIME,fit_2$CUM_RESP)
lines(fit_2$TIME,fit_2$pred, col="blue",lty=2)

##----
fit_3<-TOAST_pUnSc%>%
  filter(STATUS=="Unburnt", TYPE=="Native")
fit_c<-nls2(CUM_RESP~C*(1-exp(k*TIME)),data=fit_3,
           algorithm="port",
           start=c(C=10,k=-0.3),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
c<-summary(fit_c)
fit_3$pred <- predict(fit_c)
plot(fit_3$TIME,fit_3$CUM_RESP)
lines(fit_3$TIME,fit_3$pred, col="blue",lty=2)

##----
fit_4<-TOAST_pUnSc%>%
  filter(STATUS=="Unburnt", TYPE=="Pine")
fit_d<-nls2(CUM_RESP~C*(1-exp(k*TIME)),data=fit_4,
           algorithm="port",
           start=c(C=12,k=-0.27),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
d<-summary(fit_d)
fit_4$pred <- predict(fit_d)
plot(fit_4$TIME,fit_4$CUM_RESP)
lines(fit_4$TIME,fit_4$pred, col="blue",lty=2)


TOAST_pUnSc_ave<-data.frame(id=c(rep("P UnSc",4)),
                            C=c(a$parameters[1,1],b$parameters[1,1],c$parameters[1,1],d$parameters[1,1]),
                            k=c(a$parameters[2,1],b$parameters[2,1],c$parameters[2,1],d$parameters[2,1]))

##--------------------Native Scorched litter
TOAST_nSc<-TOAST%>%
  filter(LITTER=="N Sc")%>%
  group_by(STATUS,TYPE, TIME)

##----
fit_1<-TOAST_nSc%>%
  filter(STATUS=="Burnt", TYPE=="Native")
fit_a<-nls(CUM_RESP~C*(1-exp(k*TIME)),data=fit_1,
           algorithm="port",
           start=c(C=10,k=-1),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
a<-summary(fit_a)
fit_1$pred <- predict(fit_a)
plot(fit_1$TIME,fit_1$CUM_RESP)
lines(fit_1$TIME,fit_1$pred, col="blue",lty=2)

##----
fit_2<-TOAST_nSc%>%
  filter(STATUS=="Burnt", TYPE=="Pine")
fit_b<-nls2(CUM_RESP~C*(1-exp(k*TIME)),data=fit_2,
            algorithm="port",
            start=c(C=12,k=-0.35),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
b<-summary(fit_b)
fit_2$pred <- predict(fit_b)
plot(fit_2$TIME,fit_2$CUM_RESP)
lines(fit_2$TIME,fit_2$pred, col="blue",lty=2)

##----
fit_3<-TOAST_nSc%>%
  filter(STATUS=="Unburnt", TYPE=="Native")
fit_c<-nls2(CUM_RESP~C*(1-exp(k*TIME)),data=fit_3,
           algorithm="port",
           start=c(C=9,k=-0.25),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
c<-summary(fit_c)
fit_3$pred <- predict(fit_c)
plot(fit_3$TIME,fit_3$CUM_RESP)
lines(fit_3$TIME,fit_3$pred, col="blue",lty=2)

##----
fit_4<-TOAST_nSc%>%
  filter(STATUS=="Unburnt", TYPE=="Pine")
fit_d<-nls2(CUM_RESP~C*(1-exp(k*TIME)),data=fit_4,
           algorithm="port",
           start=c(C=11,k=-0.4),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
d<-summary(fit_d)
fit_4$pred <- predict(fit_d)
plot(fit_4$TIME,fit_4$CUM_RESP)
lines(fit_4$TIME,fit_4$pred, col="blue",lty=2)


TOAST_nSc_ave<-data.frame(id=c(rep("N Sc",4)),
                          C=c(a$parameters[1,1],b$parameters[1,1],c$parameters[1,1],d$parameters[1,1]),
                          k=c(a$parameters[2,1],b$parameters[2,1],c$parameters[2,1],d$parameters[2,1]))

##--------------------Native Uncorched litter
TOAST_nUnSc<-TOAST%>%
  filter(LITTER=="N UnSc")%>%
  group_by(STATUS,TYPE, TIME)

##----
fit_1<-TOAST_nUnSc%>%
  filter(STATUS=="Burnt", TYPE=="Native")
fit_a<-nls(CUM_RESP~C*(1-exp(k*TIME)),data=fit_1,
           algorithm="port",
           start=c(C=10,k=-1),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
a<-summary(fit_a)
fit_1$pred <- predict(fit_a)
plot(fit_1$TIME,fit_1$CUM_RESP)
lines(fit_1$TIME,fit_1$pred, col="blue",lty=2)

##----
fit_2<-TOAST_nUnSc%>%
  filter(STATUS=="Burnt", TYPE=="Pine")
fit_b<-nls(CUM_RESP~C*(1-exp(k*TIME)),data=fit_2,
           algorithm="port",
           start=c(C=10,k=-1),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
b<-summary(fit_b)
fit_2$pred <- predict(fit_b)
plot(fit_2$TIME,fit_2$CUM_RESP)
lines(fit_2$TIME,fit_2$pred, col="blue",lty=2)

##----
fit_3<-TOAST_nUnSc%>%
  filter(STATUS=="Unburnt", TYPE=="Native")
fit_c<-nls2(CUM_RESP~C*(1-exp(k*TIME)),data=fit_3,
           algorithm="port",
           start=c(C=10,k=-0.3),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
c<-summary(fit_c)
fit_3$pred <- predict(fit_c)
plot(fit_3$TIME,fit_3$CUM_RESP)
lines(fit_3$TIME,fit_3$pred, col="blue",lty=2)

##----
fit_4<-TOAST_nUnSc%>%
  filter(STATUS=="Unburnt", TYPE=="Pine")
fit_d<-nls(CUM_RESP~C*(1-exp(k*TIME)),data=fit_4,
           algorithm="port",
           start=c(C=10,k=-1),lower=c(C=0,k=-Inf), upper=c(C=Inf,k=0))
d<-summary(fit_d)
fit_4$pred <- predict(fit_d)
plot(fit_4$TIME,fit_4$CUM_RESP)
lines(fit_4$TIME,fit_4$pred, col="blue",lty=2)


TOAST_nUnSc_ave<-data.frame(id=c(rep("N UnSc",4)),
                            C=c(a$parameters[1,1],b$parameters[1,1],c$parameters[1,1],d$parameters[1,1]),
                            k=c(a$parameters[2,1],b$parameters[2,1],c$parameters[2,1],d$parameters[2,1]))

Curves<-rbind(TOAST_blank_ave,TOAST_pSc_ave,TOAST_pUnSc_ave,TOAST_nSc_ave,TOAST_nUnSc_ave)

##-----Stats

stats<-Curves%>%
  group_by(id)

anova(aov(C~id, data=stats)) #**
HSD.test(aov(C~id, data=stats), c("id"),
         console=TRUE,alpha=0.05)

aggregate(C~id, data=stats, FUN="sd")

anova(aov(k~id, data=stats)) #ns
HSD.test(aov(k~id, data=stats), c("id"),
         console=TRUE,alpha=0.05)

aggregate(k~id, data=stats, FUN="sd")

####------------Plot
fun_control <- function(x) {
  8.067852*(1-exp(-0.2238511*x))
}

fun_NUnSc <- function(x) {
  9.336568*(1-exp(-0.2926802*x))
}

fun_NSc <- function(x) {
  15.865345*(1-exp(-0.2815181*x))
}

fun_PUnSc <- function(x) {
  13.363417*(1-exp(-0.2842943*x))
}

fun_PSc <- function(x) {
  28.000682*(1-exp(-0.2615819*x))
}

TOAST_lines<-ggplot(data.frame(x=c(0,6)),aes(x=x))+
  stat_function(fun=fun_control, colour="grey")+
  stat_function(fun=fun_NUnSc, colour="lightblue")+
  stat_function(fun=fun_NSc, colour="blue")+
  stat_function(fun=fun_PUnSc, colour="orange")+
  stat_function(fun=fun_PSc, colour="red")+
  scale_y_continuous(limits=c(0,40))+
  ggsave("TOAST_a7.svg", width=200, height=100, units="mm")

TOAST%>%
  ggplot(aes(x=TIME,y=CUM_RESP))+
  geom_point(aes(colour=LITTER))+
  scale_colour_manual(values=c("grey", "blue", "lightblue", "red", "orange"))+
  xlab("Time (days)") +
  ylab(expression(paste("Cumulative Respiration (",mu,"gCO2-C/gOC)"))) +
  scale_y_continuous(limits=c(0,40))+
  ggsave("TOAST_b7.svg", width=200, height=100, units="mm")