#####################################################
# Marco FW Gauger
# 2021/02/23
#
# This is the script to obtain table 2 of the manuscript submitted to PeerJ Gauger et al 2021 Diel influences on bottlenose dolphin acoustic detection in a coastal lagoon in the southwestern Gulf of California
#####################################################
ifelse("mgcv" %in% rownames(installed.packages()) == FALSE, 
       install.packages("mgcv", dependencies = T),
       ("mgcv is already installed"))
ifelse("bbmle" %in% rownames(installed.packages()) == FALSE, 
       install.packages("bbmle", dependencies = T),
       ("bbmle is already installed"))
ifelse("chron" %in% rownames(installed.packages()) == FALSE, 
       install.packages("chron", dependencies = T),
       ("chron is already installed"))

require(mgcv)
Sys.setlocale("LC_ALL","English")
Sys.setenv(TZ='UTC')

########
analysed_data_hour <- read.csv("Supplemental Data S1.csv", sep=";")
analysed_data_hour$hour[1]
analysed_data_hour$h[1]
analysed_data_hour$datehour_lt <- as.POSIXct(paste0(analysed_data_hour$datehour,":00"), format="%d.%m.%Y %H:%M:%S", tz="MST")
chron::hours(analysed_data_hour$datehour_lt)[1]

final <-
  gam(dp10m_ELAP1~ s(lunar_phase, bs="cc", k=6) +
        te(derivate_tide, tide, bs=c("tp","tp"),
           k=c(10,10)) + 
        s(SST, k=10, bs="tp", m=2) +
        s(h, bs="cc") + 
        effort +
        s(depth, k=4) +
        s(distance_BLAP, k=4),
      data=analysed_data_hour, method="REML",family="poisson",
      knots=list(h = c(0,24)))
summary(final)
AIC(final)

final_minus_1 <-
  gam(dp10m_ELAP1~ #s(lunar_phase, bs="cc", k=6) +
        te(derivate_tide, tide, bs=c("tp","tp"),
           k=c(10,10)) + 
        s(SST, k=10, bs="tp") +
        s(h, bs="cc") + 
        effort +
        s(depth, k=4) +
        s(distance_BLAP, k=4),
      data=analysed_data_hour, method="REML",family="poisson",
      knots=list(h = c(0,24)),
      sp=final$sp[-1])
final_minus_2 <-
  gam(dp10m_ELAP1~ s(lunar_phase, bs="cc", k=6) +
        #te(derivate_tide, tide, bs=c("tp","tp"),
        #   k=c(10,10)) + 
        s(SST, k=10, bs="tp") +
        s(h, bs="cc") + 
        effort +
        s(depth, k=4) +
        s(distance_BLAP, k=4),
      data=analysed_data_hour, method="REML",family="poisson",
      knots=list(h = c(0,24)),
      sp=final$sp[-c(2,3)])
final_minus_3 <-
  gam(dp10m_ELAP1~ s(lunar_phase, bs="cc", k=6) +
        te(derivate_tide, tide, bs=c("tp","tp"),
           k=c(10,10)) + 
        #s(SST, bs="tp") +
        s(h, k=10, bs="cc") + 
        effort +
        s(depth, k=4) +
        s(distance_BLAP, k=4),
      data=analysed_data_hour, method="REML",family="poisson",
      knots=list(h = c(0,24)),
      sp=final$sp[-3])
final_minus_4 <-
  gam(dp10m_ELAP1~ s(lunar_phase, bs="cc", k=6) +
        te(derivate_tide, tide, bs=c("tp","tp"),
           k=c(10,10)) + 
        s(SST, bs="tp") +
        #s(h, k=10, bs="cc") + 
        effort +
        s(depth, k=4) +
        s(distance_BLAP, k=4),
      data=analysed_data_hour, method="REML",family="poisson",
      knots=list(h = c(0,24)),
      sp=final$sp[-4])
final_minus_5 <-
  gam(dp10m_ELAP1~ s(lunar_phase, bs="cc", k=6) +
        te(derivate_tide, tide, bs=c("tp","tp"),
           k=c(10,10)) + 
        s(SST, k=10, bs="tp") +
        s(h, bs="cc") + 
        #effort +
        s(depth, k=4) +
        s(distance_BLAP, k=4),
      data=analysed_data_hour, method="REML",family="poisson",
      knots=list(h = c(0,24)),
      sp=final$sp)
final_minus_6 <-
  gam(dp10m_ELAP1~ s(lunar_phase, bs="cc", k=6) +
        te(derivate_tide, tide, bs=c("tp","tp"),
           k=c(10,10)) + 
        s(SST, k=10, bs="tp") +
        s(h, bs="cc") + 
        effort +
        #s(depth, k=4) +
        s(distance_BLAP, k=4),
      data=analysed_data_hour, method="REML",family="poisson",
      knots=list(h = c(0,24)),
      sp=final$sp[-6])
final_minus_7 <-
  gam(dp10m_ELAP1~ s(lunar_phase, bs="cc", k=6) +
        te(derivate_tide, tide, bs=c("tp","tp"),
           k=c(10,10)) + 
        s(SST, k=10, bs="tp") +
        s(h, bs="cc") + 
        effort +
        s(depth, k=4),# +
        #s(distance_BLAP, k=4),
      data=analysed_data_hour, method="REML",family="poisson",
      knots=list(h = c(0,24)),
      sp=final$sp[-7])

AIC(final, final_minus_1, final_minus_2, final_minus_3, final_minus_4, final_minus_5, final_minus_6, final_minus_7)

bbmle::AICtab(final, final_minus_1, final_minus_2, final_minus_3, final_minus_4, final_minus_5, final_minus_6, final_minus_7, weights=T)

(dev_expl = round(cbind(c(summary(final)$dev.expl,
summary(final_minus_1)$dev.expl,
summary(final_minus_2)$dev.expl,
summary(final_minus_3)$dev.expl,
summary(final_minus_4)$dev.expl,
summary(final_minus_6)$dev.expl,
summary(final_minus_7)$dev.expl,
summary(final_minus_5)$dev.expl)),digits=3))

rbind(summary(final)$p.table[1,c(1,4)],summary(final)$s.table[,c(1,4)],summary(final)$p.table[2,c(1,4)])
dev_expl
c("",1-dev_expl[2:8]/dev_expl[1])
write.table(dev_expl, "clipboard", sep="\t")
write.table(c("",1-dev_expl[2:8]/dev_expl[1]), "clipboard", sep="\t")