#Colin Donihue
#Analysis of morphology and performance data of Podarcis erhardii collected on the Greek Island, Naxos, summer of 2014. 
#Primary comparison - wall / no wall effects on morphology and performance.


#Be sure to load:
library(nlme)
library(car)

NaxosWalls <- read.csv(“DonihueNaxosWalls.csv")

#Perform a normality test on the distribution of jumps - a count variable.
shapiro.test(NaxosWalls$Jumps)
# Result:  Shapiro-Wilk normality test
# W = 0.9435, p-value = 1.18e-09

#subset the data to only be males or females:
maledata<- subset(NaxosWalls, Sex=="M")
femaledata<- subset(NaxosWalls, Sex=="F")

#standardize the covariates
maledata$SVLstanM<-maledata$SVL - mean(maledata$SVL)
femaledata$SVLstanF<-femaledata$SVL - mean(femaledata$SVL)

###################
## Male Morphometrics
###################

model <- lme(SVL~Wall, random=~1|Site, data=maledata) #SVL
summary(model); Anova(model, type="2")
model <- lme(Humerus~Wall, random=~1|Site, data=maledata) #Humerus
summary(model); Anova(model, type="2")
model <- lme(Antebrachium~Wall, random=~1|Site, data=maledata) #Antebrachium
summary(model); Anova(model, type="2")
model <- lme(Wrist_To_Tip~Wall, random=~1|Site, data=maledata) #Wrist_To_Tip
summary(model); Anova(model, type="2")
model <- lme(Longest_Finger~Wall, random=~1|Site, data=maledata) #Longest_Finger
summary(model); Anova(model, type="2")
model <- lme(Femur~Wall, random=~1|Site, data=maledata) #Femur
summary(model); Anova(model, type="2")
model <- lme(Tibia~Wall, random=~1|Site, data=maledata) #Tibia
summary(model); Anova(model, type="2")
model <- lme(Ankle_To_Tip~Wall, random=~1|Site, data=maledata) #Ankle_To_Tip
summary(model); Anova(model, type="2")
model <- lme(Longest_Toe~Wall, random=~1|Site, data=maledata) #Longest_Toe
summary(model); Anova(model, type="2")
model <- lme(Total_Leg_Length~Wall, random=~1|Site, data=maledata) #Total_Leg_Length
summary(model); Anova(model, type="2")

###################
## Female Morphometrics
###################

model <- lme(SVL~Wall, random=~1|Site, data=femaledata) #SVL
summary(model); Anova(model, type="2")
model <- lme(Humerus~Wall, random=~1|Site, data=femaledata) #Humerus
summary(model); Anova(model, type="2")
model <- lme(Antebrachium~Wall, random=~1|Site, data=femaledata) #Antebrachium
summary(model); Anova(model, type="2")
model <- lme(Wrist_To_Tip~Wall, random=~1|Site, data=femaledata) #Wrist_To_Tip
summary(model); Anova(model, type="2")
model <- lme(Longest_Finger~Wall, random=~1|Site, data=femaledata) #Longest_Finger
summary(model); Anova(model, type="2")
model <- lme(Femur~Wall, random=~1|Site, data=femaledata) #Femur
summary(model); Anova(model, type="2")
model <- lme(Tibia~Wall, random=~1|Site, data=femaledata) #Tibia
summary(model); Anova(model, type="2")
model <- lme(Ankle_To_Tip~Wall, random=~1|Site, data=femaledata) #Ankle_To_Tip
summary(model); Anova(model, type="2")
model <- lme(Longest_Toe~Wall, random=~1|Site, data=femaledata) #Longest_Toe
summary(model); Anova(model, type="2")
model <- lme(Total_Leg_Length~Wall, random=~1|Site, data=femaledata) #Total_Leg_Length
summary(model); Anova(model, type="2")

###################
## Relative Male Morphometrics
###################

model <- lme(Humerus~Wall+SVLstanM, random=~1|Site, data=maledata) #Humerus
summary(model); Anova(model, type="2")
model <- lme(Antebrachium~Wall+SVLstanM, random=~1|Site, data=maledata) #Antebrachium
summary(model); Anova(model, type="2")
model <- lme(Wrist_To_Tip~Wall+SVLstanM, random=~1|Site, data=maledata) #Wrist_To_Tip
summary(model); Anova(model, type="2")
model <- lme(Longest_Finger~Wall+SVLstanM, random=~1|Site, data=maledata) #Longest_Finger
summary(model); Anova(model, type="2")
model <- lme(Total_Arm_Length~Wall+SVLstanM, random=~1|Site, data=maledata) #Total_Arm_Length
summary(model); Anova(model, type="2")
model <- lme(Femur~Wall+SVLstanM, random=~1|Site, data=maledata) #Femur
summary(model); Anova(model, type="2")
model <- lme(Tibia~Wall+SVLstanM, random=~1|Site, data=maledata) #Tibia
summary(model); Anova(model, type="2")
model <- lme(Ankle_To_Tip~Wall+SVLstanM, random=~1|Site, data=maledata) #Ankle_To_Tip
summary(model); Anova(model, type="2")
model <- lme(Longest_Toe~Wall+SVLstanM, random=~1|Site, data=maledata) #Longest_Toe
summary(model); Anova(model, type="2")
model <- lme(Total_Leg_Length~Wall+SVLstanM, random=~1|Site, data=maledata) #Total_Leg_Length
summary(model); Anova(model, type="2")

###################
## Relative Female Morphometrics
###################

model <- lme(Humerus~Wall+SVLstanF, random=~1|Site, data=femaledata) #Humerus
summary(model); Anova(model, type="2")
model <- lme(Antebrachium~Wall+SVLstanF, random=~1|Site, data=femaledata) #Antebrachium
summary(model); Anova(model, type="2")
model <- lme(Wrist_To_Tip~Wall+SVLstanF, random=~1|Site, data=femaledata) #Wrist_To_Tip
summary(model); Anova(model, type="2")
model <- lme(Longest_Finger~Wall+SVLstanF, random=~1|Site, data=femaledata) #Longest_Finger
summary(model); Anova(model, type="2")
model <- lme(Total_Arm_Length~Wall+SVLstanF, random=~1|Site, data=femaledata) #Total_Arm_Length
summary(model); Anova(model, type="2")
model <- lme(Femur~Wall+SVLstanF, random=~1|Site, data=femaledata) #Femur
summary(model); Anova(model, type="2")
model <- lme(Tibia~Wall+SVLstanF, random=~1|Site, data=femaledata) #Tibia
summary(model); Anova(model, type="2")
model <- lme(Ankle_To_Tip~Wall+SVLstanF, random=~1|Site, data=femaledata) #Ankle_To_Tip
summary(model); Anova(model, type="2")
model <- lme(Longest_Toe~Wall+SVLstanF, random=~1|Site, data=femaledata) #Longest_Toe
summary(model); Anova(model, type="2")
model <- lme(Total_Leg_Length~Wall+SVLstanF, random=~1|Site, data=femaledata) #Total_leg_Length
summary(model); Anova(model, type="2")

###################
## Performance Measures
###################

#Some sprint runs failed, resulting in NAs in the data. This code will subset the data to remove rows with NAs. 
MaleVelRockDF<-maledata[!is.na(maledata$Max_Vel_Rock),]
MaleAccelRockDF<-maledata[!is.na(maledata$Max_Accel_Rock),]
MaleVelSandDF<-maledata[!is.na(maledata$Max_Vel_Sand),]
MaleAccelSandDF<-maledata[!is.na(maledata$Max_Accel_Sand),]
MaleJumpsDF<-maledata[!is.na(maledata$Jumps),]

FemaleVelRockDF<-femaledata[!is.na(femaledata$Max_Vel_Rock),]
FemaleAccelRockDF<-femaledata[!is.na(femaledata$Max_Accel_Rock),]
FemaleVelSandDF<-femaledata[!is.na(femaledata$Max_Vel_Sand),]
FemaleAccelSandDF<-femaledata[!is.na(femaledata$Max_Accel_Sand),]
FemaleJumpsDF<-femaledata[!is.na(femaledata$Jumps),]

#Standardizing sprint temperature covariate within each dataframe 
MaleVelRockDF$TempstanM<-MaleVelRockDF$Sprint_Temp - mean(MaleVelRockDF$Sprint_Temp)
MaleAccelRockDF$TempstanM<-MaleAccelRockDF$Sprint_Temp - mean(MaleAccelRockDF$Sprint_Temp)
MaleVelSandDF$TempstanM<-MaleVelSandDF$Sprint_Temp - mean(MaleVelSandDF$Sprint_Temp)
MaleAccelSandDF$TempstanM<-MaleAccelSandDF$Sprint_Temp - mean(MaleAccelSandDF$Sprint_Temp)
MaleJumpsDF$TempstanM<-MaleJumpsDF$Sprint_Temp - mean(MaleJumpsDF$Sprint_Temp)

FemaleVelRockDF$TempstanF<-FemaleVelRockDF$Sprint_Temp - mean(FemaleVelRockDF$Sprint_Temp)
FemaleAccelRockDF$TempstanF<-FemaleAccelRockDF$Sprint_Temp - mean(FemaleAccelRockDF$Sprint_Temp)
FemaleVelSandDF$TempstanF<-FemaleVelSandDF$Sprint_Temp - mean(FemaleVelSandDF$Sprint_Temp)
FemaleAccelSandDF$TempstanF<-FemaleAccelSandDF$Sprint_Temp - mean(FemaleAccelSandDF$Sprint_Temp)
FemaleJumpsDF$TempstanF<-FemaleJumpsDF$Sprint_Temp - mean(FemaleJumpsDF$Sprint_Temp)

###################
## Male Performance
###################

model<-lme(Max_Vel_Rock~Wall, random=list(~1|Site, ~1|TempstanM), data=MaleVelRockDF) #Max_Vel_Rock
summary(model); Anova(model, type="2")
model<-lme(Max_Accel_Rock~Wall, random=list(~1|Site, ~1|TempstanM), data=MaleAccelRockDF) #Max_Accel_Rock
summary(model); Anova(model, type="2")
model<-lme(Max_Vel_Sand~Wall, random=list(~1|Site, ~1|TempstanM), data=MaleVelSandDF) #Max_Vel_Sand
summary(model); Anova(model, type="2")
model<-lme(Max_Accel_Sand~Wall, random=list(~1|Site, ~1|TempstanM), data=MaleAccelSandDF) #Max_Accel_Sand
summary(model); Anova(model, type="2")
model<-lme(Jumps~Wall, random=list(~1|Site, ~1|TempstanM), data=MaleJumpsDF) #Jumps
summary(model); Anova(model, type="2")

###################
## Female Performance
###################

model<-lme(Max_Vel_Rock~Wall, random=list(~1|Site, ~1|TempstanF), data=FemaleVelRockDF) #Max_Vel_Rock
summary(model); Anova(model, type="2")
model<-lme(Max_Accel_Rock~Wall, random=list(~1|Site, ~1|TempstanF), data=FemaleAccelRockDF) #Max_Accel_Rock
summary(model); Anova(model, type="2")
model<-lme(Max_Vel_Sand~Wall, random=list(~1|Site, ~1|TempstanF), data=FemaleVelSandDF) #Max_Vel_Sand
summary(model); Anova(model, type="2")
model<-lme(Max_Accel_Sand~Wall, random=list(~1|Site, ~1|TempstanF), data=FemaleAccelSandDF) #Max_Accel_Sand
summary(model); Anova(model, type="2")
model<-lme(Jumps~Wall, random=list(~1|Site, ~1|TempstanF), data=FemaleJumpsDF) #Jumps
summary(model); Anova(model, type="2")

###################
## Relative Male Performance
###################

model<-lme(Max_Vel_Rock~Wall+SVLstanM, random=list(~1|Site, ~1|TempstanM), data=MaleVelRockDF) #Max_Vel_Rock
summary(model); Anova(model, type="2")
model<-lme(Max_Accel_Rock~Wall+SVLstanM, random=list(~1|Site, ~1|TempstanM), data=MaleAccelRockDF) #Max_Accel_Rock
summary(model); Anova(model, type="2")
model<-lme(Max_Vel_Sand~Wall+SVLstanM, random=list(~1|Site, ~1|TempstanM), data=MaleVelSandDF) #Max_Vel_Sand
summary(model); Anova(model, type="2")
model<-lme(Max_Accel_Sand~Wall+SVLstanM, random=list(~1|Site, ~1|TempstanM), data=MaleAccelSandDF) #Max_Accel_Sand
summary(model); Anova(model, type="2")
model<-lme(Jumps~Wall+SVLstanM, random=list(~1|Site, ~1|TempstanM), data=MaleJumpsDF) #Jumps
summary(model); Anova(model, type="2")

###################
## Relative Female Performance
###################

model<-lme(Max_Vel_Rock~Wall+SVLstanF, random=list(~1|Site, ~1|TempstanF), data=FemaleVelRockDF) #Max_Vel_Rock
summary(model); Anova(model, type="2")
model<-lme(Max_Accel_Rock~Wall+SVLstanF, random=list(~1|Site, ~1|TempstanF), data=FemaleAccelRockDF) #Max_Accel_Rock
summary(model); Anova(model, type="2")
model<-lme(Max_Vel_Sand~Wall+SVLstanF, random=list(~1|Site, ~1|TempstanF), data=FemaleVelSandDF) #Max_Vel_Sand
summary(model); Anova(model, type="2")
model<-lme(Max_Accel_Sand~Wall+SVLstanF, random=list(~1|Site, ~1|TempstanF), data=FemaleAccelSandDF) #Max_Accel_Sand
summary(model); Anova(model, type="2")
model<-lme(Jumps~Wall+SVLstanF, random=list(~1|Site, ~1|TempstanF), data=FemaleJumpsDF) #Jumps
summary(model); Anova(model, type="2")