#R code to fit Michaelis-Menton Curves to song discrimination data
#This code opens a single randomization dataset (1000 were used in the paper) in order to demonstrate
#how to fit the Michaelis-Menton curves using our modelling approach.
##############################################################################################
##############################################################################################
##############################################################################################
###MODELS FORCED THROUGH ORIGIN (e.g. D vs genetic distance)

#STEP0: set working directory to input file
   setwd("/storage/SAMBA/home/0_PROJECTS/Playback_Experiments/Manuscript/9July2019_PLOSBIOLOGY_reviews/28July2019")



#STEP1: setup dataset
   library(qpcR)
   DATA <- read.table("DATASET.txt", header=TRUE) # a single randomized dataset 

   DATA_ALL <-  DATA[,c("LAT", "GEN_DIST", "GEN_DIST", "GEN_DIST", "D", "Contrast_ID", "Territory")]
   colnames(DATA_ALL) <- c("LAT", "TERRITORY", "GROUP", "TIME", "D", "Contrast_ID", "Territory")

   DATA_ALL[DATA[,"Territory"] == "ABSCENCE" & DATA[,"LAT"] == "TEMPERATE", "TERRITORY"] <- "ABSCENCE_TEMP"
   DATA_ALL[DATA[,"Territory"] == "ABSCENCE" & DATA[,"LAT"] == "AMAZON", "TERRITORY"] <- "ABSCENCE_AMAZ"
   DATA_ALL[DATA[,"Territory"] == "PRESENCE" & DATA[,"LAT"] == "TEMPERATE", "TERRITORY"] <- "PRESENCE_TEMP"
   DATA_ALL[DATA[,"Territory"] == "PRESENCE" & DATA[,"LAT"] == "AMAZON", "TERRITORY"] <- "PRESENCE_AMAZ"

   DATA_ALL[DATA[,"GROUP"] == "LEARNED" & DATA[,"LAT"] == "TEMPERATE", "GROUP"] <- "LEARNED_TEMP"
   DATA_ALL[DATA[,"GROUP"] == "LEARNED" & DATA[,"LAT"] == "AMAZON", "GROUP"] <- "LEARNED_AMAZ"
   DATA_ALL[DATA[,"GROUP"] != "LEARNED" & DATA[,"LAT"] == "TEMPERATE", "GROUP"] <- "INNATE_TEMP"
   DATA_ALL[DATA[,"GROUP"] != "LEARNED" & DATA[,"LAT"] == "AMAZON", "GROUP"] <- "INNATE_AMAZ"

   mm_i <- as.data.frame(DATA_ALL)


#STEP2: setup OUTPUT table to place output into
   OUTPUT <- matrix(data = NA, nrow = 1, ncol = 40)
   colnames(OUTPUT) <- c("ALL_AIC", "ALL_beta", "ALL_z", "LAT_AIC","LAT_beta_AMAZON", "LAT_beta_HIGH", "LAT_z", "GROUP_AIC", "GROUP_beta_LEARN", "GROUP_beta_INNATE", "GROUP_z", "LAT_GROUP_AIC", "LAT_GROUP_beta_AMAZON_LEARN", "LAT_GROUP_beta_TEMP_LEARN", "LAT_GROUP_beta_AMAZON_INNATE", "LAT_GROUP_beta_TEMP_INNATE", "LAT_GROUP_z", "TERRITORY_AIC", "TERRITORY_beta_NON", "TERRITORY_beta_NOTALLYEAR", "TERRITORY_beta_ALLYEAR", "TERRITORY_z", "LAT_TERRITORY_AIC", "LAT_TERRITORY_beta_AMAZON_NOTALLYEAR", "LAT_TERRITORY_beta_TEMP_NOTALLYEAR",  "LAT_TERRITORY_beta_AMAZON_ALLYEAR", "LAT_TERRITORY_beta_TEMP_ALLYEAR","LAT_TERRITORY_z", "beta_2trop_1temp_TEMP","beta_NOTALLYEAR_2trop_1temp_AMAZ", "beta_ALLYEAR_2trop_1temp_AMAZ", "2trop_1temp_AIC", "2trop_1temp_z", "Akaike_Weight_ALL", "Akaike_Weight_LAT", "Akaike_Weight_GROUP", "Akaike_Weight_LAT_GROUP",  "Akaike_Weight_TERRITORY", "Akaike_Weight_LAT_TERRITORY", "Akaike_Weight_2trop_1temp")


