# Summer 2020 Arthropods Research
library(vegan)
morphospecies=read.csv(file="Morphospecies.csv")
predictors=read.csv(file="Arthropods.csv")
summary(predictors)
nmds1<-metaMDS(morphospecies,distance="bray",trace=FALSE)
nmds1
names1<-as.character(predictors$Type)
names1
pl<-plot(nmds1,type="none")

#Correlation Coefficients
correlations.responses=read.csv(file="Correlation.Predictors.csv")
mydata.rcor=cor(correlations.responses)
mydata.rcor
write.csv(mydata.rcor,'data1.csv')

#Display names of orders and sites
points(nmds1,display="sp",col="dark gray",cex=0.8)


text(nmds1,"sp",col="black",cex=1)
text(nmds1,"sites",labels=names1,cex=1) #shows words on r-plot

#Envfit fits an environmental variables onto an ordiation
fit<-envfit(nmds1~Total.Flowers+Purple.Flowers+White.Flowers+Humidity+Temperature+Impervious+PAV+Type,predictors,perm=4999,col="dark gray")
plot(fit, col="black")

#ADONIS TEST - PERMANOVA
stat.test<-adonis(morphospecies~Total.Flowers+Purple.Flowers+White.Flowers+Humidity+Temperature+Impervious+Type+Perimeter.Area,predictors,permutations=2000,method="bray")
stat.test

# Fit of significant environmental variables
fit<-envfit(nmds1~Impervious+White.Flowers+Humidity+Type,predictors,perm=4999)
pl<-plot(nmds1,type="none")
plot(fit, col="black")

#Calculating Diversity & Evenness
diversity(orders,index="shannon")
H<-diversity(orders) #Shannon Diversity Index
J<-H/log(specnumber(orders)) # Pielou's Evenness
predictors$Diversity<-c(H) #adds diversity to the predictors dataframe
predictors$Evenness<-c(J) #adds evenness to the predictors dataframe
H

# GLM Analysis - Abundance Candidate Models
library(AER)
library(MuMIn)
library(MASS)

## List of Candidate Models for Total Arthropod Abundance
abundance.model1=glm.nb(Abundance~Total.Flowers,data=predictors)
summary(abundance.model1)
AICc(abundance.model1)
with(summary(abundance.model1),1-deviance/null.deviance)

abundance.model2=glm.nb(Abundance~Purple.Flowers,data=predictors)
summary(abundance.model2)
AICc(abundance.model2)
with(summary(abundance.model2),1-deviance/null.deviance)

abundance.model3=glm.nb(Abundance~White.Flowers,data=predictors)
summary(abundance.model3)
AICc(abundance.model3)
with(summary(abundance.model3),1-deviance/null.deviance)

abundance.model4=glm.nb(Abundance~Humidity,data=predictors)
summary(abundance.model4)
AICc(abundance.model4)
with(summary(abundance.model4),1-deviance/null.deviance)

abundance.model5=glm.nb(Abundance~Temperature,data=predictors)
summary(abundance.model5)
AICc(abundance.model5)
with(summary(abundance.model5),1-deviance/null.deviance)

abundance.model6=glm.nb(Abundance~Impervious,data=predictors)
summary(abundance.model6)
AICc(abundance.model6)
with(summary(abundance.model6),1-deviance/null.deviance)

abundance.model7=glm.nb(Abundance~Perimeter.Area,data=predictors)
summary(abundance.model7)
AICc(abundance.model7)
with(summary(abundance.model7),1-deviance/null.deviance)


## List of Candidate Models for Arthropod Diversity
diversity.model1=glm(Diversity~Total.Flowers,data=predictors)
summary(diversity.model1)
AICc(diversity.model1)
with(summary(diversity.model1),1-deviance/null.deviance)

diversity.model2=glm(Diversity~Purple.Flowers,data=predictors)
summary(diversity.model2)
AICc(diversity.model2)
with(summary(diversity.model2),1-deviance/null.deviance)

diversity.model3=glm(Diversity~White.Flowers,data=predictors)
summary(diversity.model3)
AICc(diversity.model3)
with(summary(diversity.model3),1-deviance/null.deviance)

diversity.model4=glm(Diversity~Humidity,data=predictors)
summary(diversity.model4)
AICc(diversity.model4)
with(summary(diversity.model4),1-deviance/null.deviance)

diversity.model5=glm(Diversity~Temperature,data=predictors)
summary(diversity.model5)
AICc(diversity.model5)
with(summary(diversity.model5),1-deviance/null.deviance)

