setwd("/XXX/")

#####
library("data.table")
library("tidyverse")
library("ggplot2")
library("survival")
library("coxme")
library("multcomp")
library("car")
library("survminer")
library("ggplotify")
library("emmeans") #libraries
#####

#import data
data_socind <- fread("FigS5_data.csv")

#fitsurv_infection_model <- survfit(Surv(Time, Event) ~ Infection+Injury, data = data_infection_model)

data_socind$Infection <- factor(data_socind$Infection)
data_socind$Incubation <- factor(data_socind$Incubation)
data_socind$Incubator <- factor(data_socind$Incubator)
data_socind$Replicate <- factor(data_socind$Replicate)

#Split dataset by pathogen
data_socind %>%
  filter(Infection %in% "Ecc15") -> data_socind_ecc 
data_socind %>%
  filter(Infection %in% "Bb") -> data_socind_bb


#Ecc15 analysis

null_model_ecc <- coxme ( Surv (time = Time, event = Event) ~ 1                 + ( 1 | Experimentalist) + ( 1 | Replicate) + ( 1 | ColonyNumber) + ( 1 | Incubator) , data = data_socind_ecc)
full_model_ecc <- coxme ( Surv (time = Time, event = Event) ~ 1 + Incubation + ( 1 | Experimentalist) + ( 1 | Replicate) + ( 1 | ColonyNumber) + ( 1 | Incubator), data = data_socind_ecc)

# anova(null_model   ,  full_model )
summary(full_model_ecc)

Anova(full_model_ecc, type="II")

#Bb analysis

null_model_bb <- coxme ( Surv (time = Time, event = Event) ~ 1                 + ( 1 | Experimentalist) + ( 1 | Replicate) + ( 1 | ColonyNumber) + ( 1 | Incubator) , data = data_socind_bb)
full_model_bb <- coxme ( Surv (time = Time, event = Event) ~ 1 + Incubation + ( 1 | Experimentalist) + ( 1 | Replicate) + ( 1 | ColonyNumber) + ( 1 | Incubator), data = data_socind_bb)

# anova(null_model   ,  full_model )
summary(full_model_bb)

Anova(full_model_bb, type="II")

#Plotting

surv_fitEcc <- survfit(Surv(Time, Event) ~ Group_letter, data = data_socind_ecc)

#green : #63AF87
#blue : #37799E
#yellow : #EFB565
#red : #CC5B38

ggsurvplot(surv_fitEcc,
           data = data_socind_ecc, conf.int = TRUE, legend = "top", 
           palette = c("#63AF87", "#37799E"), 
           linetype = c("dotted", "solid"), 
           legend.title="", 
           legend.labs=c("Social", "Individual"))+ 
  guides(colour = guide_legend(nrow = 1))  -> plotSocindecc

print(plotSocindecc)

plotSocindecc <- plotSocindecc$plot

ggsave("surv_socind_Ecc15.jpg", plotSocindecc, width = 8, height = 8, units = "cm")


surv_fitBb <- survfit(Surv(Time, Event) ~ Group_letter, data = data_socind_bb)

#green : #63AF87
#blue : #37799E
#yellow : #EFB565
#red : #CC5B38


ggsurvplot(surv_fitBb,
           data = data_socind_bb, conf.int = TRUE, legend = "top", 
           palette = c("#EFB565", "#CC5B38"), 
           linetype = c("dotted", "solid"), 
           legend.title="", 
           legend.labs=c("Social", "Individual"))+ 
  guides(colour = guide_legend(nrow = 1))  -> plotSocindBb

print(plotSocindBb)

plotSocindBb <- plotSocindBb$plot

plotSocindBb + 
  ggplot2::annotate("text", x = 140, y = 0.8, 
                    label = "*", size = 8) -> plotSocindBb

ggsave("surv_socind_Beauveria.jpg", plotSocindBb, width = 8, height = 8, units = "cm")