#STEP3: run model fits
      START = 2

      model.ALL <- nls(D ~ 
              (1 * TIME)/(beta + TIME), data = mm_i,
                 start = list(beta = START),  nls.control(maxiter = 500))

      model.LAT <- nls(D ~ 
              as.numeric(LAT == "TEMPERATE")*(1 * TIME)/(beta_TEMP + TIME)
            + as.numeric(LAT == "AMAZON")*(1 * TIME)/(beta_AMAZ + TIME), data = mm_i, 
                 start = list(beta_TEMP = START, 
                              beta_AMAZ = START), lower = c(0,0), upper = c(20,20), algorithm = "port",   nls.control(maxiter = 500))

      model.GROUP <- nls(D ~ 
              as.numeric(GROUP == "LEARNED_TEMP" | GROUP == "LEARNED_AMAZ")*(1 * TIME)/(beta_LEARNED + TIME)
            + as.numeric(GROUP == "INNATE_TEMP" | GROUP == "INNATE_AMAZ")*(1 * TIME)/(beta_INNATE + TIME), data = mm_i, 
                 start = list(beta_LEARNED = START, 
                              beta_INNATE = START), lower = c(0,0), upper = c(20,20), algorithm = "port",   nls.control(maxiter = 500))

      model.LAT_GROUP <- nls(D ~ 
              as.numeric(GROUP == "LEARNED_TEMP")*(1 * TIME)/(beta_LEARNED_TEMP + TIME)
            + as.numeric(GROUP == "LEARNED_AMAZ")*(1 * TIME)/(beta_LEARNED_AMAZ + TIME)
            + as.numeric(GROUP == "INNATE_TEMP")*(1 * TIME)/(beta_INNATE_TEMP + TIME)
            + as.numeric(GROUP == "INNATE_AMAZ")*(1 * TIME)/(beta_INNATE_AMAZ + TIME), data = mm_i, 
                 start = list(beta_LEARNED_TEMP = START, 
                              beta_LEARNED_AMAZ = START,
                              beta_INNATE_TEMP = START,
                              beta_INNATE_AMAZ = START),
                 lower = c(0,0,0,0), upper = c(20,20,20,20), algorithm = "port",  nls.control(maxiter = 500))




      model.TERRITORY <- nls(D ~ 
              as.numeric(TERRITORY == "ABSCENCE_TEMP" | TERRITORY == "ABSCENCE_AMAZ")*(1 * TIME)/(beta_NOTALLYEAR_TERRITORIAL + TIME)            
            + as.numeric(TERRITORY == "PRESENCE_TEMP" | TERRITORY == "PRESENCE_AMAZ")*(1 * TIME)/(beta_ALLYEAR_TERRITORIAL + TIME), data = mm_i, 
                 start = list(beta_NOTALLYEAR_TERRITORIAL = START, 
                              beta_ALLYEAR_TERRITORIAL = START), lower = c(0,0), upper = c(20,20), algorithm = "port",   nls.control(maxiter = 500))


      model.2trop_1temp <- nls(D ~ 
              as.numeric(TERRITORY == "ABSCENCE_TEMP" | TERRITORY == "PRESENCE_TEMP")*(1 * TIME)/(beta_2trop_1temp_TEMP + TIME)
            + as.numeric(TERRITORY == "ABSCENCE_AMAZ")*(1 * TIME)/(beta_NOTALLYEAR_2trop_1temp_AMAZ + TIME)
            + as.numeric(TERRITORY == "PRESENCE_AMAZ")*(1 * TIME)/(beta_ALLYEAR_2trop_1temp_AMAZ + TIME), data = mm_i, 
                 start = list(beta_2trop_1temp_TEMP = START, 
                              beta_NOTALLYEAR_2trop_1temp_AMAZ = START,
                              beta_ALLYEAR_2trop_1temp_AMAZ = START), lower = c(0,0,0), upper = c(20,20,20), algorithm = "port",   nls.control(maxiter = 500))

         OUTPUT[1,"ALL_AIC"] <- AIC(model.ALL)
         OUTPUT[1,"ALL_beta"] <- coef(model.ALL)
	
         OUTPUT[1,"LAT_AIC"] <-  AIC(model.LAT)
         OUTPUT[1,"LAT_beta_HIGH"] <- coef(model.LAT)[[1]] 
         OUTPUT[1,"LAT_beta_AMAZON"] <- coef(model.LAT)[[2]]

         OUTPUT[1,"GROUP_AIC"] <-  AIC(model.GROUP)
         OUTPUT[1,"GROUP_beta_LEARN"] <-  coef(model.GROUP)[[1]]
         OUTPUT[1,"GROUP_beta_INNATE"] <-  coef(model.GROUP)[[2]]

         OUTPUT[1,"LAT_GROUP_AIC"] <-  AIC(model.LAT_GROUP)
         OUTPUT[1,"LAT_GROUP_beta_AMAZON_LEARN"] <-  coef(model.LAT_GROUP)[[2]]
         OUTPUT[1,"LAT_GROUP_beta_TEMP_LEARN"] <-  coef(model.LAT_GROUP)[[1]]
         OUTPUT[1,"LAT_GROUP_beta_AMAZON_INNATE"] <-  coef(model.LAT_GROUP)[[4]]
         OUTPUT[1,"LAT_GROUP_beta_TEMP_INNATE"] <-  coef(model.LAT_GROUP)[[3]]

         OUTPUT[1,"TERRITORY_AIC"] <-  AIC(model.TERRITORY)
         OUTPUT[1,"TERRITORY_beta_NOTALLYEAR"] <-  coef(model.TERRITORY)[[1]]
         OUTPUT[1,"TERRITORY_beta_ALLYEAR"] <-  coef(model.TERRITORY)[[2]]

         OUTPUT[1,"2trop_1temp_AIC"] <-  AIC(model.2trop_1temp)
         OUTPUT[1,"beta_NOTALLYEAR_2trop_1temp_AMAZ"] <-  coef(model.2trop_1temp)[[2]]
         OUTPUT[1,"beta_2trop_1temp_TEMP"] <-  coef(model.2trop_1temp)[[1]]
         OUTPUT[1,"beta_ALLYEAR_2trop_1temp_AMAZ"] <-  coef(model.2trop_1temp)[[3]]
	}

      OUTPUT[1,c("Akaike_Weight_ALL", "Akaike_Weight_LAT", "Akaike_Weight_GROUP", "Akaike_Weight_LAT_GROUP", "Akaike_Weight_TERRITORY", "Akaike_Weight_2trop_1temp")] <-  akaike.weights(c(AIC(model.ALL), AIC(model.LAT), AIC(model.GROUP), AIC(model.LAT_GROUP), AIC(model.TERRITORY), AIC(model.2trop_1temp)))$weights


