Schefferville - transitions; Figures and Tables for publication

Tammy L. Elliott

Date: April 1, 2015

R version 3.1.0

# set global chunk options: 
library(knitr)
opts_chunk$set(cache=FALSE, fig.align='center')

All figures and tables are for abundance-weighted angiosperm data

Figure 1 - Plots of five cluster PCoAs and positions on grid

#Plot this relationship
#dev.new(width=11.8, height=8)
par(mfrow=c(3,3))

vasc.env.5<-cmdscale(site.env.dist, k = 54, eig = FALSE, add = FALSE, x.ret = FALSE)
plot(vasc.env.5, type="n", xlab="", ylab="", main="", cex=1, pch=16, col="black", xaxt='n', yaxt="n")
ordispider(vasc.env.5, vasc.env.clust.5$cluster, col="grey40", label=TRUE)
ordiellipse(vasc.env.5, vasc.env.clust.5$cluster, col="grey70")
box(lwd=1)

#Now do a K-means clustering analysis of these distances
ccas.env<-cascadeKM((site.env.dist),2,15)
plot(ccas.env,sortq=TRUE)

cols <- myColorRamp(c("black", "gray25", "gray60", "gray85", "white"), as.numeric(vasc.env.clust.5$cluster))
colfunc <- colorRampPalette(c("black", "gray25", "gray60", "gray85", "white"))
#par(mar=c(5.1,4.1,4.1,2.1))
plot(site.trunc$Latitude..WGS.84.datum., site.trunc$Longitude..WGS.84.datum., xlab="Latitude", ylab="Longitude",
bg=cols, pch=21, col="black",cex=1.75, cex.lab=1.5, cex.axis=1, main="", cex.main=1.25)
#xl <- 1
#yb <- -2.3
#xr <- 1.5
#yt <- 2.7
legend("bottomright", c("Group 1","Group 2","Group 3", "Group 4", "Group 5"), pch=21,
    pt.bg=c("black", "gray25", "gray60", "gray85", "white"),col="black",  cex=0.75,  bty="n")

abd.angio.abd.5<-cmdscale(angio.dist.abd, k = 100, eig = FALSE, add = FALSE, x.ret = FALSE)
plot(abd.angio.abd.5, type="n", xlab="", ylab="", main="", cex=1, pch=16, col="black", cex.main=0.9, xaxt='n', yaxt="n")
ordispider(abd.angio.abd.5, abd.angio.abd.k.clust.5$cluster, col="grey40", label=TRUE)
ordiellipse(abd.angio.abd.5, abd.angio.abd.k.clust.5$cluster, col="grey70")
box(lwd=1)

# Suggestions for optimal K-means clusters 
ccas.dist<-cascadeKM((angio.dist.abd),2,15)
plot(ccas.dist,sortq=TRUE)

myColorRamp <- function(colors, values) {
    v <- (values - min(values))/diff(range(values))
    x <- colorRamp(colors)(v)
    rgb(x[,1], x[,2], x[,3], maxColorValue = 255)
}

cols <- myColorRamp(c("black", "gray25", "gray60", "gray85", "white"), as.numeric(abd.angio.abd.k.clust.5$cluster))
colfunc <- colorRampPalette(c("black", "gray25", "gray60", "gray85", "white"))
par(mar=c(5.1,4.1,4.1,2.1))
plot(site.trunc$Latitude..WGS.84.datum., site.trunc$Longitude..WGS.84.datum., xlab="Latitude", ylab="Longitude",
bg=cols, pch=21, col="black",cex=1.75, cex.lab=1.5, cex.axis=1, main="", cex.main=1.25)
xl <- 1
yb <- -2.3
xr <- 1.5
yt <- 2.7

legend("bottomright", c("Group 1","Group 2","Group 3", "Group 4", "Group 5"), pch=21,
    pt.bg=c("black", "gray25", "gray60", "gray85", "white"),col="black",  cex=0.75,  bty="n")


angio.abd.k.5<-cmdscale(abd.angio.dist, k = 75, eig = FALSE, add = FALSE, x.ret = FALSE)
plot(angio.abd.k.5, type="n", xlab="", ylab="", main="", cex=1, pch=16, col="black", cex.main=0.9, xaxt='n', yaxt="n")
ordispider(angio.abd.k.5, abd.angio.dist.clust.five$cluster, col="grey40", label=TRUE)
ordiellipse(angio.abd.k.5, abd.angio.dist.clust.five$cluster, col="grey70")
box(lwd=1)

ccas.dist.angio<-cascadeKM((abd.angio.dist),2,15)
plot(ccas.dist.angio,sortq=TRUE)

