#Import data
setwd("C:/Users/awvater/Desktop")
Morts <- read.csv("C:/Users/awvater/Desktop/Figure_2data")
View(Morts)
Morts = subset(Morts,Morts$Day>0) # eliminate abalone not seen after start

#Install
install.packages('survival')

#Pull up libraries
library(survival)

#create vector "Survival"
Morts$Survival = Surv(Morts$Day,Morts$Event)

#Create survfit model
Surv1 = survfit(Morts$Survival ~ Morts$Species + Morts$Group + Morts$Temp)
summary(Surv1)

#Plot
plot(Surv1,col=1:6)
legend("bottomleft", c("Ca. Xc - exposed, Ambient, Red",
                       "Ca. Xc - exposed, Heated, Red", 
                       "Naive, Heated, Red", 
                       "Ca. Xc - exposed, Ambient, White",
                       "Ca. Xc - exposed, Heated, White",
                       "Naive, Heated, White"), fill = (col = 1:6))


#cox models
#exposed, heated survdiff between White vs Red
m_temp = subset(Morts,Morts$Group=='Exposed')
m_temp = subset(Morts,Morts$Temp == 'Heated')

mtemp = survdiff(Surv(Day, Event) ~ Species, data=m_temp)
summary(mtemp)
mtemp
cox6 = coxph(Surv(Day, Event) ~ Species , data=m_temp)
summary(cox6)

#interactions
cox1 = coxph(Surv(Day, Event) ~ Species + Temp + Temp*Species  , data=Morts[Morts$Group=='Exposed',])
summary(cox1)

cox2 = coxph(Surv(Day, Event) ~ Species + Group + Species:Group  , data=Morts[Morts$Temp=='Heated',])
summary(cox2)



#cox3 = coxph(Surv(Day, Event) ~ Species + Temp  , data=Morts[Morts$Group=='Exposed',])
#summary(cox3)

#cox4 = coxph(Surv(Day, Event) ~ Species + Temp + Group  , data=Morts)
#summary(cox4)

#cox5 = coxph(Surv(Day, Event) ~ Species + Temp + Temp:Species  , data=Morts)
#summary(cox5)

#cox6 = coxph(Surv(Day, Event) ~ Species , data=Morts[Morts$Temp=='Ambient',])
#summary(cox6)


### Why are these here? They seem to do nothing... 
# set "Ambient" as referent group, look at exposure effects
Morts$Temp = relevel(Morts$Temp,ref = 'Ambient')
# set unexposed as referent group, look at temp effects
Morts$Group = relevel(Morts$Group,ref = 'Unexposed')


#get rows with subset
mexposure = survdiff(Surv(Day, Event) ~ Group + Temp+Species, data=Morts)
mspecies = survdiff(Surv(Day, Event) ~ Group + Temp+Species, data=Morts)

#other plots 
install.packages('devtools')
install_github('sinhrks/ggfortify')

install.packages('ggfortify')

?autoplot
autoplot(Surv1)
ggsurv(Surv1)