#STEP4: CONSISE SUMMARY
#Akaike Weight
   #MODEL 1
   OUTPUT[,"Akaike_Weight_ALL"]

   #MODEL 2
   OUTPUT[,"Akaike_Weight_LAT"]

   #MODEL 3
   OUTPUT[,"Akaike_Weight_GROUP"]

   #MODEL 4
   OUTPUT[,"Akaike_Weight_LAT_GROUP"]

   #MODEL 5
   OUTPUT[,"Akaike_Weight_TERRITORY"]

   #MODEL 6
   OUTPUT[,"Akaike_Weight_2trop_1temp"]


#beta
   #MODEL 1
   OUTPUT[,"ALL_beta"]

   #MODEL 2
   OUTPUT[,"LAT_beta_HIGH"]
   OUTPUT[,"LAT_beta_AMAZON"]

   #MODEL 3

   OUTPUT[,"GROUP_beta_LEARN"]
   OUTPUT[,"GROUP_beta_INNATE"]

   #MODEL 4
   OUTPUT[,"LAT_GROUP_beta_TEMP_LEARN"]
   OUTPUT[,"LAT_GROUP_beta_TEMP_INNATE"]
   OUTPUT[,"LAT_GROUP_beta_AMAZON_LEARN"]
   OUTPUT[,"LAT_GROUP_beta_AMAZON_INNATE"]

   #MODEL 5
   OUTPUT[,"TERRITORY_beta_NOTALLYEAR"]
   OUTPUT[,"TERRITORY_beta_ALLYEAR"]

   #MODEL 6
   OUTPUT[,"beta_2trop_1temp_TEMP"]
   OUTPUT[,"beta_NOTALLYEAR_2trop_1temp_AMAZ"]
   OUTPUT[,"beta_ALLYEAR_2trop_1temp_AMAZ"]


