library(lmerTest)

library(car)

# set working directory

setwd()

## Analyses of parasite and host performance ####

RhiPlaData <- read.table("RhiPlaData.csv", header=TRUE, sep=";")

RhiPlaData$HostPop<-as.factor(RhiPlaData$HostPop)

RhiPlaData$ParPop<-as.factor(RhiPlaData$ParPop)

RhiPlaData$Host_Type<-as.factor(RhiPlaData$Host_Type)

RhiPlaData$naive_host<-as.factor(RhiPlaData$naive_host)

RhiPlaData$Par_yn<-as.factor(RhiPlaData$Par_yn)

# Early parasite size (and other parasite traits):

model1<-lmer(Par_leaflength~Host_Type + (1|ParPop) , data=RhiPlaData)

anova(model1)

# Early host size (and other host traits):

model2<-lmer(log_Hostsize1~naive_host*Par_yn + (1|HostPop) , data=RhiPlaData)

anova(model2)

# (no) relationship between parasite size and host size:

summary(lm(Par_leaflength~log_Hostsize1, data=RhiPlaData))

# Parasite survival since establishment:

model3<-glmer(Mort_sinceM1~Host_Type + (1|ParPop), data=RhiPlaData, family = binomial)

Anova(model3)

# Parasite biomass:

RhiPlaData$log_Parmass<-log10(RhiPlaData$Par_biomass)

model4<-lmer(log_Parmass~Host_Type + (1|ParPop) , data=RhiPlaData)

anova(model4)

# Possible effect of truncation?

# Select largest parasite from each combination of ParPop x HosPop

RhiPlaData$test <- with(RhiPlaData, ave(Par_biomass, list(HostPop, ParPop), FUN=function (x) rank(max(x, na.rm=TRUE) - x)))

RhiPlaTrunc<-subset(RhiPlaData, RhiPlaData$test == 1)

model5<-lmer(log_Parmass~Host_Type + (1|ParPop) , data=RhiPlaTrunc)

anova(model5)

# Subset: only parasite populations surviving with all host types:

RhiPlaSmall<-subset(RhiPlaData, ParPop == "2" | ParPop == "3" | ParPop == "4" | ParPop == "5")

model6<-lmer(log_Parmass~Host_Type + (1|ParPop) , data=RhiPlaSmall)

anova(model6)

# multiplicative fitness

Multifit <- read.table("Multifit.csv", header=TRUE, sep=";")

Multifit$Host_Type<-as.factor(Multifit$Host_Type)

Multifit$ParPop<-as.factor(Multifit$ParPop)

model7<-lmer(Multifit_mg~Host_Type + (1|ParPop), data=Multifit)

anova(model7)

model_dist<-lm(Multifit_mg~logPopdist+ I(logPopdist**2), data=Multifit)

summary(model_dist)

# host biomass

model8<-lmer(Host_biomass~naive_host*Par_yn + (1|HostPop) , data=RhiPlaData)

anova(model8)

# singularity => remove host pop, but keep Type III SS. For correct contrasts:

type3<-list(naive_host = contr.sum, Par_yn = contr.sum)

model8b<-lm(Host_biomass~naive_host*Par_yn , data=RhiPlaData, contrasts= type3)

Anova(model8b, type=3)

anova(model8, model8b)

# models not significantly different.

# Analysis restricted to experienced hosts

exp_hosts<-subset(RhiPlaData, Host_Type != "-1")

exp_hosts$sympatric<-ifelse(exp_hosts$Host_Type=="2", "1", "0")

exp_hosts$near_allo<-ifelse(exp_hosts$Host_Type=="1", "1", "0")

model9<-lm(Host_biomass~Par_yn + sympatric + near_allo , data=exp_hosts)

anova(model9)

# (weak) relationship between parasite size and host size:

summary(lm(log10(Host_biomass)~log_Parmass, data=RhiPlaData))

## Root Allocation Experiment ####

RAE <- read.table("RootAllocationExp.csv", header=TRUE, sep=";")

RAE$HostPop<-as.factor(RAE$HostPop)

RAE$ParPop<-as.factor(RAE$ParPop)

RAE$Par_yn<-as.factor(RAE$Par_yn)

# Parasite leaf length: simple ANOVA

m1<-lm(P_leaflength~attached, data=RAE)

anova(m1)

# including population as random factor

m2<-lmer(P_leaflength~attached + (1|HostPop) + (1|ParPop), data=RAE)

# models including Parasite Population have singularity problems

m2<-lmer(P_leaflength~attached + (1|HostPop) , data=RAE)

anova(m2)

# effect still significant (p < 0.001)

anova(m2, m1)

# model with random factor not significantly better

m3<-lmer(P_leaflength~attached*naive_host + (1|HostPop) , data=RAE)

anova(m3)

# only attachment significant (p = 0.003)

# Parasite shoot mass: simple ANOVA

m4<-lm(P_shootmass~attached, data=RAE)

anova(m4)

# Host traits

m5<-lm(H_shootmass~Par_yn + attached, data=RAE)

Anova(m5)

m6<-lmer(H_shootmass~Par_yn + attached + (1|ParPop) , data=RAE)

anova(m6)

anova(m6, m5)

# not significantly different

m7<-lmer(H_shootmass~attached*naive_host + (1|ParPop) , data=RAE)

anova(m7)

# interesting interaction: attached : naive host

m8<-lm(H_shootmass~attached*naive_host , data=RAE)

anova(m8)

anova(m7, m8)

# m8 has smaller AIC but not significantly better

# For comparison: This effect is not significant for parasite presence:

m9<-lm(H_shootmass~Par_yn*naive_host , data=RAE)

anova(m9)

# Host root mass

m10<-lm(H_rootmass~Par_yn + attached, data=RAE)

anova(m10)

# General linear models:

type3<-list(Par_yn = contr.sum)

m11<-lm(log10(H_rootmass)~log10(H_shootmass)*Par_yn, data=RAE, contrasts= type3)

Anova(m11, type=3)

m12<-lm(HLR_Length~H_rootmass*Par_yn, data=RAE, contrasts= type3)

Anova(m12, type=3)

# Effects of attachment on root traits

m13<-lm(HLR_Thickness~Par_yn, data=RAE)

anova(m13)

m14<-lm(HPR_Thickness~attached, data=RAE)

anova(m14)

# attachment significant, parasite presence not.

m15<-lm(HPR_Thickness~Par_yn + attached, data=RAE)

anova(m15)

m16<-lm(HLR_Thickness~Par_yn + attached, data=RAE)

anova(m16)

# effects on attachment?

model_attachment<-glmer(attached~naive_host + (1|ParPop), data=subset(RAE, RAE$Par_yn==1), family = binomial)

Anova(model_attachment)

model_attachment<-glm(attached~naive_host, data=subset(RAE, RAE$Par_yn==1), family = binomial)

Anova(model_attachment)

# both not significant. Attachment did not differ between naive and experienced hosts

## Host Choice Experiment ####

Agar <- read.table("HostChoiceExp.csv", header=TRUE, sep=";")

Agar$ID<-as.factor(Agar$ID)

model_choice<-lmer(sqrtLength~sterile*symp+ (1|ID)+(1|ParPop), data=Agar)

anova(model_choice, type=1)