myColorRamp <- function(colors, values) {
    v <- (values - min(values))/diff(range(values))
    x <- colorRamp(colors)(v)
    rgb(x[,1], x[,2], x[,3], maxColorValue = 255)
}

cols <- myColorRamp(c("black", "gray25", "gray60", "gray85", "white"), as.numeric(abd.angio.dist.clust.five$cluster))
colfunc <- colorRampPalette(c("black", "gray25", "gray60", "gray85", "white"))
par(mar=c(5.1,4.1,4.1,2.1))
plot(site.trunc$Latitude..WGS.84.datum., site.trunc$Longitude..WGS.84.datum., xlab="Latitude", ylab="Longitude",
bg=cols, pch=21, col="black",cex=1.75, cex.lab=1.5, cex.axis=1, main="", cex.main=1.25)
xl <- 1
yb <- -2.3
xr <- 1.5
yt <- 2.7

legend("bottomright", c("Group 1","Group 2","Group 3", "Group 4", "Group 5"), pch=21,
    pt.bg=c("black", "gray25", "gray60", "gray85", "white"),col="black",  cex=0.75,  bty="n")

textClick.bold<-function(express,col="black",font=2,cex=NULL){
par(mar = rep(0, 4),xpd=NA)
text(locator(1),express,col=col,cex=cex)
}

#textClick.bold("(A)", cex=1.75)
#textClick.bold("(B)", cex=1.75)
#textClick.bold("(C)", cex=1.75)
#textClick.bold("(D)", cex=1.75)
#textClick.bold("(E)", cex=1.75)
#textClick.bold("(F)", cex=1.75)
#textClick.bold("(G)", cex=1.75)
#textClick.bold("(H)", cex=1.75)
#textClick.bold("(I)", cex=1.75)

Compare Sums-of-squares for three different distances

ss.five<-rbind(env.dist.ss, bd.ss, ph.bd.ss)
rownames(ss.five)<-c("Gower", "Beta diversity", "Ph. beta diversity")
colnames(ss.five)<-c("Total SS", "Total Within SS", "Total Within / Total SS", "Between SS", "Between / Total SS")
ss.five

##                      Total SS Total Within SS Total Within / Total SS
## Gower                586.5044          86.893               0.1481540
## Beta diversity     35774.1975       11920.700               0.3332206
## Ph. beta diversity 35238.7403        7239.495               0.2054414
##                    Between SS Between / Total SS
## Gower                499.6114          0.8518460
## Beta diversity     23853.4979          0.6667794
## Ph. beta diversity 27999.2458          0.7945586

Partition comparison using Rand procedure

compPart.5.angio.abd<-rbind(env.dist.compPart.5,bd.abd.angio.compPart.5, phy.bd.abd.angio.compPart.5)
rownames(compPart.5.angio.abd)<-c("Env.dist",  "Bdiv.angio", "Phy.Bdiv.angio") 
colnames(compPart.5.angio.abd)<-c("Env.dist",  "Bdiv.angio", "Phy.Bdiv.angio")  
compPart.5.angio.abd

##                Env.dist Bdiv.angio Phy.Bdiv.angio
## Env.dist           1.00       0.70           0.67
## Bdiv.angio         0.70       1.00           0.75
## Phy.Bdiv.angio     0.67       0.75           1.00

Parition comparison using modified Rand procedure

compPart.5.angio.abd.crand<-rbind(env.dist.compPart.5.crand,bd.abd.angio.compPart.5.crand, phy.bd.abd.angio.compPart.5.crand)
rownames(compPart.5.angio.abd.crand)<-c("Env.dist",  "Bdiv.angio", "Phy.Bdiv.angio") 
colnames(compPart.5.angio.abd.crand)<-c("Env.dist",  "Bdiv.angio", "Phy.Bdiv.angio")  
compPart.5.angio.abd.crand

##                Env.dist Bdiv.angio Phy.Bdiv.angio
## Env.dist           1.00       0.17           0.11
## Bdiv.angio         0.17       1.00           0.25
## Phy.Bdiv.angio     0.11       0.25           1.00

Parition comparison using Nowak procedure

compPart.5.angio.abd.nowak<-rbind(env.dist.compPart.5.nowak,bd.abd.angio.compPart.5.nowak, phy.bd.abd.angio.compPart.5.nowak)
rownames(compPart.5.angio.abd.nowak)<-c("Env.dist",  "Bdiv.angio", "Phy.Bdiv.angio") 
colnames(compPart.5.angio.abd.nowak)<-c("Env.dist",  "Bdiv.angio", "Phy.Bdiv.angio")  
compPart.5.angio.abd.nowak