##############################################################################################
##############################################################################################
##############################################################################################
###MODELS NOT FORCED THROUGH ORIGIN (e.g. D vs song length divergence)

#STEP0: set working directory to input file
   setwd("/storage/SAMBA/home/0_PROJECTS/Playback_Experiments/Manuscript/9July2019_PLOSBIOLOGY_reviews/28July2019")


#STEP1: setup dataset
   library(qpcR)
   DATA <- read.table("DATASET.txt", header=TRUE)

   DATA_ALL <-  DATA[,c("LAT", "GEN_DIST", "GEN_DIST", "LENGTH", "D", "Contrast_ID", "Territory")]
   colnames(DATA_ALL) <- c("LAT", "TERRITORY", "GROUP", "SONG_LENGTH", "D", "Contrast_ID", "Territory")

   DATA_ALL[DATA[,"Territory"] == "ABSCENCE" & DATA[,"LAT"] == "TEMPERATE", "TERRITORY"] <- "ABSCENCE_TEMP"
   DATA_ALL[DATA[,"Territory"] == "ABSCENCE" & DATA[,"LAT"] == "AMAZON", "TERRITORY"] <- "ABSCENCE_AMAZ"
   DATA_ALL[DATA[,"Territory"] == "PRESENCE" & DATA[,"LAT"] == "TEMPERATE", "TERRITORY"] <- "PRESENCE_TEMP"
   DATA_ALL[DATA[,"Territory"] == "PRESENCE" & DATA[,"LAT"] == "AMAZON", "TERRITORY"] <- "PRESENCE_AMAZ"

   DATA_ALL[DATA[,"GROUP"] == "LEARNED" & DATA[,"LAT"] == "TEMPERATE", "GROUP"] <- "LEARNED_TEMP"
   DATA_ALL[DATA[,"GROUP"] == "LEARNED" & DATA[,"LAT"] == "AMAZON", "GROUP"] <- "LEARNED_AMAZ"
   DATA_ALL[DATA[,"GROUP"] != "LEARNED" & DATA[,"LAT"] == "TEMPERATE", "GROUP"] <- "INNATE_TEMP"
   DATA_ALL[DATA[,"GROUP"] != "LEARNED" & DATA[,"LAT"] == "AMAZON", "GROUP"] <- "INNATE_AMAZ"

   mm_i <- as.data.frame(DATA_ALL)


#STEP2: setup OUTPUT table to place output into
   OUTPUT <- matrix(data = NA, nrow = 1, ncol = 40)
   colnames(OUTPUT) <- c("ALL_AIC", "ALL_beta", "ALL_z", "LAT_AIC","LAT_beta_AMAZON", "LAT_beta_HIGH", "LAT_z", "GROUP_AIC", "GROUP_beta_LEARN", "GROUP_beta_INNATE", "GROUP_z", "LAT_GROUP_AIC", "LAT_GROUP_beta_AMAZON_LEARN", "LAT_GROUP_beta_TEMP_LEARN", "LAT_GROUP_beta_AMAZON_INNATE", "LAT_GROUP_beta_TEMP_INNATE", "LAT_GROUP_z", "TERRITORY_AIC", "TERRITORY_beta_NON", "TERRITORY_beta_NOTALLYEAR", "TERRITORY_beta_ALLYEAR", "TERRITORY_z", "LAT_TERRITORY_AIC", "LAT_TERRITORY_beta_AMAZON_NOTALLYEAR", "LAT_TERRITORY_beta_TEMP_NOTALLYEAR",  "LAT_TERRITORY_beta_AMAZON_ALLYEAR", "LAT_TERRITORY_beta_TEMP_ALLYEAR","LAT_TERRITORY_z", "beta_2trop_1temp_TEMP","beta_NOTALLYEAR_2trop_1temp_AMAZ", "beta_ALLYEAR_2trop_1temp_AMAZ", "2trop_1temp_AIC", "2trop_1temp_z", "Akaike_Weight_ALL", "Akaike_Weight_LAT", "Akaike_Weight_GROUP", "Akaike_Weight_LAT_GROUP",  "Akaike_Weight_TERRITORY", "Akaike_Weight_LAT_TERRITORY", "Akaike_Weight_2trop_1temp")