diversity.model6=glm(Diversity~Impervious,data=predictors)
summary(diversity.model6)
AICc(diversity.model6)
with(summary(diversity.model6),1-deviance/null.deviance)

diversity.model7=glm(Diversity~Perimeter.Area,data=predictors)
summary(diversity.model7)
AICc(diversity.model7)
with(summary(diversity.model7),1-deviance/null.deviance)

## List of Candidate Models for Arthropod Morphospecies Richness
richness.model1=glm(Richness~Total.Flowers,data=predictors,family="poisson")
summary(richness.model1)
AICc(richness.model1)
with(summary(richness.model1),1-deviance/null.deviance)

richness.model2=glm(Richness~Purple.Flowers,family="poisson",data=predictors)
summary(richness.model2)
AICc(richness.model2)
with(summary(richness.model2),1-deviance/null.deviance)

richness.model3=glm(Richness~White.Flowers,data=predictors,family="poisson")
summary(richness.model3)
AICc(richness.model3)
with(summary(richness.model3),1-deviance/null.deviance)

richness.model4=glm(Richness~Humidity,data=predictors,family="poisson")
summary(richness.model4)
AICc(richness.model4)
with(summary(richness.model4),1-deviance/null.deviance)

richness.model5=glm(Richness~Temperature,data=predictors,family="poisson")
summary(richness.model5)
AICc(richness.model5)
with(summary(richness.model5),1-deviance/null.deviance)

richness.model6=glm(Richness~Impervious,data=predictors,family="poisson")
summary(richness.model6)
AICc(richness.model6)
with(summary(richness.model6),1-deviance/null.deviance)

richness.model7=glm(Richness~Perimeter.Area,data=predictors,family="poisson")
summary(richness.model7)
AICc(richness.model7)
with(summary(richness.model7),1-deviance/null.deviance)

## List of Candidate Models for Hymenoptera Abundance
hymenoptera.abundance.model1=glm.nb(Hymenoptera.Abundance~Total.Flowers,data=predictors)
summary(hymenoptera.abundance.model1)
AICc(hymenoptera.abundance.model1)
with(summary(hymenoptera.abundance.model1),1-deviance/null.deviance)

hymenoptera.abundance.model2=glm.nb(Hymenoptera.Abundance~Purple.Flowers,data=predictors)
summary(hymenoptera.abundance.model2)
AICc(hymenoptera.abundance.model2)
with(summary(hymenoptera.abundance.model2),1-deviance/null.deviance)

hymenoptera.abundance.model3=glm.nb(Hymenoptera.Abundance~White.Flowers,data=predictors)
summary(hymenoptera.abundance.model3)
AICc(hymenoptera.abundance.model3)
with(summary(hymenoptera.abundance.model3),1-deviance/null.deviance)

hymenoptera.abundance.model4=glm.nb(Hymenoptera.Abundance~Humidity,data=predictors)
summary(hymenoptera.abundance.model4)
AICc(hymenoptera.abundance.model4)
with(summary(hymenoptera.abundance.model4),1-deviance/null.deviance)

hymenoptera.abundance.model5=glm.nb(Hymenoptera.Abundance~Temperature,data=predictors)
summary(hymenoptera.abundance.model5)
AICc(hymenoptera.abundance.model5)
with(summary(hymenoptera.abundance.model5),1-deviance/null.deviance)

hymenoptera.abundance.model6=glm.nb(Hymenoptera.Abundance~Impervious,data=predictors)
summary(hymenoptera.abundance.model6)
AICc(hymenoptera.abundance.model6)
with(summary(hymenoptera.abundance.model6),1-deviance/null.deviance)

hymenoptera.abundance.model7=glm.nb(Hymenoptera.Abundance~Perimeter.Area,data=predictors)
summary(hymenoptera.abundance.model7)
AICc(hymenoptera.abundance.model7)
with(summary(hymenoptera.abundance.model7),1-deviance/null.deviance)


## List of Candidate Models for Hemiptera Abundance
hemiptera.abundance.model1=glm.nb(Hemiptera.Abundance~Total.Flowers,data=predictors)
summary(hemiptera.abundance.model1)
AICc(hemiptera.abundance.model1)
with(summary(hemiptera.abundance.model1),1-deviance/null.deviance)

hemiptera.abundance.model2=glm.nb(Hemiptera.Abundance~Purple.Flowers,data=predictors)
summary(hemiptera.abundance.model2)
AICc(hemiptera.abundance.model2)
with(summary(hemiptera.abundance.model2),1-deviance/null.deviance)