##                Env.dist Bdiv.angio Phy.Bdiv.angio
## Env.dist           1.00       0.36           0.33
## Bdiv.angio         0.36       1.00           0.45
## Phy.Bdiv.angio     0.33       0.45           1.00

Compare beta diversity within and among partitions for angiosperms only

#dev.new(width=5.9, height=4)
par(mai=c(0.50,0.85,0.3,0.5))
plot(abd.sp.line.five.angio[1,], axes=FALSE, col="black", pch=16, cex=1.25, ylab="Beta diversity", xlab="",las=1, cex.axis=1, cex.lab=1.2,
  ylim=c(0.35,1),cex.main=0.85, bty="c")
axis(1,1:5,labels=c("Group 1", "Group 2", "Group 3", "Group 4", "Group 5"))
axis(2)
points(abd.ph.line.five.angio[1,], col="gray70", pch=16, cex=1.25)
abline(h=sp.angio.part.among$pCqN, lwd=2, lty=2)
abline(h=angio.part.among$pCqN, lwd=2, lty=2, col="gray70")
box(bty="l", lwd=3)
legend("bottomright", c("Beta diversity (among)", "Beta diversity (within)", "Phylogenetic beta diversity (among)", "Phylogenetic beta diversity (within)"), col = c("black","black", "gray70", "gray70"), cex=0.65,
      lty = c(2, 0,2,0),lwd=c(2,2,2,2), pch = c(NA, 16, NA, 16), bg = "white", bty="n")

Differences between within and among values for beta diversity and phylogenetic beta diversity

# Differences between among and within values
# Species-level differences
sp.among.within.diff<-(sp.angio.part.among$pCqN-abd.sp.part.within.1.5.angio$pCqN)+(sp.angio.part.among$pCqN-abd.sp.part.within.2.5.angio$pCqN)+(sp.angio.part.among$pCqN-abd.sp.part.within.3.5.angio$pCqN)+(sp.angio.part.among$pCqN-abd.sp.part.within.4.5.angio$pCqN)+
(sp.angio.part.among$pCqN-abd.sp.part.within.5.5.angio$pCqN) 
sp.among.within.diff

## [1] 0.2346428

# Phylogenetic differences
ph.among.within.diff<-(angio.part.among$pCqN-abd.ph.part.within.1.5.angio$pCqN)+(angio.part.among$pCqN-abd.ph.part.within.2.5.angio$pCqN)+(angio.part.among$pCqN-abd.ph.part.within.3.5.angio$pCqN)+(angio.part.among$pCqN-abd.ph.part.within.4.5.angio$pCqN)+(angio.part.among$pCqN-abd.ph.part.within.5.5.angio$pCqN)
ph.among.within.diff

## [1] 0.1524931

Species contribution to beta diversity for three different clustering methods

#dev.new(fig.width=11.8, fig.height=8)
par(mai=c(0.1,0.1,0.15,0.1))

  plot(angio.one.tree,x.lim=340,show.tip.label = FALSE,edge.width=2)

    nspecies <- length(angio.one.tree$tip.label)
    segments(rep(210,nspecies),1:nspecies,rep(210,nspecies)+(250*env.dist.ord), 1:nspecies,lwd=4, col="grey60")
    segments(rep(250,nspecies),1:nspecies,rep(250,nspecies)+(250*abd.angio.ord), 1:nspecies,lwd=4, col="grey10")
    segments(rep(300,nspecies),1:nspecies,rep(300,nspecies)+(250*ph.abd.angio.ord), 1:nspecies,lwd=4, col="grey40")

#textClick.bold("(A)", cex=1.75)
#textClick.bold("(B)", cex=1.75)
#textClick.bold("(C)", cex=1.75)

Calculate Pearson’s correlation coefficients for different SCBD

# Pearson's correlation coefficients
cor.env.dist.bd<-cor(env.dist.ord, abd.angio.ord)
cor.env.dist.ph.bd<-cor(env.dist.ord, ph.abd.angio.ord)
cor.bd.ph.bd<-cor(abd.angio.ord, ph.abd.angio.ord)

# Correlation between Gower's distance and beta diversity
cor.env.dist.bd

## [1] 0.7420528

# Correlation between Gower's distance and phylogenetic beta diversity
cor.env.dist.ph.bd

## [1] 0.7854276

# Correlation between beta diversity and phylogenetic beta diversity
cor.bd.ph.bd

## [1] 0.9414736