#STEP3: run model fits
      START = 2
      START_z = 0

      model.ALL <- nls(D ~ 
              (1 * SONG_LENGTH-z)/(beta-z + SONG_LENGTH-z), data = mm_i,
                 start = list(beta = START, z = START_z),  nls.control(maxiter = 500))

      model.LAT <- nls(D ~ 
              as.numeric(LAT == "TEMPERATE")*(1 * SONG_LENGTH-z)/(beta_TEMP-z + SONG_LENGTH-z)
            + as.numeric(LAT == "AMAZON")*(1 * SONG_LENGTH-z)/(beta_AMAZ-z + SONG_LENGTH-z), data = mm_i, 
                 start = list(beta_TEMP = START, 
                              beta_AMAZ = START, z = START_z), lower = c(-5,-5,-10), upper = c(20,20,0), algorithm = "port",   nls.control(maxiter = 500))

      model.GROUP <- nls(D ~ 
              as.numeric(GROUP == "LEARNED_TEMP" | GROUP == "LEARNED_AMAZ")*(1 * SONG_LENGTH-z)/(beta_LEARNED-z + SONG_LENGTH-z)
            + as.numeric(GROUP == "INNATE_TEMP" | GROUP == "INNATE_AMAZ")*(1 * SONG_LENGTH-z)/(beta_INNATE-z + SONG_LENGTH-z), data = mm_i, 
                 start = list(beta_LEARNED = START, 
                              beta_INNATE = START, z = START_z), lower = c(-5,-5,-10), upper = c(20,20,0), algorithm = "port",   nls.control(maxiter = 500))

      model.LAT_GROUP <- nls(D ~ 
              as.numeric(GROUP == "LEARNED_TEMP")*(1 * SONG_LENGTH-z)/(beta_LEARNED_TEMP-z + SONG_LENGTH-z)
            + as.numeric(GROUP == "LEARNED_AMAZ")*(1 * SONG_LENGTH-z)/(beta_LEARNED_AMAZ-z + SONG_LENGTH-z)
            + as.numeric(GROUP == "INNATE_TEMP")*(1 * SONG_LENGTH-z)/(beta_INNATE_TEMP-z + SONG_LENGTH-z)
            + as.numeric(GROUP == "INNATE_AMAZ")*(1 * SONG_LENGTH-z)/(beta_INNATE_AMAZ-z + SONG_LENGTH-z), data = mm_i, 
                 start = list(beta_LEARNED_TEMP = START, 
                              beta_LEARNED_AMAZ = START,
                              beta_INNATE_TEMP = START,
                              beta_INNATE_AMAZ = START, z = START_z),
                 lower = c(-5,-5,-5,-5,-10), upper = c(20,20,20,20,0), algorithm = "port",  nls.control(maxiter = 500))




      model.TERRITORY <- nls(D ~ 
              as.numeric(TERRITORY == "ABSCENCE_TEMP" | TERRITORY == "ABSCENCE_AMAZ")*(1 * SONG_LENGTH-z)/(beta_NOTALLYEAR_TERRITORIAL-z + SONG_LENGTH-z)            
            + as.numeric(TERRITORY == "PRESENCE_TEMP" | TERRITORY == "PRESENCE_AMAZ")*(1 * SONG_LENGTH-z)/(beta_ALLYEAR_TERRITORIAL-z + SONG_LENGTH-z), data = mm_i, 
                 start = list(beta_NOTALLYEAR_TERRITORIAL = START, 
                              beta_ALLYEAR_TERRITORIAL = START, z = START_z), lower = c(-5,-5,-10), upper = c(20,20,0), algorithm = "port",   nls.control(maxiter = 500))


      model.2trop_1temp <- nls(D ~ 
              as.numeric(TERRITORY == "ABSCENCE_TEMP" | TERRITORY == "PRESENCE_TEMP")*(1 * SONG_LENGTH-z)/(beta_2trop_1temp_TEMP-z + SONG_LENGTH-z)
            + as.numeric(TERRITORY == "ABSCENCE_AMAZ")*(1 * SONG_LENGTH-z)/(beta_NOTALLYEAR_2trop_1temp_AMAZ-z + SONG_LENGTH-z)
            + as.numeric(TERRITORY == "PRESENCE_AMAZ")*(1 * SONG_LENGTH-z)/(beta_ALLYEAR_2trop_1temp_AMAZ-z + SONG_LENGTH-z), data = mm_i, 
                 start = list(beta_2trop_1temp_TEMP = START, 
                              beta_NOTALLYEAR_2trop_1temp_AMAZ = START,
                              beta_ALLYEAR_2trop_1temp_AMAZ = START, z = START_z), lower = c(-5,-5,-5,-10), upper = c(20,20,20,0), algorithm = "port",   nls.control(maxiter = 500))


         OUTPUT[1,"ALL_AIC"] <- AIC(model.ALL)
         OUTPUT[1,"ALL_beta"] <- coef(model.ALL)[[1]]
         OUTPUT[1,"ALL_z"] <- coef(model.ALL)[[1]]

         OUTPUT[1,"LAT_AIC"] <-  AIC(model.LAT)
         OUTPUT[1,"LAT_beta_HIGH"] <- coef(model.LAT)[[1]] 
         OUTPUT[1,"LAT_beta_AMAZON"] <- coef(model.LAT)[[2]]
         OUTPUT[1,"LAT_z"] <- coef(model.LAT)[[3]]

         OUTPUT[1,"GROUP_AIC"] <-  AIC(model.GROUP)
         OUTPUT[1,"GROUP_beta_LEARN"] <-  coef(model.GROUP)[[1]]
         OUTPUT[1,"GROUP_beta_INNATE"] <-  coef(model.GROUP)[[2]]
         OUTPUT[1,"GROUP_z"] <-  coef(model.GROUP)[[3]]

         OUTPUT[1,"LAT_GROUP_AIC"] <-  AIC(model.LAT_GROUP)
         OUTPUT[1,"LAT_GROUP_beta_AMAZON_LEARN"] <-  coef(model.LAT_GROUP)[[2]]
         OUTPUT[1,"LAT_GROUP_beta_TEMP_LEARN"] <-  coef(model.LAT_GROUP)[[1]]
         OUTPUT[1,"LAT_GROUP_beta_AMAZON_INNATE"] <-  coef(model.LAT_GROUP)[[4]]
         OUTPUT[1,"LAT_GROUP_beta_TEMP_INNATE"] <-  coef(model.LAT_GROUP)[[3]]
         OUTPUT[1,"LAT_GROUP_z"] <-  coef(model.LAT_GROUP)[[5]]

         OUTPUT[1,"TERRITORY_AIC"] <-  AIC(model.TERRITORY)
         OUTPUT[1,"TERRITORY_beta_NOTALLYEAR"] <-  coef(model.TERRITORY)[[1]]
         OUTPUT[1,"TERRITORY_beta_ALLYEAR"] <-  coef(model.TERRITORY)[[2]]
         OUTPUT[1,"TERRITORY_z"] <- coef(model.TERRITORY)[[3]]

         OUTPUT[1,"LAT_TERRITORY_AIC"] <-  AIC(model.LAT_TERRITORY)
         OUTPUT[1,"LAT_TERRITORY_beta_AMAZON_NOTALLYEAR"] <-  coef(model.LAT_TERRITORY)[[2]]
         OUTPUT[1,"LAT_TERRITORY_beta_TEMP_NOTALLYEAR"] <-  coef(model.LAT_TERRITORY)[[1]]
         OUTPUT[1,"LAT_TERRITORY_beta_AMAZON_ALLYEAR"] <-  coef(model.LAT_TERRITORY)[[4]]
         OUTPUT[1,"LAT_TERRITORY_beta_TEMP_ALLYEAR"] <-  coef(model.LAT_TERRITORY)[[3]]
         OUTPUT[1,"LAT_TERRITORY_z"] <- coef(model.LAT_TERRITORY)[[5]]

         OUTPUT[1,"2trop_1temp_AIC"] <-  AIC(model.2trop_1temp)
         OUTPUT[1,"beta_NOTALLYEAR_2trop_1temp_AMAZ"] <-  coef(model.2trop_1temp)[[2]]
         OUTPUT[1,"beta_2trop_1temp_TEMP"] <-  coef(model.2trop_1temp)[[1]]
         OUTPUT[1,"beta_ALLYEAR_2trop_1temp_AMAZ"] <-  coef(model.2trop_1temp)[[3]]
         OUTPUT[1,"2trop_1temp_z"] <- coef(model.2trop_1temp)[[4]]

      OUTPUT[1,c("Akaike_Weight_ALL", "Akaike_Weight_LAT", "Akaike_Weight_GROUP", "Akaike_Weight_LAT_GROUP", "Akaike_Weight_TERRITORY", "Akaike_Weight_2trop_1temp")] <-  akaike.weights(c(AIC(model.ALL), AIC(model.LAT), AIC(model.GROUP), AIC(model.LAT_GROUP), AIC(model.TERRITORY), AIC(model.2trop_1temp)))$weights



