################################Spatial Autocorrelation R Code#### ############################################### ################################################################################################################## ###Table must have the form of columns with each pair-wise comparison having a genetic and euclidean distance#### Co1<-read.table(file="WItube.txt", header=T, sep="\t") #Converts to data frameW# distmat.df<-as.data.frame(Co1) #cuts the data by Euclidean distance categories, the 400 should include the maximum Euclidean distance within the analysis# #length.out corresponds the number of discrete categories## distfact<-cut(distmat.df[,3],breaks=seq(0,400,length.out=79),labels=seq(0,400,length.out=79)[-1]) #Sets the number of permutations for the spatial autocorrelation analysis#### n.permutations<-1000 #Creates a matrix to store results from spatial autocorrelation### PermutationsTab<-matrix(ncol=78,nrow=n.permutations) #loop function for spatial autorcorrelation analysis# for(p in 1: n.permutations){ distfact.R<-sample(distfact) PermutationsTab[p,]<-tapply(distmat.df[,4],distfact.R,mean,na.rm=T) } #calculates 95% confidence intervals# CIpermut<-apply(PermutationsTab,MARGIN=2,quantile,probs=c(0.025,0.975)) #########################Population Genetics in diveRsity############################################################## #loads R package diveRsity# libary(diveRsity) #All population diversity metrics including allelic richness and heterozygosities are generated using this command# Div<-divBasic(infile=“WisconsinGenotypes.txt”, outfile=NULL, gp=3, bootstraps=1000) ########################spatial Principal Components Analysis################################# ############################################################################################## #Loads adegenet if already installed# library(adegenet) #coordinate file# coord<-read.table(file="Coord.txt", header=T, sep="\t") ##read genepop file into adegenet## #Genepop file, must have .gen extension# gen<-read.genepop(file=”genepopfile.gen”) #spatial principal components analysis# spca <- spca(gen, xy=coord, ask=FALSE, type=1, nfpos=2, scannf=FALSE) #calls the spatially lagged scores from the sPCA# lagged<-spca$li ##barplot of eigenvalues## barplot(spca$eig,main="Eigenvalues of sPCA", col=rep(c("red","grey"),c(1,100))) ##screenplot of variance components## screeplot(spca) ##Gtest for global structure## myGtest <- global.rtest(gen$tab,spca$lw,nperm=1000) ##Gtest for local structure## myLtest <- local.rtest(gen$tab,spca$lw,nperm=1000) #######################Redundancy Analysis################################### ############################################################################# #Loads vegan if already installed# library(vegan) ##Load table with lagged scores and landscape data### B1<-read.table(file=“WisconsinVariables.txt”, header=TRUE, sep=“\t”) ###this file contains all explanatory and response variables### ##full model for RDA## full<-rda(PC1+PC2~Agriculture+River+Ecoregion+Condition(X+Y), data=B1) ##model selection## mod0<-rda(PC1+PC2 ~ 1) sel<-ordistep(mod0, scope = formula(full), perm.max = 10000, direction="both") #final selection## anova.cca(sel) ########################Spatially Lagged Regression############################ ############################################################################### #loads spdep if already installed# library(spdep) ##Load table with lagged scores and landscape data### B1<-read.table(file=“WI.txt”, header=TRUE, sep=“\t”) #coordinate file# coord<-read.table(file="Coord.txt", header=T, sep="\t") #neighbor list# lsg.nb <- dnearneigh (cbind (X, Y), 0, 1000) #nearest neighbor distances# nbd <- nbdists(lsg.nb, cbind (X, Y)) #create a geographic weight list# gl <- lapply(nbd, function(x) 1/x) #weight list final# lsg.listw <- nb2listw(lsg.nb, glist=gl) #null model# lsg1.null<- lagsarlm (PC1~1, data = B1, lsg.listw) summary(lsg1.null) #full model# lsg1.full<- lagsarlm (PC1~Ag+River+Eco, data = B1, lsg.listw) summary(lsg1.full)