hemiptera.abundance.model3=glm.nb(Hemiptera.Abundance~White.Flowers,data=predictors)
summary(hemiptera.abundance.model3)
AICc(hemiptera.abundance.model3)
with(summary(hemiptera.abundance.model3),1-deviance/null.deviance)

hemiptera.abundance.model4=glm.nb(Hemiptera.Abundance~Humidity,data=predictors)
summary(hemiptera.abundance.model4)
AICc(hemiptera.abundance.model4)
with(summary(hemiptera.abundance.model4),1-deviance/null.deviance)

hemiptera.abundance.model5=glm.nb(Hemiptera.Abundance~Temperature,data=predictors)
summary(hemiptera.abundance.model5)
AICc(hemiptera.abundance.model5)
with(summary(hemiptera.abundance.model5),1-deviance/null.deviance)

hemiptera.abundance.model6=glm.nb(Hemiptera.Abundance~Impervious,data=predictors)
summary(hemiptera.abundance.model6)
AICc(hemiptera.abundance.model6)
with(summary(hemiptera.abundance.model6),1-deviance/null.deviance)

hemiptera.abundance.model7=glm.nb(Hemiptera.Abundance~Perimeter.Area,data=predictors)
summary(hemiptera.abundance.model7)
AICc(hemiptera.abundance.model7)
with(summary(hemiptera.abundance.model7),1-deviance/null.deviance)


## List of Candidate Models for Diptera Abundance
diptera.abundance.model1=glm.nb(Diptera.Abundance~Total.Flowers,data=predictors)
summary(diptera.abundance.model1)
AICc(diptera.abundance.model1)
with(summary(diptera.abundance.model1),1-deviance/null.deviance)

diptera.abundance.model2=glm.nb(Diptera.Abundance~Purple.Flowers,data=predictors)
summary(diptera.abundance.model2)
AICc(diptera.abundance.model2)
with(summary(diptera.abundance.model2),1-deviance/null.deviance)

diptera.abundance.model3=glm.nb(Diptera.Abundance~White.Flowers,data=predictors)
summary(diptera.abundance.model3)
AICc(diptera.abundance.model3)
with(summary(diptera.abundance.model3),1-deviance/null.deviance)

diptera.abundance.model4=glm.nb(Diptera.Abundance~Humidity,data=predictors)
summary(diptera.abundance.model4)
AICc(diptera.abundance.model4)
with(summary(diptera.abundance.model4),1-deviance/null.deviance)

diptera.abundance.model5=glm.nb(Diptera.Abundance~Temperature,data=predictors)
summary(diptera.abundance.model5)
AICc(diptera.abundance.model5)
with(summary(diptera.abundance.model5),1-deviance/null.deviance)

diptera.abundance.model6=glm.nb(Diptera.Abundance~Impervious,data=predictors)
summary(diptera.abundance.model6)
AICc(diptera.abundance.model6)
with(summary(diptera.abundance.model6),1-deviance/null.deviance)

diptera.abundance.model7=glm.nb(Diptera.Abundance~Perimeter.Area,data=predictors)
summary(diptera.abundance.model7)
AICc(diptera.abundance.model7)
with(summary(diptera.abundance.model7),1-deviance/null.deviance)


## List of Candidate Models for Hymenoptera Richness
hymenoptera.richness.model1=glm(Hymenoptera.Richness~Total.Flowers,family="poisson",data=predictors)
summary(hymenoptera.richness.model1)
AICc(hymenoptera.richness.model1)
with(summary(hymenoptera.richness.model1),1-deviance/null.deviance)

hymenoptera.richness.model2=glm(Hymenoptera.Richness~Purple.Flowers,family="poisson",data=predictors)
summary(hymenoptera.richness.model2)
AICc(hymenoptera.richness.model2)
with(summary(hymenoptera.richness.model2),1-deviance/null.deviance)

hymenoptera.richness.model3=glm(Hymenoptera.Richness~White.Flowers,family="poisson",data=predictors)
summary(hymenoptera.richness.model3)
AICc(hymenoptera.richness.model3)
with(summary(hymenoptera.richness.model3),1-deviance/null.deviance)

hymenoptera.richness.model4=glm(Hymenoptera.Richness~Humidity,family="poisson",data=predictors)
summary(hymenoptera.richness.model4)
AICc(hymenoptera.richness.model4)
with(summary(hymenoptera.richness.model4),1-deviance/null.deviance)