#STEP4: CONSISE SUMMARY
#Akaike Weight
   #MODEL 1
   OUTPUT[,"Akaike_Weight_ALL"]

   #MODEL 2
   OUTPUT[,"Akaike_Weight_LAT"]

   #MODEL 3
   OUTPUT[,"Akaike_Weight_GROUP"]

   #MODEL 4
   OUTPUT[,"Akaike_Weight_LAT_GROUP"]

   #MODEL 5
   OUTPUT[,"Akaike_Weight_TERRITORY"]

   #MODEL 6
   OUTPUT[,"Akaike_Weight_2trop_1temp"]


#beta
   #MODEL 1
   OUTPUT[,"ALL_beta"]

   #MODEL 2
   OUTPUT[,"LAT_beta_HIGH"]
   OUTPUT[,"LAT_beta_AMAZON"]

   #MODEL 3

   OUTPUT[,"GROUP_beta_LEARN"]
   OUTPUT[,"GROUP_beta_INNATE"]

   #MODEL 4
   OUTPUT[,"LAT_GROUP_beta_TEMP_LEARN"]
   OUTPUT[,"LAT_GROUP_beta_TEMP_INNATE"]
   OUTPUT[,"LAT_GROUP_beta_AMAZON_LEARN"]
   OUTPUT[,"LAT_GROUP_beta_AMAZON_INNATE"]

   #MODEL 5
   OUTPUT[,"TERRITORY_beta_NOTALLYEAR"]
   OUTPUT[,"TERRITORY_beta_ALLYEAR"]

   #MODEL 6
   OUTPUT[,"beta_2trop_1temp_TEMP"]
   OUTPUT[,"beta_NOTALLYEAR_2trop_1temp_AMAZ"]
   OUTPUT[,"beta_ALLYEAR_2trop_1temp_AMAZ"]


#z
   #MODEL 1
   mean(OUTPUT[,"ALL_z"])

   #MODEL 2
   mean(OUTPUT[,"LAT_z"])

   #MODEL 3
   mean(OUTPUT[,"GROUP_z"])

   #MODEL 4
   mean(OUTPUT[,"LAT_GROUP_z"])

   #MODEL 5
   mean(OUTPUT[,"TERRITORY_z"])

   #MODEL 6
   mean(OUTPUT[,"2trop_1temp_z"])



















