# Day = day of experiment
# Time = time of events in hours
# Event = binary 1/dead 0/censored
# Infection and Fight = experimental conditions
# Colony series = batch number of the colony
# Colony number = unique colony identifier
# Experimentalist = experimentalist (slightly different manipulation practice)
# Group_letter are dummy variables used in the script to identify experimental conditions. 


library("data.table")
library("tidyverse")
library("ggplot2")
library("survival")
library("coxme")
library("survminer")
library("ggplotify")
library("plyr")

#import data
data_fight <- fread("Fig4D_data.csv")

data_fight$Infection <- factor(data_fight$Infection)
data_fight$Fight <- factor(data_fight$Fight, levels = c("No", "Yes"))


#Plotting

surv_fitFights <- survfit(Surv(Time, Event) ~ Group_letter, data = data_fight)

ggsurvplot(surv_fitFights,
           data = data_fight, conf.int = TRUE, legend = "top", 
           palette = c("#63AF87", "#37799E", "#EFB565", "#CC5B38"), 
           linetype = c("dotted", "solid", "dotted", "solid"), 
           legend.title="", 
           legend.labs=c("No fight - Uninfected", "No fight - Infected", "Fight - Uninfected", "Fight - Infected"))+ 
  guides(colour = guide_legend(nrow = 2)) -> plotFights

print(plotFights)


#STATS

null_model <- coxme ( Surv (time = Time, event = Event) ~ 1                + ( 1 | Replicate)  + ( 1 | Experimentalist) + (1 | ColonySeries)+ (1 | ColonyNumber) , data = data_fight)
full_model <- coxme ( Surv (time = Time, event = Event) ~ 1 + Infection + Fight + ( 1 | Replicate)  + ( 1 | Experimentalist)  , data = data_fight)
interaction_model_n <- coxme ( Surv (time = Time, event = Event) ~ 1 + Infection * Fight + ( 1 | Replicate)  + ( 1 | Experimentalist) 
                               , data = data_fight)

# AIC-based method (simply compare the AIC of each model to determine the best fit)
# necessary because there is no death in the control condition 
# and Wald tests are therefore no usable (coefficient converges to infinity)