hymenoptera.richness.model5=glm(Hymenoptera.Richness~Temperature,family="poisson",data=predictors)
summary(hymenoptera.richness.model5)
AICc(hymenoptera.richness.model5)
with(summary(hymenoptera.richness.model5),1-deviance/null.deviance)

hymenoptera.richness.model6=glm(Hymenoptera.Richness~Impervious,family="poisson",data=predictors)
summary(hymenoptera.richness.model6)
AICc(hymenoptera.richness.model6)
with(summary(hymenoptera.richness.model6),1-deviance/null.deviance)

hymenoptera.richness.model7=glm(Hymenoptera.Richness~Perimeter.Area,family="poisson",data=predictors)
summary(hymenoptera.richness.model7)
AICc(hymenoptera.richness.model7)
with(summary(hymenoptera.richness.model7),1-deviance/null.deviance)


## List of Candidate Models for Hemiptera Richness
hemiptera.richness.model1=glm(Hemiptera.Richness~Total.Flowers,family="poisson",data=predictors)
summary(hemiptera.richness.model1)
AICc(hemiptera.richness.model1)
with(summary(hemiptera.richness.model1),1-deviance/null.deviance)

hemiptera.richness.model2=glm(Hemiptera.Richness~Purple.Flowers,family="poisson",data=predictors)
summary(hemiptera.richness.model2)
AICc(hemiptera.richness.model2)
with(summary(hemiptera.richness.model2),1-deviance/null.deviance)

hemiptera.richness.model3=glm(Hemiptera.Richness~White.Flowers,family="poisson",data=predictors)
summary(hemiptera.richness.model3)
AICc(hemiptera.richness.model3)
with(summary(hemiptera.richness.model3),1-deviance/null.deviance)

hemiptera.richness.model4=glm(Hemiptera.Richness~Humidity,family="poisson",data=predictors)
summary(hemiptera.richness.model4)
AICc(hemiptera.richness.model4)
with(summary(hemiptera.richness.model4),1-deviance/null.deviance)

hemiptera.richness.model5=glm(Hemiptera.Richness~Temperature,family="poisson",data=predictors)
summary(hemiptera.richness.model5)
AICc(hemiptera.richness.model5)
with(summary(hemiptera.richness.model5),1-deviance/null.deviance)

hemiptera.richness.model6=glm(Hemiptera.Richness~Impervious,family="poisson",data=predictors)
summary(hemiptera.richness.model6)
AICc(hemiptera.richness.model6)
with(summary(hemiptera.richness.model6),1-deviance/null.deviance)

hemiptera.richness.model7=glm(Hemiptera.Richness~Perimeter.Area,family="poisson",data=predictors)
summary(hemiptera.richness.model7)
AICc(hemiptera.richness.model7)
with(summary(hemiptera.richness.model7),1-deviance/null.deviance)


## List of Candidate Models for Diptera Richness
diptera.richness.model1=glm(Diptera.Richness~Total.Flowers,family="poisson",data=predictors)
summary(diptera.richness.model1)
AICc(diptera.richness.model1)
with(summary(diptera.richness.model1),1-deviance/null.deviance)

diptera.richness.model2=glm(Diptera.Richness~Purple.Flowers,family="poisson",data=predictors)
summary(diptera.richness.model2)
AICc(diptera.richness.model2)
with(summary(diptera.richness.model2),1-deviance/null.deviance)

diptera.richness.model3=glm(Diptera.Richness~White.Flowers,family="poisson",data=predictors)
summary(diptera.richness.model3)
AICc(diptera.richness.model3)
with(summary(diptera.richness.model3),1-deviance/null.deviance)

diptera.richness.model4=glm(Diptera.Richness~Humidity,family="poisson",data=predictors)
summary(diptera.richness.model4)
AICc(diptera.richness.model4)
with(summary(diptera.richness.model4),1-deviance/null.deviance)

diptera.richness.model5=glm(Diptera.Richness~Temperature,family="poisson",data=predictors)
summary(diptera.richness.model5)
AICc(diptera.richness.model5)
with(summary(diptera.richness.model5),1-deviance/null.deviance)

diptera.richness.model6=glm(Diptera.Richness~Impervious,family="poisson",data=predictors)
summary(diptera.richness.model6)
AICc(diptera.richness.model6)
with(summary(diptera.richness.model6),1-deviance/null.deviance)

diptera.richness.model7=glm(Diptera.Richness~Perimeter.Area,family="poisson",data=predictors)
summary(diptera.richness.model7)
AICc(diptera.richness.model7)
with(summary(diptera.richness.model7),1-deviance/null.deviance)